**Overview**

The

*Black-Scholes-Merton*set of stochastic partial differential equations is the most well-known formula to price option according the price, volatility and time to expiration of the underlying security. Functional languages such as Scala allows the developer to introduce and test different volatility and time decay modeling formula.

**Note**: For the sake of readability of the implementation of algorithms, all non-essential code such as error checking, comments, exception, validation of class and method arguments, scoping qualifiers or import is omitted**Theory**

Black-Scholes approach assumes that the price of the underlying security follows a geometric Brownian distribution known as the Wiener process. If we note

*S(t)*the price of security (stock, bond, future),*V(t)*the price of the option over time,*r*the annualized compound risk-free interest rate or rate of return, and*sigma*the volatility of the price of the underlying security S then the Black-Scholes equation can be described as \[\frac{\partial V}{\partial t} + \frac{\sigma ^{2}S^{2}}{2}\frac{\partial^2 V}{\partial S^2} + r (S\frac{\partial V}{\partial S}- V) = 0\] The resolution of the partial differential equations produce the following solution for a call option: \[C = N(d)S - N(d - \sigma \sqrt{T-t})Ke^{-r(T-t))}\,\,with\\N(x)=\frac{1}{\sqrt{2 \pi }}\int_{-\infty }^{x}e^{-\frac{x^{2}}{2}}dt\,\,\,Gauss\,distribution\\d = \frac{1}{\sigma \sqrt{T-t}}[log(\frac{S}{K}) + (r+\frac{\sigma ^{2}}{2})(T-t)]\] K is the strike price for the underlying security. The price of a put option is computed as \[P=N(-d + \sigma \sqrt{T-t})Ke^{-r(T-t))} - N(-d)S\] As a reminder, a financial*call option*gives a buyer the right for a premium fee but not the obligation to acquire an agreed amount of assets*S*(securities, bonds, commodity) from a seller at a future date*T*for a certain predefined price (strike price:*K*). The seller is obligated to sell the commodity should the buyer so decide. A*put option*gives the buyer the right but not the obligation to sell an asset by a future data at a predefined (Strike) price.**Implementation**

The computation of the price of an option is CPU intensive. The Scala code snippet below implements the pricing of a call option for a security

