Monte Carlo Simulation

Please read my Disclaimer.

Function European_Call_bMCBS( _
K As Double, _
T As Double, _
S As Double, _
sig As Double, _
r As Double, _
div As Double, _
N As Long, _
M As Long) As Double
'European Call Option by Monte Carlo with Black-Scholes
'Clewlow/Strickland: Implementing Derivatives Models, page 85
Dim dt As Double
Dim nudt As Double, sigsdt As Double
Dim lnS As Double, lnSt As Double
Dim sum_CT As Double, sum_CT2 As Double
Dim e As Double
Dim St As Double, CT As Double
Dim SD As Double, SE As Double
Dim i As Long, j As Long
'Precompute constants
dt = T / N
nudt = (r - div - 0.5 * sig ^ 2#) * dt
sigsdt = sig * Sqr(dt)
lnS = Log(S)
sum_CT = 0#
sum_CT2 = 0#
For j = 1 To M
'For each simulation
lnSt = lnS
For i = 1 To N
'For each time step
' e = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
e = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
lnSt = lnSt + nudt + sigsdt * e 'Evolve the stock price
Next i
St = Exp(lnSt)
CT = St - K: If CT < 0# Then CT = 0#
sum_CT = sum_CT + CT
sum_CT2 = sum_CT + CT * CT
Next j
SD = Sqr(sum_CT * sum_CT / M - sum_CT2) * Exp(-2# * r * T) / (M - 1) 'Standard deviation
SE = SD / Sqr(M) 'Standard error
European_Call_bMCBS = sum_CT / M * Exp(-r * T)
End Function
Function European_Call_bMCBS_AVR( _
K As Double, _
T As Double, _
S As Double, _
sig As Double, _
r As Double, _
div As Double, _
N As Long, _
M As Long) As Double
'European Call Option by Monte Carlo with Black-Scholes and Antithetic Variance Reduction
'Clewlow/Strickland: Implementing Derivatives Models, page 89
Dim dt As Double
Dim nudt As Double, sigsdt As Double
Dim lnS As Double, LnSt1 As Double, LnSt2 As Double
Dim sum_CT As Double, sum_CT2 As Double
Dim e As Double
Dim St1 As Double, St2 As Double
Dim CT As Double, CT1 As Double, CT2 As Double
Dim SD As Double, SE As Double
Dim i As Long, j As Long
'Precompute constants
dt = T / N
nudt = (r - div - 0.5 * sig ^ 2#) * dt
sigsdt = sig * Sqr(dt)
lnS = Log(S)
sum_CT = 0#
sum_CT2 = 0#
For j = 1 To M
'For each simulation
LnSt1 = lnS
LnSt2 = lnS
For i = 1 To N
'For each time step
' e = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
e = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
LnSt1 = LnSt1 + nudt + sigsdt * e 'Evolve the stock price
LnSt2 = LnSt2 + nudt + sigsdt * -e 'Evolve the stock price
Next i
St1 = Exp(LnSt1)
St2 = Exp(LnSt2)
CT1 = St1 - K: If CT1 < 0# Then CT1 = 0#
CT2 = St2 - K: If CT2 < 0# Then CT2 = 0#
CT = (CT1 + CT2) / 2#
sum_CT = sum_CT + CT
sum_CT2 = sum_CT + CT * CT
Next j
SD = Sqr(sum_CT * sum_CT / M - sum_CT2) * Exp(-2# * r * T) / (M - 1) 'Standard deviation
SE = SD / Sqr(M) 'Standard error
European_Call_bMCBS_AVR = sum_CT / M * Exp(-r * T)
End Function
Function European_Call_bMCBS_DCV( _
K As Double, _
T As Double, _
S As Double, _
sig As Double, _
r As Double, _
div As Double, _
N As Long, _
M As Long) As Double
'European Call Option by Monte Carlo with Black-Scholes and Delta-based Control Variate
'Clewlow/Strickland: Implementing Derivatives Models, page 97
Dim dt As Double
Dim nudt As Double, sigsdt As Double, erddt As Double
Dim sum_CT As Double, sum_CT2 As Double
Dim e As Double
Dim St As Double, Stn As Double
Dim CT As Double
Dim SD As Double, SE As Double
Dim beta1 As Double, t1 As Double
Dim cv As Double
Dim delta1 As Double
Dim i As Long, j As Long
'Dim eps As Variant
'eps = Array(0.5391, 0.1692, 0.4583, 0.4164, 0.9223, 0.1597, 0.4011, -0.2382, 1.2152, -0.5506)
'Precompute constants
dt = T / N
nudt = (r - div - 0.5 * sig ^ 2#) * dt
sigsdt = sig * Sqr(dt)
erddt = Exp((r - div) * dt)
beta1 = -1#
'Debug.Print "dt = " & dt & ", nudt = " & nudt & ", sigsdt = " & sigsdt & ", erddt = " & erddt & ", beta1 = " & beta1
sum_CT = 0#
sum_CT2 = 0#
For j = 1 To M
'For each simulation
St = S
cv = 0#
For i = 1 To N
'For each time step
t1 = CDbl(i - 1) * dt
delta1 = Black_Scholes_Delta(St, t1, K, T, sig, r, div)
' Debug.Print "e = " & e & ", St = " & St & ", delta1 = " & delta1 & ", cv = " & cv
' e = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
e = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
' e = eps(i - 1)
Stn = St * Exp(nudt + sigsdt * e)
cv = cv + delta1 * (Stn - St * erddt) '* exp(T - (t1 + dt)) 'inflation factor *exp can be omitted
St = Stn
Next i
' Debug.Print "e = " & e & ", St = " & St & ", delta1 = " & delta1 & ", cv = " & cv
CT = St - K: If CT < 0# Then CT = 0#
CT = CT + beta1 * cv
sum_CT = sum_CT + CT
sum_CT2 = sum_CT + CT * CT
' Debug.Print "CT = " & CT & ", CT2 = " & CT * CT
Next j
SD = Sqr(sum_CT * sum_CT / M - sum_CT2) * Exp(-2# * r * T) / (M - 1) 'Standard deviation
SE = SD / Sqr(M) 'Standard error
European_Call_bMCBS_DCV = sum_CT / M * Exp(-r * T)
End Function
Function European_Call_bMCBS_ADCV( _
K As Double, _
T As Double, _
S As Double, _
sig As Double, _
r As Double, _
div As Double, _
N As Long, _
M As Long) As Double
'European Call Option by Monte Carlo with Black-Scholes and Antithetic and Delta-based Control Variates
'Clewlow/Strickland: Implementing Derivatives Models, page 100
Dim dt As Double
Dim nudt As Double, sigsdt As Double, erddt As Double
Dim sum_CT As Double, sum_CT2 As Double
Dim e As Double
Dim St1 As Double, Stn1 As Double
Dim St2 As Double, Stn2 As Double
Dim CT As Double, CT1 As Double, CT2 As Double
Dim SD As Double, SE As Double
Dim beta1 As Double, t1 As Double
Dim cv1 As Double, cv2 As Double
Dim delta1 As Double, delta2 As Double
Dim i As Long, j As Long
'Dim eps As Variant
'eps = Array(0.5391, 0.1692, 0.4583, 0.4164, 0.9223, 0.1597, 0.4011, -0.2382, 1.2152, -0.5506)
'Precompute constants
dt = T / N
nudt = (r - div - 0.5 * sig ^ 2#) * dt
sigsdt = sig * Sqr(dt)
erddt = Exp((r - div) * dt)
beta1 = -1#
'Debug.Print "dt = " & dt & ", nudt = " & nudt & ", sigsdt = " & sigsdt & ", erddt = " & erddt & ", beta1 = " & beta1
sum_CT = 0#
sum_CT2 = 0#
For j = 1 To M
'For each simulation
St1 = S
St2 = S
cv1 = 0#
cv2 = 0#
For i = 1 To N
'For each time step
t1 = CDbl(i - 1) * dt
delta1 = Black_Scholes_Delta(St1, t1, K, T, sig, r, div)
delta2 = Black_Scholes_Delta(St2, t1, K, T, sig, r, div)
' Debug.Print "e = " & e & ", St = " & St & ", delta1 = " & delta1 & ", cv = " & cv
' e = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
e = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
' e = eps(i - 1)
Stn1 = St1 * Exp(nudt + sigsdt * e)
Stn2 = St2 * Exp(nudt + sigsdt * -e)
cv1 = cv1 + delta1 * (Stn1 - St1 * erddt) '* exp(T - (t1 + dt)) 'inflation factor *exp can be omitted
cv2 = cv2 + delta2 * (Stn2 - St2 * erddt) '* exp(T - (t1 + dt)) 'inflation factor *exp can be omitted
St1 = Stn1
St2 = Stn2
Next i
' Debug.Print "e = " & e & ", St = " & St & ", delta1 = " & delta1 & ", cv = " & cv
CT1 = St1 - K: If CT1 < 0# Then CT1 = 0#
CT2 = St2 - K: If CT2 < 0# Then CT2 = 0#
CT = (CT1 + beta1 * cv1 + CT2 + beta1 * cv2) / 2#
sum_CT = sum_CT + CT
sum_CT2 = sum_CT + CT * CT
' Debug.Print "CT = " & CT & ", CT2 = " & CT * CT
Next j
SD = Sqr(sum_CT * sum_CT / M - sum_CT2) * Exp(-2# * r * T) / (M - 1) 'Standard deviation
SE = SD / Sqr(M) 'Standard error
European_Call_bMCBS_ADCV = sum_CT / M * Exp(-r * T)
End Function
Function European_Call_bMCBS_ADGCV( _
K As Double, _
T As Double, _
S As Double, _
sig As Double, _
r As Double, _
div As Double, _
N As Long, _
M As Long) As Double
'European Call Option by Monte Carlo with Black-Scholes and Antithetic, Delta-
'and Gamma-based Control Variates
'Clewlow/Strickland: Implementing Derivatives Models, page 102
Dim dt As Double
Dim nudt As Double, sigsdt As Double
Dim erddt As Double, egamma As Double
Dim sum_CT As Double, sum_CT2 As Double
Dim e As Double
Dim St1 As Double, Stn1 As Double
Dim St2 As Double, Stn2 As Double
Dim CT As Double, CT1 As Double, CT2 As Double
Dim SD As Double, SE As Double
Dim beta1 As Double, beta2 As Double, t1 As Double
Dim cv1 As Double, cv2 As Double
Dim delta1 As Double, delta2 As Double
Dim gamma1 As Double, gamma2 As Double
Dim i As Long, j As Long
'Dim eps As Variant
'eps = Array(0.5391, 0.1692, 0.4583, 0.4164, 0.9223, 0.1597, 0.4011, -0.2382, 1.2152, -0.5506)
'Precompute constants
dt = T / N
nudt = (r - div - 0.5 * sig ^ 2#) * dt
sigsdt = sig * Sqr(dt)
erddt = Exp((r - div) * dt)
egamma = Exp((2# * (r - div) + sig * sig) * dt) - 2# * erddt + 1#
beta1 = -1#
beta2 = -0.5
'Debug.Print "dt = " & dt & ", nudt = " & nudt & ", sigsdt = " & sigsdt & ", erddt = " & erddt & ", beta1 = " & beta1
sum_CT = 0#
sum_CT2 = 0#
For j = 1 To M
'For each simulation
St1 = S
St2 = S
cv1 = 0#
cv2 = 0#
For i = 1 To N
'For each time step
'Compute hedge sensitivities
t1 = CDbl(i - 1) * dt
delta1 = Black_Scholes_Delta(St1, t1, K, T, sig, r, div)
delta2 = Black_Scholes_Delta(St2, t1, K, T, sig, r, div)
gamma1 = Black_Scholes_Gamma(St1, t1, K, T, sig, r, div)
gamma2 = Black_Scholes_Gamma(St2, t1, K, T, sig, r, div)
'Evolve asset prices
' Debug.Print "e = " & e & ", St = " & St & ", delta1 = " & delta1 & ", cv = " & cv
' e = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
e = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
' e = eps(i - 1)
Stn1 = St1 * Exp(nudt + sigsdt * e)
Stn2 = St2 * Exp(nudt + sigsdt * -e)
'Accumulate control variates
cv1 = cv1 + delta1 * (Stn1 - St1 * erddt) + _
delta2 * (Stn2 - St2 * erddt)
cv2 = cv2 + gamma1 * ((Stn1 - St1) ^ 2 - St1 ^ 2 * egamma) + _
gamma2 * ((Stn2 - St2) ^ 2 - St2 ^ 2 * egamma)
St1 = Stn1
St2 = Stn2
Next i
' Debug.Print "e = " & e & ", St = " & St & ", delta1 = " & delta1 & ", cv = " & cv
CT1 = St1 - K: If CT1 < 0# Then CT1 = 0#
CT2 = St2 - K: If CT2 < 0# Then CT2 = 0#
CT = (CT1 + CT2 + beta1 * cv1 + beta2 * cv2) / 2#
sum_CT = sum_CT + CT
sum_CT2 = sum_CT + CT * CT
' Debug.Print "CT = " & CT & ", CT2 = " & CT * CT
Next j
SD = Sqr(sum_CT * sum_CT / M - sum_CT2) * Exp(-2# * r * T) / (M - 1) 'Standard deviation
SE = SD / Sqr(M) 'Standard error
European_Call_bMCBS_ADGCV = sum_CT / M * Exp(-r * T)
End Function
Function European_Spread_bMCBS( _
K As Double, _
T As Double, _
S1 As Double, _
S2 As Double, _
sig1 As Double, _
sig2 As Double, _
div1 As Double, _
div2 As Double, _
rho As Double, _
r As Double, _
N As Long, _
M As Long) As Double
'European Spread Option by Monte Carlo with Black-Scholes
'Clewlow/Strickland: Implementing Derivatives Models, page 110
Dim dt As Double
Dim nu1dt As Double, nu2dt As Double
Dim sig1sdt As Double, sig2sdt As Double
Dim sum_CT As Double, sum_CT2 As Double
Dim e1 As Double, e2 As Double
Dim z1 As Double, z2 As Double
Dim St1 As Double, St2 As Double
Dim CT As Double
Dim SD As Double, SE As Double
Dim srho As Double
Dim i As Long, j As Long
'Precompute constants
dt = T / N
nu1dt = (r - div1 - 0.5 * sig1 ^ 2#) * dt
nu2dt = (r - div2 - 0.5 * sig2 ^ 2#) * dt
sig1sdt = sig1 * Sqr(dt)
sig2sdt = sig2 * Sqr(dt)
srho = Sqr(1# - rho ^ 2#)
sum_CT = 0#
sum_CT2 = 0#
For j = 1 To M
'For each simulation
St1 = S1
St2 = S2
For i = 1 To N
'For each time step
e1 = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
e2 = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
' e1 = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
' Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
' e2 = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
' Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
z1 = e1
z2 = rho * e1 + srho * e2
St1 = St1 * Exp(nu1dt + sig1sdt * z1)
St2 = St2 * Exp(nu2dt + sig2sdt * z2)
Next i
CT = St1 - St2 - K: If CT < 0# Then CT = 0#
sum_CT = sum_CT + CT
sum_CT2 = sum_CT + CT * CT
Next j
SD = Sqr(sum_CT * sum_CT / M - sum_CT2) * Exp(-2# * r * T) / (M - 1) 'Standard deviation
SE = SD / Sqr(M) 'Standard error
European_Spread_bMCBS = sum_CT / M * Exp(-r * T)
End Function
Function European_Down_and_Out_Call_bMCBS( _
K As Double, _
T As Double, _
S As Double, _
sig As Double, _
r As Double, _
div As Double, _
H As Long, _
N As Long, _
M As Long) As Double
'European Down and Out Call Option by Monte Carlo with Black-Scholes
'Clewlow/Strickland: Implementing Derivatives Models, page 116
Dim dt As Double
Dim nudt As Double, sigsdt As Double
Dim lnS As Double, lnSt As Double
Dim sum_CT As Double, sum_CT2 As Double
Dim e As Double
Dim St As Double, CT As Double
Dim SD As Double, SE As Double
Dim i As Long, j As Long
Dim bBarrierCrossed As Boolean
#If MODULE_TEST Then
e = Rnd(-1)
Randomize 1 'maybe even additional parameter to function
#End If
'Precompute constants
dt = T / N
nudt = (r - div - sig * sig / 2#) * dt
sigsdt = sig * Sqr(dt)
sum_CT = 0#
sum_CT2 = 0#
For j = 1 To M
'For each simulation
St = S
bBarrierCrossed = False
For i = 1 To N
'For each time step
' e = Application.WorksheetFunction.NormSInv(Rnd()) 'Standard normal sample
e = Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + Rnd() + _
Rnd() + Rnd() + Rnd() + Rnd() - 6# 'Standard normal sample
St = St * Exp(nudt + sigsdt * e) 'Evolve the stock price
If St <= H Then
bBarrierCrossed = True
Exit For
End If
Next i
If bBarrierCrossed Then
CT = 0#
Else
CT = St - K
If CT < 0# Then CT = 0#
End If
sum_CT = sum_CT + CT
sum_CT2 = sum_CT + CT * CT
Next j
SD = Sqr(Abs(sum_CT * sum_CT / M - sum_CT2)) * Exp(-2# * r * T) / (M - 1) 'Standard deviation
SE = SD / Sqr(M) 'Standard error
European_Down_and_Out_Call_bMCBS = sum_CT / M * Exp(-r * T)
End Function