蒙特卡洛估值⼏种不同的计算⽅式（Python）
版权声明：本⽂为博主原创⽂章，未经博主允许不得转载。
蒙特卡洛模拟计算看涨期权
蒙特卡洛估值是⼀个计算量较⾼的⼀个算法，对于简单的问题也需要海量的计算，因此要⾼效的实现蒙特卡洛算法． 下⾯使⽤不同的⽅法实现蒙特卡洛算法
１．Scipy
２．纯Python
３．Numpy(1)
４．Numpy(2)
⽅法⼀：Scipy
欧式看涨期权的定价公式Black-Scholes-Merton(1973):
  公式１：     
, t时刻股票或者指数的⽔平
, ⾏权价格
, 期权的到期年限(距离到期⽇时间间隔)
, ⽆风险利率
，波动率(收益标准差)
# 倒⼊⽤到的库
from math import log, sqrt, exp
from scipy import stats
# 期权的定价计算，根据公式１．
def bsm_call_value(S_0, K, T, r, sigma):
S_0 = float(S_0)
d_1 = (log(S_0 / K) + (r + 0.5 *sigma **2) *T)/(sigma * sqrt(T))
d_2 = (log(S_0 / K) + (r - 0.5 *sigma **2) *T)/(sigma * sqrt(T))
C_0 = (S_0 * df(d_1, 0.0, 1.0) - K * exp(-r * T) * df(d_2, 0.0, 1.0))
return C_0
# 计算的⼀些初始值
S_0 = 100.0# 股票或指数初始的价格;
K = 105#  ⾏权价格
T = 1.0#  期权的到期年限(距离到期⽇时间间隔)
r = 0.05#  ⽆风险利率
sigma = 0.2# 波动率(收益标准差)
# 到期期权价值
%time print bsm_call_value(S_0, K, T, r, sigma)
估值结果
8.021********
CPU times: user 0 ns, sys: 0 ns, total: 0 ns
Wall time: 1.01 ms　
8.0213522作为蒙特卡洛估值的基准值，计算耗时1.01ms，等于0.00101s
⽅法⼆：Python
下⾯的计算仍基于BSM(balck-scholes-merton)，模型中的⾼风险标识(股票指数)在风险中⽴的情况下遵循以随机微分⽅程(SDE)表⽰的布朗运动．
随机微分⽅程：
是⼀个服从布朗运动的随机变量　
下⾯是蒙特卡洛估值公式２，中的参数与公式１中的定义相同：
变量 是服从正态分布的随机变量，是⼀个⾜够⼩的⼀个时间间隔．期权的到期年限T满⾜0<5<=T(Hilpisch的著作中)
蒙特卡洛数值化计算⽅法：
１．将到期⽇之前的时间间隔[0,T]，分割为多个等距的⼦间隔
２．进⾏ Ｉ 次的模拟， {１，２，…Ｉ}
　每次模拟，循环M个⼦间隔, {…T}
　(1)，每个间隔点，取⼀个随机数
　(2)，根据离散化公式３，计算
　(3)，计算T时刻，期权的内在价值 :
３．根据Ｉ次的模拟，计算期权到期价值：
Python实现
# 纯Python实现
from time import time
from math import exp, sqrt, log
from random import gauss, seed
seed(2000)
# 计算的⼀些初始值
S_0 = 100.0# 股票或指数初始的价格;
K = 105#  ⾏权价格
T = 1.0#  期权的到期年限(距离到期⽇时间间隔)
r = 0.05#  ⽆风险利率
sigma = 0.2# 波动率(收益标准差)
M = 50# number of time steps
dt = T/M      # time enterval
I = 20000# number of simulation
random python
start = time()
S = []    #
for i in range(I):
path = []    # 时间间隔上的模拟路径
for t in range(M+1):
if t==0:
path.append(S_0)
else:
z = gauss(0.0, 1.0)
S_t = path[t-1] * exp((r-0.5*sigma**2) * dt + sigma * sqrt(dt) * z)
path.append(S_t)
S.append(path)
# 计算期权现值
C_0 = exp(-r * T) *sum([max(path[-1] -K, 0) for path in S])/I
total_time = time() - start
print'European Option value %.6f'% C_0
print'total time is %.6f seconds'% total_time
估值结果
European Option value 8.159995
total time is 2.384639 seconds
前３０条模拟路径
# 选取部分模拟路径可视化
import matplotlib.pyplot as plt
%matplotlib inline
plt.figure(figsize=(10,7))
plt.xlabel('Time step')
plt.ylabel('index level')
for i in range(30):
plt.plot(S[i])
⽅法三：Numpy
使⽤numpy的⼀些数组计算，减少for循环使⽤
import numpy as np
from time import time
# 计算的⼀些初始值
S_0 = 100.0# 股票或指数初始的价格;
K = 105#  ⾏权价格
T = 1.0#  期权的到期年限(距离到期⽇时间间隔)
r = 0.05#  ⽆风险利率
sigma = 0.2# 波动率(收益标准差)
M = 50# number of time steps
dt = T/M      # time enterval
I = 20000# number of simulation
# 20000条模拟路径，每条路径５０个时间步数
S = np.zeros((M+1, I))
S[0] = S_0
np.random.seed(2000)
start = time()
for t in range(1, M+1):
z = np.random.standard_normal(I)
S[t] = S[t-1] * np.exp((r- 0.5 * sigma **2)* dt + sigma * np.sqrt(dt)*z)
C_0 = np.exp(-r * T)* np.sum(np.maximum(S[-1] - K, 0))/I
end = time()
# 估值结果
print'total time is %.6f seconds'%(end-start)
print'European Option Value %.6f'%C_0
估值结果
total time is 0.166926 seconds
European Option Value 7.993282
前２０条模拟路径
# 前２０条模拟路径
import matplotlib.pyplot as plt
%matplotlib inline
plt.figure(figsize=(10,7))
plt.xlabel('Time step')
plt.ylabel('index level')
for i in range(20):
plt.plot(S.T[i])
到期指数模拟⽔平
# 到期时所有模拟指数⽔平的频率直⽅图
%matplotlib inline
plt.hist(S[-1], bins=50)
plt.xlabel('index level')
plt.ylabel('frequency')
到期期权内在价值
# 模拟期权到期⽇的内在价值
%matplotlib inline
plt.hist(np.maximum(S[-1]-K, 0), bins=50)
plt.xlabel('option inner value')
plt.ylabel('frequency')

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。