r - 在 R 编程中修复波动率曲面图的插值

Question

下面的脚本通过 quantmod 中的函数提取 yahoo 数据，然后处理周围的数据，以使用 RGL 库形成 3D 图形，附加一个 ggplot 来显示我正在尝试使用单独的线几何图形创建表面的数据。问题是 3D 图表看起来非常难看并且被切割，因为前一个月到期的点数量有限。谁能告诉我这里发生了什么，我能做些什么来解决这个问题。我需要平滑吗？然后插入每个到期行......？ volsurface http://img15.imageshack.us/img15/7338/surface.png ggplot2_smile http://img402.imageshack.us/img402/1272/volatilitysmilegoog.png

library(RQuantLib)
library(quantmod)
library(rgl)
library(akima)
library(ggplot2)
library(plyr)

GetIV <- function(type, value,
                  underlying, strike,dividendYield, riskFreeRate, maturity, volatility,
                  timeSteps=150, gridPoints=151) {

    AmericanOptionImpliedVolatility(type, value,
                                    underlying, strike,dividendYield, riskFreeRate, maturity, volatility,
                                    timeSteps=150, gridPoints=151)$impliedVol
}


GetDelta <- function(type, underlying, strike,
                     dividendYield, riskFreeRate, maturity, volatility, 
                     timeSteps=150, gridPoints=149, engine="CrankNicolson") {

    AmericanOption(type,underlying, strike, dividendYield, riskFreeRate, maturity, volatility,
                   timeSteps=150, gridPoints=149, engine="CrankNicolson")$delta
}
# set what symbol you want vol surface for
underlying <- 'GOOG'
# set what your volatility forcast or assumption is
volforcast <- .25
# Get symbols current price
underlying.price <- getQuote(underlying,what=yahooQF("Last Trade (Price Only)"))$Last

OC <- getOptionChain(underlying, NULL)
#check data
head(OC)
lputs <- lapply(OC, FUN = function(x) x$puts[grep("[A-Z]\\d{6}[CP]\\d{8}$", rownames(x$puts)), ])
head(lputs) #check for NA values, yahoo returns all NA values sometimes
puts <- do.call('rbind', lputs )
#check data
head(puts,5)

symbols <- as.vector(unlist(lapply(lputs, rownames)))
expiries <- unlist(lapply(symbols, FUN = function(x) regmatches(x=x, regexpr('[0-9]{6}', x) )))
puts$maturity <- as.numeric((as.Date(expiries, "%y%m%d") - Sys.Date())/365)

puts$IV <- mapply(GetIV, value = puts$Ask, strike = puts$Strike, maturity = puts$maturity,
                  MoreArgs= list(type='put', underlying= underlying.price,
                                 dividendYield=0, riskFreeRate = 0.01,  
                                 volatility = volforcast), SIMPLIFY=TRUE)

puts$delta <- mapply(GetDelta, strike =  puts$Strike, volatility = puts$IV,
                     maturity = puts$maturity, MoreArgs= list(type='put', 
                                                              underlying=underlying.price, dividendYield=0, 
                                                              riskFreeRate = 0.01 ), SIMPLIFY=TRUE)

# subset out itm puts
puts <- subset(puts, delta < -.09 & delta > -.5 )

expiries.formated <- format(as.Date(levels(factor(expiries)), format = '%y%m%d'), "%B %d, %Y")

fractionofyear.levels <- levels(factor(puts$maturity))

xyz <- with(puts, interp(x=maturity, y=delta*100, z=IV*100, 
                         xo=sort(unique(maturity)), extrap=FALSE ))

with(xyz, persp3d(x,y,z, col=heat.colors(length(z))[rank(z)], xlab='maturity', 
                  ylab='delta', zlab='IV', main='IV Surface'))

putsplot <- ggplot(puts, aes(delta, IV, group = factor(maturity), color = factor(maturity))) +
    labs(x = "Delta", y = "Implied Volatilty", title="Volatility Smile", color = "GooG \nExpiration") +
    scale_colour_discrete( breaks=c(fractionofyear.levels),
                           labels=c(expiries.formated)) + 
    geom_line() +
    geom_point()

putsplot

Answer 1

akima 包正是您所需要的，但我认为您需要减少 y 轴上的插值点的数量，即 delta 变量。您现在设置的方式使用默认的 40 点网格。

# No interpolation on x-axis, but uses the default 40 point grid on the y-axis
xyz <- with(puts, interp(x=maturity, y=delta*100, z=IV*100, 
            xo=sort(unique(maturity)), extrap=FALSE ))
# By setting to use less points, it will "stretch" the surface over those points.
xyz <- with(puts, interp(x=maturity, y=delta*100, z=IV*100, 
            xo=sort(unique(maturity)), 
            yo=seq(min(delta*100), max(delta*100), length = 15), extrap=FALSE ))

您可以使用 seq 函数中的长度变量来获得不同级别的平滑度。

<小时/>

我仍然不完全明白你想要什么，但也许你想通过成熟度来平滑？看起来是这样的:

# This smooths just by x.
xyz <- with(puts, interp(x=maturity, y=delta*100, z=IV*100, 
            xo=seq(min(maturity), max(maturity), length = 5), 
            , extrap=FALSE ))

with(xyz, persp3d(x,y,z, col=heat.colors(length(z))[rank(z)], xlab='maturity', 
                  ylab='delta', zlab='IV', main='IV Surface'))