[R-SIG-Finance] Option valuation for arbitrary distribution using monte carlo simulation

[Joachim Breit]

I gave it a try. See the (ugly) code below. Comments are welcome. There 
are still a few issues to solve:

- How can one do the Monte Carlo simulation faster/more efficient?
- Can somebody explain Winston's daily risk free rate?
- Is there a better optimization algorithm for the utility function than 
genoud?

----------------------------------------------
library("fImport")
library('fOptions')
library("rgenoud")
msft = 
yahooImport(query="MSFT&a=00&b=4&c=1999&d=03&e=9&f=1999&ignore=.csv",file = 
"msft", source = "http://ichart.finance.yahoo.com/table.csv?", save = TRUE)

cl = msft at data[,"Close"]   ### unadjusted close prices
cl[10:67] = cl[10:67]/2    ### a 2:1 split on 1999-03-29

ret = cl[1:66] /cl[2:67] - 1   ### the returns
logret = log(1+ret)            ### log returns
vola = sd(logret) * sqrt(252)	###  0.4091751
bsm = GBSOption(TypeFlag = "p", S = 90, X = 80, Time = 63/252, r = 
log(1+0.08), b = log(1+0.08), sigma = vola)@price  #### 2.577398 (just 
to compare)

rf = .000317   ## don't know where this comes from in Winston's paper;
                ## should be 1.08^(1/252)-1 = 0.0003054476, right?

meanutil <- function (alpha) {
x = log(alpha * (1 + ret) + (1 - alpha) * (1 + rf))  ### portf. utility
x[which(is.nan(x))] = -100000    ### correct for ln(x<0) undefined
x = mean(x)
return(x)
}

optalpha = genoud(meanutil, nvars = 1, max = TRUE, pop.size = 1000)$par
### the optimal portfolio under utility function ln(x)

portret = optalpha * (1 + ret) + (1 - optalpha) * (1 + rf)
## the optimal portfolio returns (in the sense of: daily end wealth)

num = (1 / length(ret)) * (1 / portret)
rnprob = num/sum(num)  ### the risk neutral probabilities

ov=0           #### monte carlo simulation
for (i in 1:100000) {
ov[i] = max(80 - 90 * prod(1+ sample(ret, size=63, replace = TRUE, prob 
= rnprob)),   0  )
}

mov = 1/1.08^(63/252) * mean(ov)     ### 2.512722 = option value


Am 25.11.2011 09:48, schrieb Charles Ward:
> Luenberger discusses the task of creating risk-neutral probabilities
> from any set of observed returns. I cannot find the specific reference
> in his book (Investment Science) but there is a discussion paper by
> Winston (google) that applies his method using Excel and @Risk to a
> sample of Microsoft returns. I think it would be relatively simple to
> program it in R.
>
> Charles Ward