Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Monte Carlo pricing of Exotic Options in F#
open MathNet.Numerics.Distributions
open MathNet.Numerics.Statistics
let callPayoff strike price = max (price - strike) 0.0
let europeanPayoff payoff assetPath = assetPath |> Seq.last |> payoff
let europeanCallPayoff strike assetPath = assetPath |> europeanPayoff (callPayoff strike)
let asianArithmeticMeanCallPayoff strike (assetPath:seq<float>) = assetPath.Mean() |> callPayoff strike
let upAndOutPayoff payoff upperBarrier assetPath =
if assetPath |> Seq.exists (fun x -> x > upperBarrier) then 0.0 else payoff assetPath
let upAndOutEuropeanCallPayoff strike = upAndOutPayoff (europeanCallPayoff strike)
let doubleBarrierPayoff payoff upperBarrier lowerBarrier assetPath =
if assetPath |> Seq.exists (fun x -> x > upperBarrier || x < lowerBarrier) then 0.0 else payoff assetPath
let doubleBarrierEuropeanCallPayoff strike = doubleBarrierPayoff (europeanCallPayoff strike)
let getAssetPath S0 r deltaT sigma (normal:Normal) numSamples =
Seq.unfold (fun S -> let sNew = (S * exp(((r - (0.5 * sigma * sigma)) * deltaT) + (sigma * sqrt(deltaT) * normal.Sample())))
Some(sNew, sNew)) S0
|> Seq.take numSamples
let phi x =
let normal = new Normal()
normal.CumulativeDistribution(x)
// Exact option pricing:
let priceEuropeanCall S0 strike r T sigma =
let discountedK = strike * exp(-r * T)
let totalVolatility = sigma * sqrt(T)
let d_minus = log(S0 / discountedK) / totalVolatility
let d_plus = d_minus + (0.5 * totalVolatility)
let d_minus2 = d_minus - (0.5 * totalVolatility)
(S0 * phi(d_plus)) - (discountedK * phi(d_minus2))
let computeD S0 strike r sigma T =
let d1 = (log (S0 / strike) + ((r + (sigma * sigma / 2.0)) * T)) / (sigma * sqrt(T))
let d0 = d1 - (sigma * sqrt(T))
(d0, d1)
let priceBarrierUpAndOutCall S0 strike B r T sigma =
let normal = new Normal()
let mu = r - (0.5 * sigma * sigma)
let nu = 2.0 * (mu / (sigma * sigma))
let aux1 = (B / S0) ** (nu + 2.0)
let aux2 = (B / S0) ** nu
let d = computeD S0 strike r sigma T
let price = S0 * normal.CumulativeDistribution(snd d) - (strike * exp(-r * T) * normal.CumulativeDistribution(fst d))
let d2 = computeD S0 B r sigma T
let price2 = price - (S0 * normal.CumulativeDistribution(snd d2)) + (strike * exp(-r * T) * normal.CumulativeDistribution(fst d2))
let d3 = computeD (1.0/S0) (1.0/B) r sigma T
let price3 = price2 + (aux1 * S0 * normal.CumulativeDistribution(snd d3)) - (aux2 * strike * exp(-r * T) * normal.CumulativeDistribution(fst d3))
let d4 = computeD (1.0 / (strike * S0)) (1.0 / (B * B)) r sigma T
let price4 = price3 - (aux1 * S0 * normal.CumulativeDistribution(snd d4)) + (aux2 * strike * exp(-r * T) * normal.CumulativeDistribution(fst d4))
price4
// Monte Carlo pricing
let priceAsianArithmeticMeanMC S0 strike r T sigma numTrajectories numSamples =
let normal = new Normal(0.0, 1.0)
let deltaT = T / float numSamples
let payoffs = seq { for n in 1 .. numTrajectories do
let assetPath = getAssetPath S0 r deltaT sigma normal numSamples |> Seq.toList
yield assetPath |> asianArithmeticMeanCallPayoff strike
}
let discountFactor = exp(-r * T)
let priceMC = discountFactor * payoffs.Mean()
let stddevMC = discountFactor * payoffs.StandardDeviation() / sqrt(float numTrajectories)
(priceMC, stddevMC)
let priceDoubleBarrierCallMC S0 strike upperBarrier lowerBarrier r T sigma numTrajectories numSamples =
let normal = new Normal(0.0, 1.0)
let deltaT = T / float numSamples
let payoffs = seq { for n in 1 .. numTrajectories do
let assetPath = getAssetPath S0 r deltaT sigma normal numSamples |> Seq.toList
let price = assetPath |> doubleBarrierEuropeanCallPayoff strike upperBarrier lowerBarrier
yield price
}
let discountFactor = exp(-r * T)
let priceMC = discountFactor * payoffs.Mean()
let stddevMC = discountFactor * payoffs.StandardDeviation() / sqrt(float numTrajectories)
(priceMC, stddevMC)
let controlVariate numTrajectories discountFactor (samplesFactory : unit -> seq<float>) controlExact (payoff : seq<float> -> float) (controlPayoff) =
let samplePayoffs = seq { for n in 1 .. numTrajectories do
let assetPath = samplesFactory() |> Seq.toList
let payoff = (payoff assetPath, controlPayoff assetPath)
yield payoff
}
|> Seq.toList
let payoffs = samplePayoffs |> Seq.map (fun i -> fst i) |> Seq.toArray
let payoffsControlVariate = samplePayoffs |> Seq.map (fun i -> snd i) |> Seq.toArray
let variance = ArrayStatistics.PopulationCovariance(payoffs, payoffs)
let varianceControlVariate = ArrayStatistics.PopulationCovariance(payoffsControlVariate, payoffsControlVariate)
let covariance = ArrayStatistics.PopulationCovariance(payoffs, payoffsControlVariate)
let correlation = covariance / sqrt(varianceControlVariate * variance)
let priceMC = discountFactor * payoffs.Mean()
let stddevMC = discountFactor * payoffs.StandardDeviation() / sqrt(float numTrajectories)
let priceControlVariateMC = discountFactor * payoffsControlVariate.Mean()
let stddevControlVariateMC = discountFactor * payoffsControlVariate.StandardDeviation() / sqrt(float numTrajectories)
let price = priceMC - ((covariance / varianceControlVariate) * (priceControlVariateMC - controlExact))
let stddev = stddevMC * sqrt(1.0 - (correlation * correlation))
(price, stddev)
let priceDoubleBarrierCall_MC_CV_EU S0 strike upperBarrier lowerBarrier r T sigma numTrajectories numSamples =
let normal = new Normal(0.0, 1.0)
let discountFactor = exp(-r * T)
let deltaT = T / float numSamples;
let priceEuropeanExact = priceEuropeanCall S0 strike r T sigma
let europeanCallfn = (fun assetPath -> europeanCallPayoff strike assetPath)
let barrierCallPayoff = (fun assetPath -> doubleBarrierEuropeanCallPayoff strike upperBarrier lowerBarrier assetPath)
controlVariate numTrajectories
discountFactor
(fun () -> getAssetPath S0 r deltaT sigma normal numSamples)
priceEuropeanExact
barrierCallPayoff
europeanCallfn
let priceDoubleBarrierCall_MC_CV_UpAndOut S0 strike upperBarrier lowerBarrier r T sigma numTrajectories numSamples =
let normal = new Normal(0.0, 1.0)
let discountFactor = exp(-r * T)
let deltaT = T / float numSamples;
let priceUpAndOutExact = priceBarrierUpAndOutCall S0 strike upperBarrier r T sigma
let barrierCallPayoff = (fun assetPath -> doubleBarrierEuropeanCallPayoff strike upperBarrier lowerBarrier assetPath)
let upAndOut = (fun assetPath -> upAndOutEuropeanCallPayoff strike upperBarrier assetPath)
controlVariate numTrajectories
discountFactor
(fun () -> getAssetPath S0 r deltaT sigma normal numSamples)
priceUpAndOutExact
barrierCallPayoff
upAndOut
let asian =
let S0 = 100.0
let strike = 90.0
let r = 0.05
let T = 1.0
let sigma = 0.2
let numTrajectories = 100000
let numSamples = 12
let result = priceAsianArithmeticMeanMC S0 strike r T sigma numTrajectories numSamples
printfn "Asian arithmetic mean:%f stddev:%f" (fst result) (snd result)
let barrier =
let S0 = 100.0
let strike = 90.0
let upperBarrier = 160.0
let lowerBarrier = 75.0
let r = 0.05
let T = 0.5
let sigma = 0.4
let numTrajectories = 1000
let numSamples = 10000
let result = priceDoubleBarrierCallMC S0 strike upperBarrier lowerBarrier r T sigma numTrajectories numSamples
printfn "Monte Carlo double barrier call:%f stddev:%f" (fst result) (snd result)
let resultMCCVEU = priceDoubleBarrierCall_MC_CV_EU S0 strike upperBarrier lowerBarrier r T sigma numTrajectories numSamples
printfn "Control Variate with European Call price:%f stddev:%f" (fst resultMCCVEU) (snd resultMCCVEU)
let resultMCCVUpAndOut = priceDoubleBarrierCall_MC_CV_UpAndOut S0 strike upperBarrier lowerBarrier r T sigma numTrajectories numSamples
printfn "Control Variate with Up And Out price:%f stddev:%f" (fst resultMCCVUpAndOut) (snd resultMCCVUpAndOut)
[<EntryPoint>]
let main argv =
asian
barrier
0
@taumuon

This comment has been minimized.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.