Option valuation for arbitrary distribution using monte carlo simulation

Wed, Nov 23, 2011 6:03 AM

Dear Massimo,

thanks for the hint to Iacus' book. I was aware of this title already 
but could not get hold of it by now. Thank you for the code snippet as 
well. Unfortunately, it does not seem to help with my main problem: How 
to correct the original data such that Monte Carlo sampling from it 
gives me samples from a risk-neutral but NONNORMAL process. Maybe my 
description is a bit awkward. Also I may lack the math and statistics 
skills. I try do say it another way:

Let
r : original log return series
     ( with e.g. mean(r)  = -0.0001574849 and sd(r)  = 0.01479017 )
x:  risk-free rate of, say, 5%.
Snow : stock price now
time to maturity = 1 year = 252 trading days
St : stock price in 1 year
X: call option strike
E[] : expected value
C : the option price

I am looking for a corrected return series rcorr = r + z such that
E[ Snow * exp(  sum(  sample(rcorr, 252)     ) ) ] = 100 * exp(0.05)

because then I can find the call option price with
for (i in 1:bignumber) {
   St[i] = Snow * exp(  sum(   sample(rcorr, 252)     ) )
   op[i] = max(0, exp(-.05) * (St[i] - X)  )
}
C = mean(op)

I could of course use some search algorithm to find z. In fact thats 
what I did. But I thought this is so common a problem for finance folks 
that there must be some elegant R procedure for th task.

Your code seems to assume that we are dealing with normal log returns 
but this is not the case for my sample.

Do you see what I mean? Is this stuff dealt with in the Iacus book?

Greetings
Joachim

Hi Joachim , I'm not a prof but I recommend you the following book:
Option Pricing and Estimation of Financial Models with R
Stafano M. Iacus
http://www.amazon.com/Option-Pricing-Estimation-Financial-Models/dp/0470745843/ref=sr_1_48?ie=UTF8&qid=1321974343&sr=8-48

it'll explain a lot of thinks with code samples.
Some very base simple code, I started with this one and while was reading
the book I started to improve it

#define functions
require(package='fBasics',quietly=FALSE)

#This function to create innovation
#We keep innovation function autside of process because
#once computed leave it stored for next usage
innovation.norm<-function(pTrials=5000){
   # remember to add the seed
   #generate standard normal innovation
   z<- rnorm(n=pTrials,mean=0,sd=1)
   #genereta anthitetic innovation
   ant.z<--z
   #combine in the new one:
   z<- c(z,ant.z)
   #return innovation vector
   return(z)
}

#this function to create the price distribution
st.proc.norm<-function(pInnVect=z,mu=0,pSnow=15,pVol=0.15,pTradDaysToExp=20){
   #set the number of trading days in one year
   yearTradDays<-252
   #compute the interval time
   dt<-pTradDaysToExp/yearTradDays
   s.final<- rep(0,length(pInnVect))
   s.final<-pSnow*exp(pVol*sqrt(dt)*pInnVect)
   out<-c(s.final,pSnow)
   return(s.final)
}

#this function to build the option chain
price.chain<-function(pSFinal,pIntRate=0.02,pStrikeMin=10,pStrikeMax=15,pDStrike=0.5,pCalDaysToExp=30){
   #build strike chain
   K<-seq(from=pStrikeMin,to=pStrikeMax,by=pDStrike)
   #compute dicount factor
   calDays<-365
   dt<-pCalDaysToExp/calDays
   #prepare storage for chain
   prCall<-rep(x=0,length(K))
   prPut<-rep(x=0,length(K))
   se.call<-rep(x=0,length(K))
   se.put<-rep(x=0,length(K))

   for (i in (1:length(K))){
     #print(i)
     #print(K[i])
     #formula for pricing
     c<- ifelse(pSFinal>K[i],exp(-pIntRate*dt)*(pSFinal-K[i]), 0)
     p<- ifelse(pSFinal<K[i],exp(-pIntRate*dt)*(K[i]-pSFinal), 0)
     #price as average
     prCall[i]<-mean(c)
     prPut[i]<-mean(p)
     #standard error
     se.call[i]<-sd(c)/sqrt(length(pSFinal))
     se.put[i]<-sd(p)/sqrt(length(pSFinal))

     ## make a confidence interval estimate
     alpha<-0.05
     t.star<-qt(p=(1-alpha/2),df=length(pSFinal)-1,lower.tail=TRUE)
     prCall[i]+c(-1,1)*t.star*se.call[i]
     prPut[i]+c(-1,1)*t.star*se.put[i]
   }

   #build the data frame for option's chain

pr.chain<-data.frame(CALL=round(prCall,digits=3),Strike=K,PUT=round(prPut,digits=3))

   return(pr.chain)
}


##################################################################################
##  PRICING WITH NORMAL
inn.n<-innovation.norm(pTrials=500)
hist(inn.n)
#you must estimate volatility from past data
s.final.n<-st.proc.norm(pInnVect=inn.n,pSnow=13.40,pVol=0.25,pTradDaysToExp=20)
hist(s.final.n)
#you must estimate interest rate
price.chain(pSFinal=s.final.n,pIntRate=0.02,pStrikeMin=10,pStrikeMax=15,pDStrike=0.5,pCalDaysToExp=20)

Option valuation for arbitrary distribution using monte carlo simulation

Thread (15 messages)