*S*.1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 | import org.apache.commons.math3.analysis.function.Gaussian def callPrice( S: Double, X: Double, r: Double, sigma: Double, T: Double): Double = { val sigmaT = sigma * Math.sqrt(T) val d1 = (Math.log(S/X) + (r*T + 0.5*sigmaT*sigmaT) /sigmaT val d2 = d1 - sigmaT val gauss = new Gaussian S* gauss.value(d1 - X*Math.exp(-r*T)*gauss.value(d2) } |

However, professional traders in commodities and futures for instance are more interested simulating the impact of volatility, sigma, time decay t -T or price of the commodity on the price of the option using Monte Carlo simulation. In this use case, the computation of the entire Black-Scholes formula is not necessary. The objective is to minimize the cost of computation by implementing the algorithm into a workflow. Therefore the steps in the workflow that are not affected by the change in value of a single variable are not executed. For instance a change in the price of the underlying security, S requires the computation of a single workflow component: \[S=> log(S/X)\]. The Black-Scholes algorithm is broken down into 5 internal computational steps:.

val bs1F = (x: Double, y: Double) => Math.log(x/y) val bs2F = (x: Double, y: Double) => x*y val bs3F = (x: Double, y: Double) => 0.5*x*x*y val bs4F = (x: Double, y: Double) => x*Math.sqrt(y) val bs5F = (x: Double, y: Double, z: Double) => - x* Math.exp(-y*z)

The 5 computational steps are used to maintain intermediate results so the entire Black-Scholes formula does not have to be computed from scratch each time one of the parameters

The computational state is completely described by

*S*,*X*,*r*,*sigma*or*T*is changed.The computational state is completely described by

- Parameters (or coefficients) {S, X, r, sigma, T
- Internals or intermediate results {bs1 to bs5}

class Coefficients( val S: Double, val X: Double, val r: Double, val sigma: Double, val T: Double)

Next, we need to wrap the computation and update of intermediate results of the Black-Scholes formula.

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 | class Internals( val bs1: Double, val bs2: Double, val bs3: Double, val bs4: Double, val bs5: Double) { lazy val d1: Double = (bs1 + bs2 + bs3)/bs4 def setBs1(coefs: Coefficients): Internals = new Internals(bs1, bs2F(coefs.r, coefs.T), bs3, bs4, bs5) def setBs2Bs5(coefs: Coefficients): Internals = new Internals(bs1, bs2F(coefs.r, coefs.T), bs3, bs4, bs5F(coefs.X, coefs.r, coefs.T)) def setBs3Bs5(coefs: Coefficients): Internals = new Internals(bs1, bs2, bs3F(coefs.sigma, coefs.T), bs4, bs5F(coefs.X, coefs.r, coefs.T)) } |

As mentioned earlier, the computation state is defined by the Black-Scholes coefficients and intermediate results values.

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 | class State(val coefs: Coefficients, val internals: Internals) { def this(coefs: Coefficients) = this(coefs, initInternals(coefs)) def this(S: Double, X: Double, r: Double, sigma: Double, T: Double) = this(new Coefficients(S, X, r, sigma, T)) def setS(S: Double): State = { val newCoefs = new Coefficients(S, coefs.X, coefs.r, coefs.sigma, coefs.T) new State(newCoefs, internals.setBs1(newCoefs)) } def setr(r: Double): State = { val newCoefs = new Coefficients(coefs.S, coefs.X, r, coefs.sigma, coefs.T) new State(newCoefs, internals.setBs2Bs5(newCoefs)) } def setSigma(sigma: Double): State = { val newCoefs = new Coefficients(coefs.S, coefs.X, coefs.r, sigma, coefs.T) new State(newCoefs, internals.setBs3Bs5(newCoefs)) } lazy val call: Double = { import org.apache.commons.math3.analysis.function.Gaussian val gauss = new Gaussian coefs.S * gauss.value(internals.d1) - internals.bs5*gauss.value(internals.d1 -internals.bs4) } } |

The method

The price of the call

The update of the state of the Black-Scholes computation is rather straight forward.

*setS*,*setr*,*setSigma*update one of the Black-Scholes coefficients that automatically triggers a selective recomputation of some of the intermediate results.The price of the call

*call*is a lazy value to it is computed only once (immutable state)The update of the state of the Black-Scholes computation is rather straight forward.

1 2 3 4 5 6 7 8 9 | // Initial computation/state of Black-Scholes val state0 = new State(0.4, 0.6, 0.1, 2.0, 1.4) // Compute the price of a call val callValue = state0.call // Update the value of sigma val state1 = state0.setSigma(0.14) val newCallValue = state1.call |

**Monte-Carlo Simulation**

The Black-Scholes formula can be validated through a Monte-Carlo simulation, using a Gaussian distribution to represent the stochastic component of the formula.

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 | def monteCarlo(state: State, maxIters: Int) : Option[Double] = { val coefs = state.coefs val sigmaT = coefs.sigma * Math.sqrt(coefs.T) val alpha1 = Math.exp(coefs.r*coefs.T) val alpha2 = coefs.S*Math.exp(- 0.5*sigmaT*sigmaT) val theta = alpha1*alpha2 val delta = Math.exp( sigmaT) val totalValues = Range(0, maxIters).scanLeft(0.0)((_, n) => theta*Math.pow(delta, Random.nextGaussian) ).takeWhile( _ > coefs.X) if( totalValues.size < maxIters) Some(alpha1 *totalValues.last / maxIters) else None } |

**References**

*Real options analysis: Tools and techniques for valuing strategic investments - Johnathan Mun - John Wiley & Sons 2002*

github.com/prnicolas

I like the introduction to the BS .. easy to follow. I have a couple of questions however.

ReplyDelete1. Is there an alternative implementation to boost performance? Parallel collections in Scala?

2. What is the purpose of the Monte Carlo simulation? Generate stochastic data? The last section of the post is a bit confusing..

No Sweat!