modification Carr-Madan

This commit is contained in:
François Pelletier 2014-05-27 23:06:48 -04:00
parent 92b291e2b4
commit 2dc91e644c
2 changed files with 31 additions and 32 deletions

View file

@ -9,6 +9,7 @@
#' Call price using the Carr-Madan damping parameter and FFT
#' @param strikeprice Vector of strike prices, relative to a unit stock price
#' @param char.fn Characteristic function of the log-price process
#' @param param Characteristic function parameters
#' @param eval.time Evaluation time
#' @param expiry.time Expiry time
#' @param rate Continuously compounded interest rate (force of interest)
@ -18,11 +19,9 @@
#' @return A European call option price vector
#' @export callCarrMadan
#' @author Francois Pelletier
callCarrMadan <- function(strikeprice,char.fn,eval.time,expiry.time,rate,alpha,
...,fft.control=list(N=2^14,eta=.1))
callCarrMadan <- function(strikeprice,char.fn,param,eval.time,expiry.time,rate,alpha,
...,fft.control=list(N=2^10,eta=.1))
{
# Determine moneyness
moneyness <- strikeprice < 1
# Discretization step for Fourier transform
lambda <- lambda <- (2*pi) / (fft.control$N*fft.control$eta)
# Evaluation points of the damped characteristic function of the call option log-price
@ -32,26 +31,25 @@ callCarrMadan <- function(strikeprice,char.fn,eval.time,expiry.time,rate,alpha,
# Vector of indices
jvec <- 1:fft.control$N
# Simpson's hypothesis
simpsonh <- ((dampedcfcallCarrMadan(u,char.fn,eval.time,expiry.time,rate,alpha,...,moneyness)*
exp(1i*u*b)*fft.control$eta)/3)*
(3+(-1)^jvec+((jvec-1)==0))
simpsonh_money <- ((dampedcfcallCarrMadan(u,char.fn,param,eval.time,expiry.time,rate,alpha,moneyness=TRUE)*
exp(1i*u*b)*fft.control$eta)/3)*(3+(-1)^jvec+((jvec-1)==0))
simpsonh_nomoney <- ((dampedcfcallCarrMadan(u,char.fn,param,eval.time,expiry.time,rate,alpha,moneyness=FALSE)*
exp(1i*u*b)*fft.control$eta)/3)*(3+(-1)^jvec+((jvec-1)==0))
# Log-price vector
ku <- seq(-b,(fft.control$N-1)*lambda-b,lambda)
# Log-price of the call option vector
if(moneyness)
{
callvec <- Re((exp(-alpha*ku)*fft(simpsonh))/pi)
}
else
{
callvec <- fft(simpsonh)/(sinh(alpha*ku)*pi)
}
# Price vector
Ku <- exp(ku)
# Log-price of the call option vector
callvec_money <- Re((exp(-alpha*ku)*fft(simpsonh_money))/pi)
callvec_nomoney <- Re(fft(simpsonh_nomoney)/(sinh(alpha*ku)*pi))
callvec <- callvec_money * (Ku < 1) + callvec_nomoney * (Ku >= 1)
# Index to select subset of prices in the strikeprice vector
Kindex <- Ku>=(min(strikeprice)-1) & Ku<=(max(strikeprice)+1)
Kindex <- Ku>=(min(strikeprice)-1) & Ku<=(max(strikeprice)+1) & abs(callvec) < 10^9 & abs(callvec) > 10^-9
# We use a smooth spline to get the prices for the strikeprice vector
sp0 <- smooth.spline(x=Ku[indice],y=callvec[indice])
sp0 <- smooth.spline(x=round(Ku[Kindex],digits=10),y=round(callvec[Kindex],digits=10))
predict(sp0,strikeprice)$y
}

View file

@ -7,36 +7,37 @@
#' Damped characteristic function of the call option log-price
#' @param u Transform variate
#' @param u Transform variate (vector)
#' @param char.fn Characteristic function of the log-price process
#' @param param Characteristic function parameters
#' @param eval.time Evaluation time
#' @param expiry.time Expiry time
#' @param rate Continuously compounded interest rate (force of interest)
#' @param alpha Damping parameter
#' @param ... Parameters of the characteristic function
#' @param moneyness Boolean for moneyness of call option
#' (TRUE if strike price is lower than stock price)
#' @return Characteristic function value
#' @export dampedcfcallCarrMadan
#' @author Francois Pelletier
dampedcfcallCarrMadan <- function(u,char.fn,eval.time,expiry.time,rate,alpha,...,moneyness=TRUE)
dampedcfcallCarrMadan <- function(u,char.fn,param,eval.time,expiry.time,rate,alpha,moneyness,...)
{
if(moneyness)
auxiliairyf <- function(u,char.fn,param,eval.time,expiry.time,rate,alpha,...)
{
exp(-rate*(expiry.time-eval.time))*
(1/(1+1i*u)-exp(rate*(expiry.time-eval.time))/
(1i*u)-char.fn(u-1i,param,eval.time,expiry.time,...)/(u^2-1i*u))
}
if (moneyness)
{
return(exp(-rate*(expiry.time-eval.time))*
char.fn(u-1i*(alpha+1),expiry.time-eval.time,...) /
(alpha^2+alpha-u^2+1i*u*(2*alpha+1)))
char.fn(u-1i*(alpha+1),param,eval.time,expiry.time,...) /
(alpha^2+alpha-u^2+1i*u*(2*alpha+1)))
}
else
{
auxiliairyf <- function(u,char.fn,eval.time,expiry.time,rate,alpha,...)
{
exp(-rate*(expiry.time-eval.time))*
(1/(1+1i*u)-exp(rate*(expiry.time-eval.time))/
(1i*u)-char.fn(u-1i,expiry.time-eval.time,...)/(u^2-1i*u))
}
return((auxiliairyf(u-1i*alpha,char.fn,eval.time,expiry.time,rate,alpha,...)-
auxiliairyf(u+1i*alpha,char.fn,eval.time,expiry.time,rate,alpha,...))/2)
return((auxiliairyf(u-1i*alpha,char.fn,param,eval.time,expiry.time,rate,alpha,...)-
auxiliairyf(u+1i*alpha,char.fn,param,eval.time,expiry.time,rate,alpha,...))/2)
}
}