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

[Patrick Burns]

Joachim,

I'm guessing that you could get
some better speed by doing the
'sample' call all at once, do it
with 'logret' instead of 'ret'
and sum them efficiently.  Something
like:

sampmat <- array(sample(logret,
    size=63 * 1e5, replace=TRUE,
    prob=rnprob), c(63, 1e5))

samplret <- colSums(sampmat)

# then transform to simple returns
# and use 'pmax' instead of 'max'

If the matrix is too large for your
memory, you could do it in a few
blocks.

On 23/12/2011 14:59, Joachim Breit wrote:
> 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
>
> _______________________________________________
> R-SIG-Finance at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-finance
> -- Subscriber-posting only. If you want to post, subscribe first.
> -- Also note that this is not the r-help list where general R questions
> should go.
>

-- 
Patrick Burns
patrick at burns-stat.com
http://www.burns-stat.com
http://www.portfolioprobe.com/blog
twitter: @portfolioprobe