Edited to include VBA code for comparison
Also, we know the analytical value of the simple Call option, which is 8.021, towards which the Monte-Carlo should converge, which makes the comparison easier.
Excel VBA gives 8.067 based on averaging 5 Monte-Carlo simulations (7.989, 8.187, 8.045, 8.034, 8.075)
Python gives 7.973 based on 5 MCs (7.913, 7.915, 8.203, 7.739, 8.095) and a larger Variance!
The VBA code is using a rather bad way to produce samples from Standard Normal...
I am running a super simple code in Python to price European Call Option via Monte Carlo, and I am surprised at how "bad" the convergence is with 10,000 "simulated paths". Usually, when running a Monte-Carlo for this simple problem in C++ or even VBA, I get better convergence.
I show the code below (the code is taken from Textbook "Python for Finance" and I run in in Visual Studio Code under Python 3.7.7, 64-bit version): I get the following results, as an example: Run 1 = 7.913, Run 2 = 7.915, Run 3 = 8.203, Run 4 = 7.739, Run 5 = 8.095,
Results such as the above, that differ by so much, would be unacceptable. How can the convergence be improved??? (Obviously by running more paths, but as I said: for 10,000 paths, the result should already have converged much better):
#MonteCarlo valuation of European Call Option import math import numpy as np #Parameter Values S_0 = 100. # initial value K = 105. # strike T = 1.0 # time to maturity r = 0.05 # short rate (constant) sigma = 0.2 # vol nr_simulations = 10000 #Valuation Algo: # Notice the vectorization below, instead of a loop z = np.random.standard_normal(nr_simulations) # Notice that the S_T below is a VECTOR! S_T = S_0 * np.exp((r-0.5*sigma**2)*T+math.sqrt(T)*sigma*z) #Call option pay-off at maturity (Vector!) C_T = np.maximum((S_T-K),0) # C_0 is a scalar C_0 = math.exp(-r*T)*np.average(C_T) print('Value of the European Call is: ', C_0)
I also include VBA code, which produces slightly better results (in my opinion): with the VBA code below, I get 7.989, 8.187, 8.045, 8.034, 8.075.
Option Explicit Sub monteCarlo() ' variable declaration ' stock initial & final values, option pay-off at maturity Dim stockInitial, stockFinal, optionFinal As Double ' r = rate, sigma = volatility, strike = strike price Dim r, sigma, strike As Double 'maturity of the option Dim maturity As Double ' instatiate variables stockInitial = 100# r = 0.05 maturity = 1# sigma = 0.2 strike = 105# ' normal is Standard Normal Dim normal As Double ' randomNr is randomly generated nr via "rnd()" function, between 0 & 1 Dim randomNr As Double ' variable for storing the final result value Dim result As Double Dim i, j As Long, monteCarlo As Long monteCarlo = 10000 For j = 1 To 5 result = 0# For i = 1 To monteCarlo ' get random nr between 0 and 1 randomNr = Rnd() 'max(Rnd(), 0.000000001) ' standard Normal normal = Application.WorksheetFunction.Norm_S_Inv(randomNr) stockFinal = stockInitial * Exp((r - (0.5 * (sigma ^ 2)))*maturity + (sigma * Sqr(maturity) * normal)) optionFinal = max((stockFinal - strike), 0) result = result + optionFinal Next i result = result / monteCarlo result = result * Exp(-r * maturity) Worksheets("sheet1").Cells(j, 1) = result Next j MsgBox "Done" End Sub Function max(ByVal number1 As Double, ByVal number2 As Double) If number1 > number2 Then max = number1 Else max = number2 End If End Function