In [1]:
Copied!
import matplotlib.pyplot as plt
import numpy as np
import qmcpy as qp
import matplotlib.pyplot as plt
import numpy as np
import qmcpy as qp
In [2]:
Copied!
seed = 7
seed = 7
Comparison of multilevel (Quasi-)Monte Carlo for an Asian option problem¶
Compute the exact value of the Asian option with single level QMC, for an increasing number of time steps:
In [4]:
Copied!
for level in range(5):
aco = qp.FinancialOption(qp.Sobol(2*2**level, seed=seed), option="ASIAN", volatility=.2, start_price=100, strike_price=100, interest_rate=.05)
approx_solution, data = qp.CubQMCSobolG(aco, abs_tol=1e-4).integrate()
print("Asian Option true value (%d time steps): %.5f (to within 1e-4)"%(2*2**level, approx_solution))
for level in range(5):
aco = qp.FinancialOption(qp.Sobol(2*2**level, seed=seed), option="ASIAN", volatility=.2, start_price=100, strike_price=100, interest_rate=.05)
approx_solution, data = qp.CubQMCSobolG(aco, abs_tol=1e-4).integrate()
print("Asian Option true value (%d time steps): %.5f (to within 1e-4)"%(2*2**level, approx_solution))
Asian Option true value (2 time steps): 5.63593 (to within 1e-4) Asian Option true value (4 time steps): 5.73171 (to within 1e-4) Asian Option true value (8 time steps): 5.75526 (to within 1e-4) Asian Option true value (16 time steps): 5.76112 (to within 1e-4) Asian Option true value (32 time steps): 5.76260 (to within 1e-4)
This function compares 4 different algorithms: Multilevel Monte Carlo (CubMLMC), Multilevel Quasi-Monte Carlo (CubMLQMC), continuation Multilevel Monte Carlo (CubMLMCCont) and Multilevel Quasi-Monte Carlo (CubMLQMCCont):
In [5]:
Copied!
def eval_option(option_mc, option_qmc, abs_tol):
stopping_criteria = {
"MLMC" : qp.CubMLMC(option_mc, abs_tol=abs_tol, levels_max=15),
"continuation MLMC" : qp.CubMLMCCont(option_mc, abs_tol=abs_tol, levels_max=15),
"MLQMC" : qp.CubMLQMC(option_qmc, abs_tol=abs_tol, levels_max=15),
"continuation MLQMC" : qp.CubMLQMCCont(option_qmc, abs_tol=abs_tol, levels_max=15)
}
levels = []
times = []
for name, stopper in stopping_criteria.items():
sol, data = stopper.integrate()
levels.append(data.levels)
times.append(data.time_integrate)
print("\t%-20s solution %-10.4f number of levels %-6d time %.3f"%(name, sol, levels[-1], times[-1]))
return levels, times
def eval_option(option_mc, option_qmc, abs_tol):
stopping_criteria = {
"MLMC" : qp.CubMLMC(option_mc, abs_tol=abs_tol, levels_max=15),
"continuation MLMC" : qp.CubMLMCCont(option_mc, abs_tol=abs_tol, levels_max=15),
"MLQMC" : qp.CubMLQMC(option_qmc, abs_tol=abs_tol, levels_max=15),
"continuation MLQMC" : qp.CubMLQMCCont(option_qmc, abs_tol=abs_tol, levels_max=15)
}
levels = []
times = []
for name, stopper in stopping_criteria.items():
sol, data = stopper.integrate()
levels.append(data.levels)
times.append(data.time_integrate)
print("\t%-20s solution %-10.4f number of levels %-6d time %.3f"%(name, sol, levels[-1], times[-1]))
return levels, times
Define the Multilevel Asian options:
In [8]:
Copied!
option_mc = qp.FinancialOption(qp.IIDStdUniform(seed=seed), option="asian")
option_qmc = qp.FinancialOption(qp.Lattice(seed=seed, replications=8), option="asian")
option_mc = qp.FinancialOption(qp.IIDStdUniform(seed=seed), option="asian")
option_qmc = qp.FinancialOption(qp.Lattice(seed=seed, replications=8), option="asian")
Run and compare each of the 4 algorithms for the Asian option problem:
In [9]:
Copied!
eval_option(option_mc, option_qmc, abs_tol=5e-3);
eval_option(option_mc, option_qmc, abs_tol=5e-3);
MLMC solution 1.7859 number of levels 5 time 14.886 continuation MLMC solution 1.7825 number of levels 3 time 9.339 MLQMC solution 1.7818 number of levels 3 time 0.492 continuation MLQMC solution 1.7812 number of levels 3 time 0.180
Repeat this comparison for a sequence of decreasing tolerances, with 5 different random seeds each. This will allow us to visualize the asymptotic cost complexity of each method.
In [10]:
Copied!
repetitions = 5
tolerances = 5*np.logspace(-1, -3, num=5)
levels = {}
times = {}
for t in range(len(tolerances)):
for r in range(repetitions):
print("tolerance = %10.4e, repetition = %d/%d"%(tolerances[t], r + 1, repetitions))
levels[t, r], times[t, r] = eval_option(option_mc, option_qmc, tolerances[t])
repetitions = 5
tolerances = 5*np.logspace(-1, -3, num=5)
levels = {}
times = {}
for t in range(len(tolerances)):
for r in range(repetitions):
print("tolerance = %10.4e, repetition = %d/%d"%(tolerances[t], r + 1, repetitions))
levels[t, r], times[t, r] = eval_option(option_mc, option_qmc, tolerances[t])
tolerance = 5.0000e-01, repetition = 1/5 MLMC solution 1.7574 number of levels 3 time 0.012 continuation MLMC solution 1.7765 number of levels 3 time 0.014 MLQMC solution 1.7840 number of levels 3 time 0.012 continuation MLQMC solution 1.7897 number of levels 3 time 0.012 tolerance = 5.0000e-01, repetition = 2/5 MLMC solution 1.8907 number of levels 3 time 0.005 continuation MLMC solution 1.5579 number of levels 3 time 0.010 MLQMC solution 1.7822 number of levels 3 time 0.009 continuation MLQMC solution 1.7790 number of levels 3 time 0.011 tolerance = 5.0000e-01, repetition = 3/5 MLMC solution 1.5368 number of levels 3 time 0.005 continuation MLMC solution 1.5291 number of levels 3 time 0.009 MLQMC solution 1.8200 number of levels 3 time 0.009 continuation MLQMC solution 1.8018 number of levels 3 time 0.010 tolerance = 5.0000e-01, repetition = 4/5 MLMC solution 1.6869 number of levels 3 time 0.005 continuation MLMC solution 1.8590 number of levels 3 time 0.009 MLQMC solution 1.8063 number of levels 3 time 0.008 continuation MLQMC solution 1.8205 number of levels 3 time 0.010 tolerance = 5.0000e-01, repetition = 5/5 MLMC solution 1.5989 number of levels 3 time 0.006 continuation MLMC solution 1.7292 number of levels 3 time 0.009 MLQMC solution 1.7672 number of levels 3 time 0.007 continuation MLQMC solution 1.7521 number of levels 3 time 0.008 tolerance = 1.5811e-01, repetition = 1/5 MLMC solution 1.6774 number of levels 3 time 0.010 continuation MLMC solution 1.8315 number of levels 3 time 0.016 MLQMC solution 1.7883 number of levels 3 time 0.010 continuation MLQMC solution 1.7917 number of levels 3 time 0.011 tolerance = 1.5811e-01, repetition = 2/5 MLMC solution 1.7921 number of levels 3 time 0.012 continuation MLMC solution 1.6978 number of levels 3 time 0.017 MLQMC solution 1.7939 number of levels 3 time 0.015 continuation MLQMC solution 1.7888 number of levels 3 time 0.012 tolerance = 1.5811e-01, repetition = 3/5 MLMC solution 1.7723 number of levels 3 time 0.013 continuation MLMC solution 1.8004 number of levels 4 time 0.021 MLQMC solution 1.7896 number of levels 3 time 0.011 continuation MLQMC solution 1.7780 number of levels 3 time 0.009 tolerance = 1.5811e-01, repetition = 4/5 MLMC solution 1.7671 number of levels 3 time 0.011 continuation MLMC solution 1.7734 number of levels 4 time 0.015 MLQMC solution 1.7888 number of levels 3 time 0.012 continuation MLQMC solution 1.7693 number of levels 3 time 0.008 tolerance = 1.5811e-01, repetition = 5/5 MLMC solution 1.7919 number of levels 3 time 0.012 continuation MLMC solution 1.8158 number of levels 3 time 0.010 MLQMC solution 1.7866 number of levels 3 time 0.011 continuation MLQMC solution 1.7808 number of levels 3 time 0.012 tolerance = 5.0000e-02, repetition = 1/5 MLMC solution 1.7964 number of levels 3 time 0.127 continuation MLMC solution 1.8056 number of levels 4 time 0.065 MLQMC solution 1.7835 number of levels 3 time 0.042 continuation MLQMC solution 1.7734 number of levels 3 time 0.029 tolerance = 5.0000e-02, repetition = 2/5 MLMC solution 1.7826 number of levels 3 time 0.096 continuation MLMC solution 1.7682 number of levels 3 time 0.100 MLQMC solution 1.7884 number of levels 3 time 0.043 continuation MLQMC solution 1.7705 number of levels 3 time 0.029 tolerance = 5.0000e-02, repetition = 3/5 MLMC solution 1.7604 number of levels 3 time 0.072 continuation MLMC solution 1.8240 number of levels 3 time 0.064 MLQMC solution 1.7786 number of levels 3 time 0.050 continuation MLQMC solution 1.7769 number of levels 3 time 0.025 tolerance = 5.0000e-02, repetition = 4/5 MLMC solution 1.7693 number of levels 3 time 0.094 continuation MLMC solution 1.8409 number of levels 4 time 0.178 MLQMC solution 1.7839 number of levels 3 time 0.041 continuation MLQMC solution 1.7835 number of levels 3 time 0.039 tolerance = 5.0000e-02, repetition = 5/5 MLMC solution 1.7808 number of levels 3 time 0.090 continuation MLMC solution 1.7848 number of levels 3 time 0.062 MLQMC solution 1.7768 number of levels 3 time 0.027 continuation MLQMC solution 1.7758 number of levels 3 time 0.024 tolerance = 1.5811e-02, repetition = 1/5 MLMC solution 1.7800 number of levels 3 time 1.225 continuation MLMC solution 1.7838 number of levels 3 time 0.410 MLQMC solution 1.7835 number of levels 3 time 0.099 continuation MLQMC solution 1.7781 number of levels 3 time 0.059 tolerance = 1.5811e-02, repetition = 2/5 MLMC solution 1.7941 number of levels 4 time 1.002 continuation MLMC solution 1.7848 number of levels 3 time 0.380 MLQMC solution 1.7783 number of levels 3 time 0.096 continuation MLQMC solution 1.7760 number of levels 3 time 0.056 tolerance = 1.5811e-02, repetition = 3/5 MLMC solution 1.7787 number of levels 3 time 0.809 continuation MLMC solution 1.7846 number of levels 3 time 0.791 MLQMC solution 1.7820 number of levels 3 time 0.132 continuation MLQMC solution 1.7776 number of levels 3 time 0.059 tolerance = 1.5811e-02, repetition = 4/5 MLMC solution 1.7806 number of levels 4 time 1.074 continuation MLMC solution 1.7810 number of levels 3 time 0.482 MLQMC solution 1.7806 number of levels 3 time 0.151 continuation MLQMC solution 1.7823 number of levels 3 time 0.051 tolerance = 1.5811e-02, repetition = 5/5 MLMC solution 1.7914 number of levels 4 time 0.952 continuation MLMC solution 1.7858 number of levels 3 time 0.769 MLQMC solution 1.7830 number of levels 3 time 0.111 continuation MLQMC solution 1.7804 number of levels 3 time 0.072 tolerance = 5.0000e-03, repetition = 1/5 MLMC solution 1.7827 number of levels 5 time 18.777 continuation MLMC solution 1.7823 number of levels 3 time 4.222 MLQMC solution 1.7817 number of levels 3 time 0.291 continuation MLQMC solution 1.7814 number of levels 3 time 0.425 tolerance = 5.0000e-03, repetition = 2/5 MLMC solution 1.7844 number of levels 4 time 10.280 continuation MLMC solution 1.7807 number of levels 3 time 8.755 MLQMC solution 1.7811 number of levels 3 time 0.355 continuation MLQMC solution 1.7812 number of levels 3 time 0.460 tolerance = 5.0000e-03, repetition = 3/5 MLMC solution 1.7819 number of levels 4 time 16.041 continuation MLMC solution 1.7825 number of levels 3 time 6.892 MLQMC solution 1.7813 number of levels 3 time 0.272 continuation MLQMC solution 1.7814 number of levels 3 time 0.369 tolerance = 5.0000e-03, repetition = 4/5 MLMC solution 1.7843 number of levels 4 time 12.328 continuation MLMC solution 1.7804 number of levels 3 time 4.089 MLQMC solution 1.7805 number of levels 3 time 0.498 continuation MLQMC solution 1.7823 number of levels 3 time 0.162 tolerance = 5.0000e-03, repetition = 5/5 MLMC solution 1.7822 number of levels 4 time 10.548 continuation MLMC solution 1.7805 number of levels 3 time 9.212 MLQMC solution 1.7818 number of levels 3 time 0.460 continuation MLQMC solution 1.7805 number of levels 3 time 0.508
Compute and plot the asymptotic cost complexity.
In [11]:
Copied!
avg_time = {}
for method in range(4):
avg_time[method] = [np.mean([times[t, r][method] for r in range(repetitions)]) for t in range(len(tolerances))]
avg_time = {}
for method in range(4):
avg_time[method] = [np.mean([times[t, r][method] for r in range(repetitions)]) for t in range(len(tolerances))]
In [12]:
Copied!
plt.figure(figsize=(10,7))
plt.plot(tolerances, avg_time[0], label="MLMC")
plt.plot(tolerances, avg_time[1], label="continuation MLMC")
plt.plot(tolerances, avg_time[2], label="MLQMC")
plt.plot(tolerances, avg_time[3], label="continuation MLQMC")
plt.xscale("log")
plt.yscale("log")
plt.xlabel("requested absolute tolerance")
plt.ylabel("average run time in seconds")
plt.legend();
plt.figure(figsize=(10,7))
plt.plot(tolerances, avg_time[0], label="MLMC")
plt.plot(tolerances, avg_time[1], label="continuation MLMC")
plt.plot(tolerances, avg_time[2], label="MLQMC")
plt.plot(tolerances, avg_time[3], label="continuation MLQMC")
plt.xscale("log")
plt.yscale("log")
plt.xlabel("requested absolute tolerance")
plt.ylabel("average run time in seconds")
plt.legend();
In [13]:
Copied!
max_levels = {}
for method in range(4):
levels_rep = np.array([levels[len(tolerances)-1, r][method] for r in range(repetitions)])
max_levels[method] = [np.count_nonzero(levels_rep == level)/repetitions for level in range(15)]
max_levels = {}
for method in range(4):
levels_rep = np.array([levels[len(tolerances)-1, r][method] for r in range(repetitions)])
max_levels[method] = [np.count_nonzero(levels_rep == level)/repetitions for level in range(15)]
In [14]:
Copied!
plt.figure(figsize=(14,10))
plt.subplot(2,2,1); plt.bar(range(15), max_levels[0], label="MLMC old", color="C0"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend()
plt.subplot(2,2,2); plt.bar(range(15), max_levels[1], label="MLMC new", color="C1"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend()
plt.subplot(2,2,3); plt.bar(range(15), max_levels[2], label="MLQMC old", color="C2"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend()
plt.subplot(2,2,4); plt.bar(range(15), max_levels[3], label="MLQMC new", color="C3"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend();
plt.figure(figsize=(14,10))
plt.subplot(2,2,1); plt.bar(range(15), max_levels[0], label="MLMC old", color="C0"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend()
plt.subplot(2,2,2); plt.bar(range(15), max_levels[1], label="MLMC new", color="C1"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend()
plt.subplot(2,2,3); plt.bar(range(15), max_levels[2], label="MLQMC old", color="C2"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend()
plt.subplot(2,2,4); plt.bar(range(15), max_levels[3], label="MLQMC new", color="C3"); plt.xlabel("max level"); plt.ylabel("fraction of runs"); plt.ylim(0, 1); plt.legend();
In [ ]:
Copied!