This is the second post in a series about Python for Finance. We will pick up the example of the first post and show how to incorporate variance reduction techniques for better accuracy and increased speed.
In the first post, we presented a rather concise and nevertheless quite fast implementation of the Longstaff-Schwartz algorithm (LSM) to value American options by Monte Carlo simulation (MCS). Apart from pure Python, the script presented uses extensively the array library NumPy. The enhanced Python example of this post will again make use of the capabilities of this important library.
One source of error when implementing MCS algorithms stems from the (pseudo-)random numbers generated. If you need a set of standard normally distributed random numbers you would expect to get a set of numbers with a mean of zero and a standard deviation of one. However, when generating only a relatively small set of (pseudo-)random numbers, you might experience the following:
In : rn=standard_normal(1000) In : mean(rn) Out: -0.019913964971170758 In : std(rn) Out: 1.0221048379067881 In :
As you can see, the first and second moments deviate from the desired values. However, with two lines of code you can easily correct the first and second moment:
In : rn=rn-mean(rn) In : rn=rn/std(rn) In : mean(rn) Out: -2.7644553313166398e-17 In : std(rn) Out: 0.99999999999999922 In :
In the script presented below, there are similar lines of code which correct the moments of the set of pseudo-random numbers. This variance reduction technique is called Moment Matching. In addition, a technique called Antithetic Variates is implemented where for each pseudo-random number the negative value is also used. For these and other variance reduction techniques refer to Glasserman, Paul (2004): Monte Carlo Methods in Financial Engineering. Springer, New York.
The script values the first (non-trivial) example of the seminal paper Longstaff-Schwartz (2001). For details refer to the paper itself (see publisher’s Web site) or to the recent book “Derivatives Analytics with Python” by Dr. Yves J. Hilpisch, MD of Visixion GmbH.
# # Least-Squares Monte Carlo for American Put # with Antithetic Paths and Moment Matching # (c) Visixion GmbH -- For Illustration Purposes Only. # from pylab import * from time import time t0=time() ## Parameters M=50;I=50000 # Simulation S0=36.;K=40.;r=0.06;vol=0.2;T=1.0 # Option dt=T/M;df=exp(-r*dt) # Time & Discounting ## Random Numbers ran=standard_normal((M,I/2)) rand=concatenate((ran,-ran),1) # Antithetic Variates rand=rand-mean(rand) # Matching of 1st Moment rand=rand/std(rand) # Matching of 2nd Moment ## LSM Algorithm S=ones((M+1,I),'d') S[1:,:]=exp(cumsum((r-vol**2/2)*dt+ vol*rand*sqrt(dt), axis=0)) V=maximum(K-S*S0,0) for t in range(M-1,0,-1): # Dynamic Optimization/Exercise Decisions V[t,:]=where(V[t,:]>polyval(polyfit(S[t,:], V[t+1,:]*df,3),S[t,:]), V[t,:],V[t+1,:]*df) ## Output print "Value %8.3f" %(sum(V[1,:]*df)/I) # LSM estimator print "Time %8.3f" %(time()-t0)
This script shows a better convergence of the LSM estimator and may therefore reach similar accuracy levels as the first version in much less time (or better accuracy in the same time). If you want to benchmark the scripts: the true value of the American put option is 4.486 when using a binomial model with 500 time steps.
Go to Python for Finance (I)
Go to Python for Finance (III)