of random numbers will lead to (slightly) different estimates. 20,000 paths per simulation
can also be considered a bit too low for robust Monte Carlo estimates (leading, however,
to high valuation speeds). By contrast, the binomial option pricing model with 1,000 time
steps is rather precise but also takes much longer in this case.
Again, we can try NumPy vectorization techniques to come up with equally precise but
faster results from the binomial approach. The binomial_np function might seem a bit
cryptic at first sight; however, when you step through the individual construction steps and
inspect the results, it becomes clear what happens behind the (NumPy) scenes:
In [ 68 ]: def binomial_np(strike):
”’ Binomial option pricing with NumPy.
Parameters
==========
strike : float
strike price of the European call option
”’
# Index Levels with NumPy
mu = np.arange(M + 1 )
mu = np.resize(mu, (M + 1 , M + 1 ))
md = np.transpose(mu)
mu = u ** (mu - md)
md = d ** md
S = S0 * mu * md
# Valuation Loop
pv = np.maximum(S - strike, 0 )
z = 0
for t in range(M - 1 , - 1 , - 1 ): # backward iteration
pv[ 0 :M - z, t] = (q * pv[ 0 :M - z, t + 1 ]
+ ( 1 – q) * pv[ 1 :M - z + 1 , t + 1 ]) * df
z += 1
return pv[ 0 , 0 ]
Let us briefly take a look behind the scenes. For simplicity and readability, consider only
M=4 time steps. The first step:
In [ 69 ]: M = 4 # four time steps only
mu = np.arange(M + 1 )
mu
Out[69]: array([0, 1, 2, 3, 4])
The second step of the construction:
In [ 70 ]: mu = np.resize(mu, (M + 1 , M + 1 ))
mu
Out[70]: array([[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4]])
The third one:
In [ 71 ]: md = np.transpose(mu)
md
Out[71]: array([[0, 0, 0, 0, 0],
[1, 1, 1, 1, 1],
[2, 2, 2, 2, 2],
[3, 3, 3, 3, 3],
[4, 4, 4, 4, 4]])
The fourth and fifth steps:
In [ 72 ]: mu = u ** (mu - md)
mu.round( 3 )
Out[72]: array([[ 1. , 1.006, 1.013, 1.019, 1.026],
[ 0.994, 1. , 1.006, 1.013, 1.019],
