Skip to content

optim function multi variables -newbie

6 messages · Uwe Ligges, michalseneca, Hans W Borchers

michalseneca
# 
Hi ,
can i ask if somebody know the syntax for optim function with multiple
variables:

I have function  f  (params,y,g,h,j) returnig z
y<-C(0.2,0.5,0.6)
where 
params<-c(x1,x2)
I define g,h,j and then I define params, then use
optim to structure it ..however if i use syntax

optim(params , f,y=y,g=g,h=h,j=j) to optimize params , but i got result
equalling to my params ..which is kind of wrong. 
I and if i use 
optim(params , f(params,y=y,g=g,h=h,j=j)) I get 

Error in function (par)  : could not find function "fn"

how ever none of them seems to be producing meaningful results

Please I am newbie and could not find a proper answer to that..

Help me.!!

Thanks
Uwe Ligges
#
On 27.01.2011 14:08, michalseneca wrote:
The syntax above is correct, hence you need to describe the problem you 
try to solve, ideally including the full code for f and some example data.

Uwe Ligges
>
michalseneca
# 
Hi,
Basically what I am trying is to rewrite matlab code into R ...This is code
for famous SABR model Calibration.
I did most of the code up to optimalization issue.
In the attachment you can find matlab and my unfinished R code..This is very
helpful code...Basically I am trying to achieve the same results.


Thank you very much for your help..This might help me much with my thesis


Best regards 
Michal 


Michal http://r.789695.n4.nabble.com/file/n3242469/optimtest optimtest 
http://r.789695.n4.nabble.com/file/n3242469/Hagan_Figure_33.zip
Hagan_Figure_33.zip
Hans W Borchers
# 
michalseneca <michalseneca <at> gmail.com> writes:
Apparently, you are resetting the input 'params' within the function
EstimateRhoAndVol() to c(0.5, 0.5) regardless of what the actual input is.
No wonder that 'optim' returns the starting values: it sees a constant
function.

Besides that, you have to be more careful with NAs in the SABR model.
Did you try out the function? --- Are you testing the list?

--Hans Werner
michalseneca
# 
I tried to modify the code,thanks for noticing...now i get that the function
cannot be evaluated at initial parameters.However I do not know what does
that mean..Should I try to modify parameters. I am still not sure of syntax
of function. I cannot get that optimize work correctly. I am trying to test
that..however I cannot get to a result..In matlab however that works
correctly...Did anybody of you tried to succesfully run it ?

Michal
Hans W Borchers
# 
michalseneca <michalseneca <at> gmail.com> writes:
I went through your functions and the Matlab code, below I will append my
versions. Matlab and R now return the same values, at least for a few test
cases. One obstacle, for instance, was that solve.polynom() expects the
coefficients in reverse order.

Here is an example of how to apply optim(). Whether the result is
resonable you'll have to find out for yourself:

    MktVol <- c(20.8, 19.8, 18.1, 16.1, 15.1, 14.5, 14.0, 14.3, 15.0,
                15.9, 16.2, 16.4, 16.6, 17.3, 19.1) / 100;
    MktStrike <- c(4.3, 4.8, 5.3, 5.8, 6.3, 6.8, 7.3, 7.8, 8.3, 8.8,
                   9.3, 9.8, 10.3, 10.8, 11.3) / 100;
    F <- 0.073
    ATMVol <- 0.14
    T <- 1
    b<-1
    parsbox<-c(0.7,0.5)

    f <- function(x) EstimateRhoAndVol(x, MktStrike, MktVol, ATMVol, F, T, b)
    optim(parsbox, f)
    # $par
    # [1] -0.06374008  0.60383201

Of course, there are some more open points, for instance will 'findAlpha'
always find a real, positive cubic root? etc.

--Hans Werner

### --- snip ---------------------------

library(polynom)

