基于二叉树、蒙特卡洛模拟、BS方程的期权定价模型
点击蓝字
一
摘要：
本模型基于二叉树、蒙特卡洛模拟法、BS方程，对美式看涨期权(call)、美式看跌期权(put)、欧式看涨期权(Call)、欧式看跌期权(Put)进行了期权定价。这里我们使用了期权定价的二叉树模型和Black-Scholes-Merton定价公式来对期权定价，并使用蒙特卡洛模拟算法来模拟股票价格中的 （服从标准正态分布）的抽样，将得到的期权价值结果与使用二叉树、BSM定价公式得到的结果进行对比，并多次重复实验以减小误差，从而得到蒙特卡洛模拟后的期权价格。
【关键词】蒙特卡洛模拟、期权定价模型、二叉树模型、BSM方程
二
变量名称：
三
定价原理
1. 基于二叉树的期权定价模型
(1) 思想：假定到期且只有两种可能，而且涨跌幅均为10%的假设都很粗略。修改为：在T分为很多小的时间间隔Δt，而在每一个Δt，股票价格变化由S到Su或Sd。如果价格上涨概率为p，那么下跌的概率为1-p。
(2) u,p,d的确定：
我们假定市场为风险中性，也就是说股票预期收益率μ等于无风险利率r
这里我们假设：
又有：
并且股票的价格符合布朗运动：
得：
因此我们可以确定p、u、d的值：
令，
得：
（若这里有连续复利的股息：我们应该用r-q替换r）
在行权日T时刻，期权的价格为：
但考虑到折现的因素，在第t日，期权的价格为距离到期日T的预期收益折现值（若考虑连续股息率，需将r替换为r-q）
2. 基于蒙特卡洛模拟算法的期权定价模型
鉴于美式期权会在执行期到达前行权，因此常使用蒙特卡罗算法来估计欧式期权的价值以及标的资产的到期价格。
使用蒙卡模拟进行期权定价，需要在风险中性环境中（风险中性定价原则认为①标的资产的预期收益率=无风险利率，②从而对未来的现金流进行无风险利率折现；仅仅是概率空间的变换，没有提到投资者是风险中性的，因此不用对投资者的风险喜好进行假设）。
若考虑会有连续股息率为q，
则
在行权日T时刻，期权的价格为：
但考虑到折现的因素，在第t日，期权的价格为距离到期日T的预期收益折现值（若考虑连续股息率，需将r替换为r-q）
通过伊藤过程和维纳过程：我们得到ST的公式：
求得：
这里 ε 是从期望为0，标准差为1.0的正态分布中的抽样，接下来对 ε 进行蒙特卡罗模拟并代入S(T)方程中，即可得到结果。
对 ε 进行1000次蒙特卡罗模拟后得到 ε 的均值，将该均值作为 ε 的估计。按照样本的标准差，同时得到期权在T时刻的价格（欧式看涨期则为CT；欧式看跌期则为PT）
针对S0=50，K=52，r=0.04, q=0.01, T=2, n=2的欧式看涨期权
1） 蒙特卡洛模拟原理：
由蒙特卡洛算法的理论基础分析，蒙特卡洛算法是由随机变量X的样本X1,X2,…,XN的算术平均值：
作为我们所求解的近似值。而X1,X2,…,XN独立同分布，且具有有限期望值（E(X)＜∞）
由大数定律可知：
即：随机变量X的抽样样本均值XN，当样本量N足够大时，依概率收敛于它的期望值E(X)。
由伯努利大数定理，可以用k/n作为蒙特卡洛模拟得到的概率p的估计值，所以可将N个随机变量的值g(x1),g(x2),...,g(xN)的算术平均值：
作为概率p的估计值。
2） 求解误差：
考虑给定的置信度为1-α（0<α<1），设估计精度为 ε ，即：
成立。
记n次独立实验中A出现的次数为Kn，则Kn~B(n,p)，根据T-L中心极限定理知：
令，
即：
则标准正态分布的上侧分位数
解得：
由此可以看出，n 与 ε2 呈反比，即统计样本数量越大，误差 ε 越小。
3） 方差分析
若将误差 ε 定义为：
（α为置信度，N为样本数量，σ为标准差），则σ越小， ε 越小。（此处求得的误差 ε 为概率误差，且均方差σ是未知的）可以通过样本求得σ的估计值
来代替σ。
当实验不同次数后，可以得到不同实验次数N对应的σ值：
当N很小（0-100）时，σ的值偏差极大。
3. 基于BS方程的期权定价模型
我们可以使用下面的方程对欧式期权进行定价分析：
四
使用方法及效果展示
1. 按照提示数字输入参数：
(1) 若是欧式期权：
欧式期权输入参数示例
(2) 若是美式期权：
美式期权输入参数示例
2. 得到结果
(1) 欧式期权结果：
欧式期权期权价值预测结果
(2) 美式期权：
美式期权期权价值预测结果
3. 二叉树分析
（1） 美式期权：
美式看跌期权二叉树
注：图中行权代表将在此处美式期权会提前行权
该期权为四步二叉树，也可以为N步(N>1)
美式看涨期权二叉树
（2） 欧式期权
欧式看跌期权二叉树
注：图中MC为由蒙特卡洛模拟算出的期权价值
图中BS为由BS方程算出的期权价值
欧式看涨期权二叉树
参考文献
[1] John C.Hull. Options,Futures,and Other Derivatives Eighth Edition [M].北京：清华大学出版社，2020.
[2]毛啸峰. 期权定价数值方法的有效性分析[D].西南大学,2020.
五
定价程序解释
1. 二叉树计算Python代码：
1.# 欧式看跌#2.def europput(c_values, prices, K, steps, P, u,d,node,num):3. for i in range(1, steps + 1):4. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格5. c_values[i] = np.maximum(K - prices[i], 0) # 计算最后一列的期权价值6.7. node[num][0] = steps+18. node[num][1] = prices[i]9. node[num][2] = c_values[i]10. num = num + 1………………26.# 二叉树模型对欧式看涨期权定价27.def europcall(c_values, prices, K, steps, P, u,d,node,num):28. for i in range(1, steps + 1):29. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格30. c_values[i] = np.maximum(prices[i] - K,0)# 计算最后一列的期权价值31.32. node[num][0] = steps + 133. node[num][1] = prices[i]34. node[num][2] = c_values[i]35. num = num + 1………………45.# 美式看跌#46.def americaput(c_values, prices, K, steps, P, u,d,node,num):47. for i in range(1, steps + 1):48. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格49. c_values[i] = np.maximum(K - prices[i], 0)50.51. node[num][0] = steps+152. node[num][1] = prices[i]53. node[num][2] = c_values[i]54. num = num + 155.56. # print(prices[i])57. # print(c_values[i])58. for j in range(steps, 0, -1):59. for i in range(0, j):60. prices[i] = prices[i + 1] * d………………75.# 美式看涨#76.def americacall(c_values, prices, K, steps, P, u,d,node,num):77. for i in range(1, steps + 1):78. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格79. c_values[i] = np.maximum(prices[i] - K, 0)80.81. node[num][0] = steps+182. node[num][1] = prices[i]83. node[num][2] = c_values[i]84. num = num + 185.86. # print(prices[i])87. # print(c_values[i])88. for j in range(steps, 0, -1):89. for i in range(0, j):90. prices[i] = prices[i + 1] * d………………
2. 蒙特卡洛模拟算法Python代码：
1.def MTKL(t,r,q,sigma,S,K,kind):2. miu=r-q #miu是风险中性世界下的期望收益率3. ans_ST=[]4. ans_ST.append(S)5. ans_EU=[]6. for i in range(1000):7. epsilon = np.random.normal(0, 1) #模拟epsilon8. ST=S*np.exp((miu-sigma*sigma/2)*t + sigma*epsilon* math.sqrt(t))9. ans_ST.append(ST)………………
3. BS方程计算欧式期权代码：
1.def BS(S,K,T,r,v,q,kind):2. d1 = (np.log(S / K) + (r + v ** 2 / 2.) * T) / (v * np.sqrt(T))3. d2 = d1 - v * np.sqrt(T)4. if kind == 3:………………
4. 绘制二叉树、BS和MC值Python代码：
1.def show_graph(node,A,steps,u,d,kind,len,mtkl=0,bs=0):2.3. npNode=np.array(node)4. x_step=npNode[:,0]5. y_price=npNode[:,1]6.7. #x=(0,steps+1,1)8. #y=[]9. #min=A*d**steps10. #for i in range(2*steps+1):11. # y.append(min*u**i)12. plt.rcParams['font.family']='Kaiti'13. plt.rcParams['font.size'] = 2014. plt.plot(x_step,y_price,'*')15. plt.xlabel("时间阶段")16. plt.ylabel("标的资产价格（单位：美元）")17.18. #plt.legend(x_step,color='blue', linewidth = 2.5, linestyle = "-", label ='cosine')19. #plt.legend(y_price,color='red', linewidth = 2.5, linestyle = "-", label ='sine')20. #plt.legend(loc='upper left’)21.22. for i in range(len): #i从1到2到3,i是横坐标x23. for j in range(len): #j代表的是npNode数组的索引24. if(i+1==npNode[j][0]): #i与step重合的时候s25. x_src=npNode[j][0]26. y_src=npNode[j][1]27. font=npNode[j][1]28. for jj in range(len):29. if(i+2==npNode[jj][0]): #i与step+1重合的时候30. if(round(npNode[jj][1])==round(font*d) or round(npNode[jj][1])==round(font*u)):31. x_dest = npNode[jj][0]32. y_dest = npNode[jj][1]33. plt.annotate('',xy=(x_dest,y_dest),xytext=(x_src,y_src),arrowprops=dict(facecolor='black',shrink=0.1,width=0.3))34.35. for i in range(len):36. #print(i)37. show_1 = str(round(npNode[i][1],4))38. show_2 = str(round(npNode[i][2],4))39. #print(npNode[i][0])………………
5. 主函数代码：
1.if __name__ == '__main__':2.3. S = float(input("请输入标的资产（股票）的初始价格（美元）: "))4. K = float(input("请输入该期权的执行价格（美元）:"))5. r = float(input("请输入一年期无风险利率:"))6. q = float(input("请输入股息收益率（若无分红，股息收益率为0）:"))7. sigma = float(input("请输入标的资产（股票）价格波动率:")) # 好像是自己来求？8. t = float(input("请输入总时间（年）:"))9. steps = int(input("请输入该二叉树的步数（步）:"))10. kind = int(input("请输入期权属性：1：美式看涨 2：美式看跌 3：欧式看涨 4：欧式看跌"))11.12. MTKL_kind=-113. BS_kind =-114. if kind== 3 or kind == 4:15. MTKL_kind=int(input("是否选择附加使用蒙特卡洛模拟求期权价值：5：是 6：否"))16. BS_kind = int(input("是否选择附加使用BS定价模型求期权价值：7：是 8：否"))17.18. len=int((steps+2)*(steps+1)/2) #总结点数19. node = [[0 for col in range(3)] for row in range(len)] # 一个n*3的矩阵，-1代表未初始化20. num=021. #node[0]:steps node[1]:price node[2]:value#22.23. u = np.exp(sigma * np.sqrt(t / steps)) # 时间间隔为△t=t/steps24. d = 1 / u25. P = (np.exp((r - q) * t / steps) - d) / (u - d)26. prices = np.zeros(steps + 1) # 生成最后一列的股票价格空数组27. c_values = np.zeros(steps + 1) # 生成最后一列的期权价值空数组28. prices[0] = S * d ** steps # 最后一行最后一列的股票价格29. c_values[0] = np.maximum(K - prices[0], 0) # 最后一行最后一列的期权价值30.31. node[num][0]=steps+132. node[num][1]=prices[0]33. node[num][2]=c_values[0]34. num=num+1………………
6. 程序总代码：
1.import numpy as np2.from networkx import DiGraph3.import networkx as nx4.import matplotlib.pyplot as plt5.import math6.from scipy.stats import norm7.import numpy as np8.n = norm.pdf # 概率密度函数9.N = norm.cdf # 累计概率函数10.11.#G: DiGraph = nx.DiGraph() # 有向有权图12.#num = 013.14.# 欧式看跌#15.def europput(c_values, prices, K, steps, P, u,d,node,num):16. for i in range(1, steps + 1):17. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格18. c_values[i] = np.maximum(K - prices[i], 0) # 计算最后一列的期权价值19.20. node[num][0] = steps+121. node[num][1] = prices[i]22. node[num][2] = c_values[i]23. num = num + 1………………39.# 二叉树模型对欧式看涨期权定价40.def europcall(c_values, prices, K, steps, P, u,d,node,num):41. for i in range(1, steps + 1):42. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格43. c_values[i] = np.maximum(prices[i] - K,0)# 计算最后一列的期权价值44.45. node[num][0] = steps + 146. node[num][1] = prices[i]47. node[num][2] = c_values[i]48. num = num + 1………………58.# 美式看跌#59.def americaput(c_values, prices, K, steps, P, u,d,node,num):60. for i in range(1, steps + 1):61. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格62. c_values[i] = np.maximum(K - prices[i], 0)63.64. node[num][0] = steps+165. node[num][1] = prices[i]66. node[num][2] = c_values[i]67. num = num + 1………………88.# 美式看涨#89.def americacall(c_values, prices, K, steps, P, u,d,node,num):90. for i in range(1, steps + 1):91. prices[i] = prices[i - 1] * (u ** 2) # 计算最后一列的股票价格92. c_values[i] = np.maximum(prices[i] - K, 0)93.94. node[num][0] = steps+195. node[num][1] = prices[i]96. node[num][2] = c_values[i]97. num = num + 1…………………119.def show_graph(node,A,steps,u,d,kind,len,mtkl=0,bs=0):120.121. npNode=np.array(node)122. x_step=npNode[:,0]123. y_price=npNode[:,1]124.125. #x=(0,steps+1,1)126. #y=[]127. #min=A*d**steps128. #for i in range(2*steps+1):129. # y.append(min*u**i)130. plt.rcParams['font.family']='Kaiti'131. plt.rcParams['font.size'] = 20132. plt.plot(x_step,y_price,'*')133. plt.xlabel("时间阶段")134. plt.ylabel("标的资产价格（单位：美元）")135.136. #plt.legend(x_step,color='blue', linewidth = 2.5, linestyle = "-", label ='cosine')137. #plt.legend(y_price,color='red', linewidth = 2.5, linestyle = "-", label ='sine')138. #plt.legend(loc='upper left’)139.140. for i in range(len): #i从1到2到3,i是横坐标x141. for j in range(len): #j代表的是npNode数组的索引142. if(i+1==npNode[j][0]): #i与step重合的时候s143. x_src=npNode[j][0]144. y_src=npNode[j][1]145. font=npNode[j][1]146. for jj in range(len):………………193.def MTKL(t,r,q,sigma,S,K,kind):194. miu=r-q #miu是风险中性世界下的期望收益率195. ans_ST=[]196. ans_ST.append(S)197. ans_EU=[]198. for i in range(1000):199. epsilon = np.random.normal(0, 1) #模拟epsilon200. ST=S*np.exp((miu-sigma*sigma/2)*t + sigma*epsilon* math.sqrt(t))201. ans_ST.append(ST)………………211.def BS(S,K,T,r,v,q,kind):212. d1 = (np.log(S / K) + (r + v ** 2 / 2.) * T) / (v * np.sqrt(T))213. d2 = d1 - v * np.sqrt(T)214. if kind == 3:215. V = S * np.exp(-q * T) * N(d1) - K * np.exp(-r * T) * N(d2)216. if kind == 4 :217. V = K * np.exp(-r * T) * N(-d2) - S * np.exp(-q * T) * N(-d1)218. return V219.220.if __name__ == '__main__':221.222. S = float(input("请输入标的资产（股票）的初始价格（美元）: "))223. K = float(input("请输入该期权的执行价格（美元）:"))224. r = float(input("请输入一年期无风险利率:"))225. q = float(input("请输入股息收益率（若无分红，股息收益率为0）:"))226. sigma = float(input("请输入标的资产（股票）价格波动率:")) # 好像是自己来求？227. t = float(input("请输入总时间（年）:"))228. steps = int(input("请输入该二叉树的步数（步）:"))229. kind = int(input("请输入期权属性：1：美式看涨 2：美式看跌 3：欧式看涨 4：欧式看跌"))230.231. MTKL_kind=-1232. BS_kind =-1233. if kind== 3 or kind == 4:234. MTKL_kind=int(input("是否选择附加使用蒙特卡洛模拟求期权价值：5：是 6：否"))235. BS_kind = int(input("是否选择附加使用BS定价模型求期权价值：7：是 8：否"))236.237. len=int((steps+2)*(steps+1)/2) #总结点数238. node = [[0 for col in range(3)] for row in range(len)] # 一个n*3的矩阵，-1代表未初始化239. num=0240. #node[0]:steps node[1]:price node[2]:value#241.242. u = np.exp(sigma * np.sqrt(t / steps)) # 时间间隔为△t=t/steps243. d = 1 / u244. P = (np.exp((r - q) * t / steps) - d) / (u - d)245. prices = np.zeros(steps + 1) # 生成最后一列的股票价格空数组246. c_values = np.zeros(steps + 1) # 生成最后一列的期权价值空数组………………
注：应作者要求只展示部分代码，省略号处表示代码未完。
APC编辑部科普组
往期高赞好文：
1.
2.
3.
4.
5.
A.P.C.总部招生群:1027526691
格物社•科普平台：605923025
我们是传播科普的学生大家庭，
很高兴你能看到这里~
责任编辑：一毫秒的永恒
APC科普编辑部
👈扫码加入总部招生群
你的每个在看，我都当了喜欢