R: Calculating IV using Black-Scholes and bisection method, loop refusing to work

Question

I have my Black-Scholes function and my bisection model for call options with data from a CSV. It appears to be getting stuck in the inner loop because it stays above the tolerance. My Black-Scholes does calculate accurately and I am using the average of bid and ask for the market price instead of the actual price of the option. After working on this for hours, maybe I am just missing something obvious.

The link to the CSV is here: http://s000.tinyupload.com/?file_id=06213890949979926112

########################################################################
#Black-Scholes-Merton Call
bsmCall <- function(S, K, M, sig, r) {
  yrTime=(M/252)
  d1 <- (log(S/K)+(r+(sig^2/2))*(yrTime))/(sig*(sqrt(yrTime)))
  d2 <- d1-sig*(sqrt(yrTime))
  C <- (S*(pnorm(d1)))-((pnorm(d2))*K*(exp(-r*yrTime)))
  return(C)
}
########################################################################

myData = read.csv("09-26-16.csv", stringsAsFactors=FALSE)    #DATA
myData <- myData[,2:24]   #omit first column

####### start bisection method of CALLS and put IV in database #######
i <- 1    # reset counter
tol <- 0.000001   #tolerance

while(i <= nrow(myData)) {
  if((myData[i,5] != 0) & (myData[i,6] != 0)) {
    volLower <- .0001    #will need to reset with each iteration
    volUpper <- 1         #will need to reset with each iteration
    volMid <- (volLower + volUpper) / 2   #will need to reset with each iteration

    while(abs(bsmCall(as.numeric(as.character(myData[i,17])),as.numeric(as.character(myData[i,1])),as.numeric(as.character(myData[i,22])),volMid,as.numeric(as.character(myData[i,23])))-(as.numeric(as.character(myData[i,5])))) >= tol) {
      if((bsmCall(as.numeric(as.character(myData[i,17])),as.numeric(as.character(myData[i,1])),as.numeric(as.character(myData[i,22])),volMid,as.numeric(as.character(myData[i,23])))-(as.numeric(as.character(myData[i,5])))) < 0) {
        volLower <- volMid
        volMid <- (volUpper + volMid)/2
      } else {
        volUpper <- volMid
        volMid <- (volLower + volMid)/2
      }
    }
    myData[i,8] <- volMid
  } else { myData[i,8] <- 0 }
  i=i+1
}

Answer 1

The problem is here:

while(abs(bsmCall(as.numeric(as.character(myData[i,17])),
                  as.numeric(as.character(myData[i,1])),
                  as.numeric(as.character(myData[i,22])),
                  volMid,
                  as.numeric(as.character(myData[i,23])))-(as.numeric(as.character(myData[i,5])))) >= tol)

You're using a while loop on a condition that, if true, is always true. It's an infinite loop. On your first row of data this problem is encountered.

How to fix this error is specific to your use case, but if you just change while to if you'll see the loop complete immediately.

You asked about the bisection method. There are a few in packages and here's another from here:

bisect <- function(fn, lower, upper, tol=1.e-07, ...) {
f.lo <- fn(lower, ...)
f.hi <- fn(upper, ...)
feval <- 2

if (f.lo * f.hi > 0) stop("Root is not bracketed in the specified interval
\n")
chg <- upper - lower

while (abs(chg) > tol) {
        x.new <- (lower + upper) / 2
        f.new <- fn(x.new, ...)
        if (abs(f.new) <= tol) break
        if (f.lo * f.new < 0) upper <- x.new
        if (f.hi * f.new < 0) lower <- x.new
        chg <- upper - lower
        feval <- feval + 1
}
list(x = x.new, value = f.new, fevals=feval)
}

# An example
fn1 <- function(x, a) {
exp(-x) - a*x
}

bisect(fn1, 0, 2, a=1)

bisect(fn1, 0, 2, a=2)

Recursive version:

bisectMatt <- function(fn, lo, hi, tol = 1e-7, ...) {

    flo <- fn(lo, ...)
    fhi <- fn(hi, ...)

    if(flo * fhi > 0)
        stop("root is not bracketed by lo and hi")

    mid <- (lo + hi) / 2
    fmid <- fn(mid, ...)
    if(abs(fmid) <= tol || abs(hi-lo) <= tol)
        return(mid)


    if(fmid * fhi > 0)
        return(bisectMatt(fn, lo, mid, tol, ...))

    return(bisectMatt(fn, mid, hi, tol, ...))
}

Answer 2

Jeez, this is my 3rd edit so far...

Lets reconstruct the while loop when i=1 and print the volMid - the of the only part of the while condition that is updating after each iteration

i <- 1
volLower <- .0001    #will need to reset with each iteration
volUpper <- 1         #will need to reset with each iteration
volMid <- (volLower + volUpper) / 2   #will need to reset with each iteration

j <- 1
while(abs(bsmCall(myData[i,17], myData[i,1], myData[i,22],volMid,myData[i,23])-myData[i,5]) >= tol & j < 30) {
  if(bsmCall(myData[i,17], myData[i,1], myData[i,22],volMid,myData[i,23])-myData[i,5] < 0) {
volLower <- volMid
volMid <- (volUpper + volMid)/2
  } else {
    print("pos")
    volUpper <- volMid
    volMid <- (volLower + volMid)/2
  }
  j <- j + 1
  print(volMid)
}

Result:

#[1] 0.750025
#[1] 0.8750125
#[1] 0.9375062
#[1] 0.9687531
#[1] 0.9843766
#[1] 0.9921883
#[1] 0.9960941
#[1] 0.9980471
#[1] 0.9990235
#[1] 0.9995118
#[1] 0.9997559
#[1] 0.9998779
#[1] 0.999939
#[1] 0.9999695
#[1] 0.9999847
#[1] 0.9999924
#[1] 0.9999962
#[1] 0.9999981
#[1] 0.999999
#[1] 0.9999995
#[1] 0.9999998
#[1] 0.9999999
#[1] 0.9999999
#[1] 1
#[1] 1
#[1] 1
#[1] 1
#[1] 1
#[1] 1

volMid converges to 1 after less than 30 iterations, and from there on out, it's stuck.