EstimateRhoAndVol <- function(params, MktStrike, MktVol, ATMVol, F, T, b)
{
# -------------------------------------------------------------------------
# Returns the following SABR parameters:
# r = rho
# v = vol-of-vol
# Uses ATM volatility to estimate alpha
# Required inputs:
# MktStrike = Vector of Strikes
# MktVol    = Vector of corresponding volatilities
# ATMVol = ATM volatility
# F = spot price
# T = maturity
# b = beta parameter
# -------------------------------------------------------------------------
  r <- params[1]
  v <- params[2]
  a <- findAlpha(F, F, T, ATMVol, b, r, v)
  N <- length(MktVol)

  ModelVol <- numeric(N)
  error <- numeric(N)
  # Define the model volatility and the squared error terms
  for (i in 1:N) {
    ModelVol[i] = SABRvol(a, b, r, v, F, MktStrike[i], T)
    error[i] = (ModelVol[i] - MktVol[i])^2
  }

  # Return the SSE
  y <- sum(error, na.rm=TRUE)

  # Impose the constraint that -1 <= rho <= +1
  # via a penalty on the objective function
  if (abs(r) > 1)
    y <- Inf
  return(y)
}


SABRvol <- function(a, b, r, v, F, K, T)
{
# ---------------------------------------------
# Returns the SABR volatility.
# Required inputs:
# a = alpha parameter
# b = beta parameter
# r = rho parameter
# v = vol of vol parameter
# F = spot price
# K = strike price
# T = maturity
# ---------------------------------------------

  # Separate into ATM case (simplest case) and
  # Non-ATM case (most general case)

  if (abs(F-K) <= 0.001) {  # ATM vol
	Term1 <- a/F^(1-b)
	Term2 <- ((1-b)^2/24*a^2/F^(2-2*b) + r*b*a*v/4/F^(1-b) + (2-3*r^2)*v^2/24)
	y <- Term1 * (1 + Term2*T)

  } else {  # Non-ATM vol
	FK <- F*K
	z <- v/a*(FK)^((1-b)/2)*log(F/K)
	x <- log((sqrt(1 - 2*r*z + z^2) + z - r)/(1-r))
	Term1 <- a / FK^((1-b)/2) / (1 + (1-b)^2/24*log(F/K)^2 + 
	         (1-b)^4/1920*log(F/K)^4)
	if (abs(x-z) < 1e-10) Term2 <- 1
	else                  Term2 <- z / x
	Term3 <- 1 + ((1-b)^2/24*a^2/FK^(1-b) + 
	             r*b*v*a/4/FK^((1-b)/2) + (2-3*r^2)/24*v^2)*T
	y <- Term1 * Term2 * Term3
  }
  return(y)
}


findAlpha <- function(F, K, T, ATMvol, b, r, v)
{
# -------------------------------------------------------------------------
# Finds the SABR "alpha" parameter by solving the cubic equation described
# in West (2005) "Calibration of the SABR Model in Illiquid Markets"
# The function can return multiple roots.  In that case, the program 
# eliminates negative roots, and select the smallest root among the 
# positive roots that remain.
# Required inputs:
# F = spot price
# K = strike price
# T = maturity
# ATMvol = ATM market volatility
# b = beta parameter
# r = rho parameter
# v = vol of vol parameter
# -------------------------------------------------------------------------

  # Find the coefficients of the cubic equation for alpha
  C0 <- -ATMvol * F^(1-b)
  C1 <- (1 + (2-3*r^2) * v^2 * T / 24)
  C2 <- r * b * v * T / 4 / F^(1-b)
  C3 <- (1-b)^2 * T / 24 / F^(2-2*b)

  # Return the roots of the cubic equation (multiple roots)
  AlphaVector <- solve(as.polynomial(c(C0, C1, C2, C3)))

  # Find and return the smallest positive root
  index <- which(Im(AlphaVector) == 0 & Re(AlphaVector) > 0)
  Alpha <- AlphaVector[index]
  if (length(Alpha) == 0) Alpha <- 0
  return(min(Re(Alpha)))
}