一、模型思想
首先我们先弄清楚累计期权是指有多个交割日的期权，对于存续期为1个月的每日观察的累计期权而言，每个交易日都必须进行一次期权的行权，客户可能会产生亏损或收益也可能既不亏损也没收益。一个月下来的总的收益就是每一个交易日期权亏损的总和，并且在期权方面客户和金融机构是一个零和博弈，客户产生的收益即为金融机构的亏损，客户的亏损就是金融机构的收益。那么金融机构可以卖一个累计期权给客户的同时，对标的物进行一个反向操作，就可以进行实现delta对冲。
该模型是基于单一累计期权，进行一个月的模拟对冲，对未来可能出现的每天的价格进行预测，并且模拟次数达到一万次，查看对单一障碍期权做对冲可以实现的收益和收益率。理论上，当模拟次数足够多并且资产价格预测是服从随机游走时，对该期权进行对冲交易的收益的均值应该在0附近，收益率也几乎为0，也就是说正常情况下，100%的（不留敞口的）delta对冲是无法获得超额收益的。
此外，该模型也致力于预测当标的物价格的真实波动率变化时会怎样影响对冲收益和收益率。因为对一个月的标的物价格进行模拟，一共模拟了一万次，每一次模拟都要在服从随机游走的前提下生成22个资产价格（假设一个月有22个交易日）。那么一万次模拟的标的物价格数据就有各种各样不同的资产价格的真实波动率，那么我们统计出波动率在不同的区间，我们去做对冲可以实现怎样的收益和收益率？
二、代码实现
import numpy as np
from scipy.stats import norm
import pandas as pd
from tqdm import tqdm
def barrier_option_price(s, x, r, T, q, barrier, sigma, option_type, is_call):
"""
使用Black-Scholes模型计算障碍期权的价格。参数：
s: float, 标的资产当前价格
x: float, 期权的行权价
r: float, 无风险利率
T: float, 期权有效期的剩余时间，以年为单位
q: float, 标的资产的连续红利收益率
barrier: float, 期权的障碍水平
sigma: float, 标的资产的波动率
option_type: str, 期权类型（例如："UpAndIn"、"UpAndOut"、"DownAndIn"、"DownAndOut"）
is_call: bool, True表示看涨期权，False表示看跌期权返回：
float, 障碍期权的价格
"""def rubin_A(pv_S, x1, pv_X, Sig_Root_T, Phi):
return (Phi * (pv_S * norm.cdf(Phi * x1) - pv_X * norm.cdf(Phi * (x1 - Sig_Root_T))))def rubin_B(pv_S, pv_X, barrier, lambda0, Y, Sig_Root_T, Phi, eta):
lambda0 = min(lambda0, 500) # 确保lambda不超过500，避免数值不稳定
term1 = pv_S * (barrier ** (2 * lambda0)) * norm.cdf(eta * Y)
term2 = pv_X * (barrier ** (2 * lambda0 - 2)) * norm.cdf(eta * (Y - Sig_Root_T))
return (Phi * (term1 - term2))
# 计算中间变量
pv_X = x * np.exp(-r * T)
pv_S = s * np.exp(-q * T)
if T <= 0: # 检验剩余时间大于0
T = 0.00001
Sig_Root_T = sigma * np.sqrt(T)
D1 = np.log(pv_S / pv_X) / Sig_Root_T + 0.5 * Sig_Root_T
D2 = D1 - Sig_Root_T
mu = r - q - 0.5 * sigma * sigma
lambda0 = 1 + (mu / (sigma * sigma))
x1 = np.log(s / barrier) / Sig_Root_T + lambda0 * Sig_Root_T
Y = np.log((barrier**2/(s * x))) / Sig_Root_T +lambda0 * Sig_Root_T
y1 = np.log(barrier / s) / Sig_Root_T + lambda0 * Sig_Root_T
# 对期权进行分类分类顺序：
#1、是否看涨
#2、障碍价格是否大于执行价
#3、期权类型（向上敲入、向上敲出、向下敲入、向下敲出）
#4、标的物初始价格s是否大于障碍价格
price = 999 # 设置初始变量记录barrier_price
if is_call:
if (barrier >= x):
if option_type == "UpAndIn":
if s < barrier:
price = rubin_A(pv_S, x1, pv_X, Sig_Root_T, 1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, 1, -1) + rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, 1, -1)
else:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, 1)
elif option_type == "UpAndOut":
if s < barrier:
price = rubin_A(pv_S, x1, pv_X, Sig_Root_T, 1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, 1, -1) + rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, 1, -1)
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, 1) - price
else:
price = 0
elif option_type == "DownAndIn":
if s > barrier:
price = rubin_A(pv_S, x1, pv_X, Sig_Root_T, 1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, 1, 1)
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, 1) - price
else:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, 1)
elif option_type == "DownAndOut" :
if s > barrier:
price = rubin_A(pv_S, x1, pv_X, Sig_Root_T, 1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, 1, 1)
else:
price = 0
else: # 即当barrier < x
if option_type == "UpAndIn":
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, 1)
elif option_type == "UpAndOut":
price = 0
elif option_type == "DownAndIn":
if s > barrier:
price = rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, 1, 1)
else:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, 1)
elif option_type == "DownAndOut" :
if s > barrier:
price = rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, 1, 1)
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, 1) - price
else:
price = 0
else:
if (barrier >= x):
if option_type == "UpAndIn":
if s < barrier:
price = rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, -1, -1)
else:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, -1)
elif option_type == "UpAndOut":
if s < barrier:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, -1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, -1, -1)
else:
price = 0
elif option_type == "DownAndIn":
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, -1)
elif option_type == "DownAndOut" :
price = 0
else: # 即当barrier < x
if option_type == "UpAndIn":
if s < barrier:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, -1) - rubin_A(pv_S, x1, pv_X, Sig_Root_T, -1) + rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, -1, -1)
else:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, -1)
elif option_type == "UpAndOut":
if s < barrier:
price = rubin_A(pv_S, x1, pv_X, Sig_Root_T, -1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, -1, -1)
else:
price = 0
elif option_type == "DownAndIn":
if s > barrier:
price = rubin_A(pv_S, x1, pv_X, Sig_Root_T, -1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, -1, 1) + rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, -1, 1)
else:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, -1)
elif option_type == "DownAndOut" :
if s > barrier:
price = rubin_A(pv_S, D1, pv_X, Sig_Root_T, -1) - rubin_A(pv_S, x1, pv_X, Sig_Root_T, -1) + rubin_B(pv_S, pv_X, barrier / s, lambda0, Y, Sig_Root_T, -1, 1) - rubin_B(pv_S, pv_X, barrier / s, lambda0, y1, Sig_Root_T, -1, 1)
else:
price = 0
return price
def barrier_option_delta(s, x, r, T, q, barrier, sigma, option_type, is_call):
# todo 加其他的希腊值
"""
使用Black-Scholes模型计算障碍期权的delta。参数：
# 参数同上返回：
float, 障碍期权的delta
"""
price_up = barrier_option_price(s * 1.0001, x, r, T, q, barrier, sigma, option_type, is_call)
price_down = barrier_option_price(s * 0.9999, x, r, T, q, barrier, sigma, option_type, is_call)
delta = (price_up - price_down) / (s * 0.0002)return delta
# 标的物为期货的价格模拟
def simulation(S, T_remaining, r, sigma, I, pricelimitrate, steps):
"""
:param S: 初始价格
:param r: 无风险收益率
:param T_remaining: 剩余自然日（年）
:param sigma: 资产收益率的波动率
:param pricelimitrate: 涨跌幅限制比率
:param I: 路径次数(eg:I=300000,即模拟30万次)
:param steps: 每一次模拟的路径步数（应该与观察的次数保持一致）
"""delta_t = float(T_remaining) / steps
spath = np.zeros((steps + 1, I)) # 生成元素全为0的数组，行数为step+1，列数为I模拟次数
spath[0] = S # 设置数组第一行的的值全为1for t in range(1, steps + 1):
z = np.random.standard_normal(I) # 对于每一步都生成一个一维随机数组，且每一个元素都服从标准正态分布
# print(z)
middle1 = spath[t - 1, :] * np.exp(
(- 0.5 * sigma ** 2) * delta_t + sigma * np.sqrt(delta_t) * z) # 两个形状一样的一维数组，对应位置元素相乘
uplimit = spath[t - 1] * (1 + pricelimitrate) # 洗掉涨跌停板
lowlimit = spath[t - 1] * (1 - pricelimitrate)# 考虑涨跌幅度的限制，确定生成的价格在合理的范围内
temp = np.where(uplimit < middle1, uplimit, middle1) # 比较uplimit和middle1，取两者的最小，得到新的tem数组
temp = np.where(lowlimit > middle1, lowlimit, temp) # 比较lowlimit和temp，取两者的最大，得到新的tem
spath[t, :] = temp # 将新的tem赋值给第t行的价格数组spathreturn spath # 得到(step+1)*I的价格矩阵
# 基于蒙特卡洛的得到的标的资产价格矩阵计算payoff的现值
def cashflow_perday_observe(price_path, I, T, X, barrier, steps, r):
payoff = np.zeros(I)
for i in range(I):
knocked_out = False # 添加一个变量来跟踪是否发生敲出
for x in range(1, steps + 1):
if price_path[x, i] >= barrier:
knocked_out = True # 如果发生敲出，设置标志为True
payoff[i] = 0
break
if not knocked_out:
# 只有在没有敲出时，观察最后一天的payoff
payoff[i] = max(price_path[steps, i] - X, 0) * np.exp(-r * T)
return payoff.mean()
# 障碍期权公式法
S = 5000
X = 5000
r = 0.035
T=1/12 # 年的单位
q = 0
barrier = 5600
sigma = 0.2
option_type = 'UpAndOut'
is_call = True
barrier_price = barrier_option_price(S,X, r, T, q, barrier, sigma, option_type, is_call)
print(barrier_price)
打印累积期权的价格
T_remaining = 1/12
I = 10000
pricelimitrate = 0.1
steps = 21# 模拟的资产价格矩阵
spath = simulation(S, T_remaining, r, sigma, I, pricelimitrate, steps)
spath
输出模拟的资产价格矩阵
# 每日观察MC方法
# 相同参数下每日观察MC方法的期权价格
price = cashflow_perday_observe(spath, I, T, X, barrier, steps, r)
price
这一步是输出每日观察情况下的累积期权价格，对于敲出期权，每日观察更可能发生敲出，从而使得期权多头利润下降，故期权权利金低于到期日才观察的权利金。
# 到期日才观察MC方法
# 相同参数下MC方法到期日才观察的期权价格
def cashflow_expire(spath, I, T, X, barrier, steps, r):
return1 = np.zeros(I)
for i in range(I):
if spath[steps,i] > barrier: # 如果到期日的价格在X和barrier之间
return1[i] = 0
else: # 如果到期日的价格在执行价格X以下
return1[i] = (max(spath[steps,i] - X, 0)) * np.exp(-r * T)
return return1.mean()price1 = cashflow_expire(spath, I, T, X, barrier, steps, r)
print(price1)
这一步是输出到期日才观察规则下的期权价格
# 初始化数组记录每一次模拟、每一个时间节点的delta
delta_spath = np.zeros_like(spath)
for j in tqdm(range(I)): # 循环列
for i in range(steps + 1): # 循环行
delta_spath[i, j] = barrier_option_delta(spath[i, j], X, r, T, q, barrier, sigma, option_type, is_call)delta_spath
这一步是输出得到每一次模拟的delta矩阵
# 计算每一次模拟的资产价格收益率
returns = (spath[1:] - spath[:-1]) / spath[:-1]# 计算每一次模拟的收益率的波动率（日波动率）并年化
volatility_per_column = np.std(returns, axis=0) * np.sqrt(252)volatility_per_column
这一步是计算每一次模拟的收益率的波动率（日波动率）并年化，然后打印出来
volatility_per_column.mean()
一万次模拟的波动率的均值应该要和我们最初设定的波动率0.2较为接近，我算出来的是
0.19261179508827073，大家可以自己尝试一下
接下来开始对冲：对冲收益=期货端收益+期权端收益
delta_spath1 = delta_spath.copy()
# 这里对得到的delta矩阵的数据进行了处理：除了第0行的delta不变，其他各行的delta都要减去前一行的delta，如果是正值就表示要在期货市场上做多
# 如果是负值就表示要在期货市场上做空
for i in range(len(delta_spath1)):
if i > 0:
delta_spath1[i, :] = delta_spath1[i, :] - delta_spath1[i - 1, :]print(delta_spath1)
这一步输出新的delta矩阵，为后面的期货对冲多空方向和具体数量做准备
# 初始化数组用来记录每次模拟的期货端收益
returns_futures = np.zeros(I)
for i in range(I):
futures_payoff = 0
for j in range(steps):
# 这里可以这样写是因为不论是期货空头的卖方还是多头的卖方，都可以理解成买方，
# 比如多头的卖方可以理解成空头的买方，空头的卖方理解成多头的买方，这样就可以统一的与交割日当天的期货价格进行交割，
# 从而实现计算每一次对冲调仓的收益
futures_payoff += delta_spath1[j,i] * (spath[-1,i] - spath[j,i])returns_futures[i] = futures_payoff
returns_futures
这一步可以得到每一次模拟完期货端可以实现盈利的多少
根据我的计算，在一万次模拟对冲中，有5216次模拟在期货端是实现盈利的，并且对冲后期货端的收益的均值为0.33，几乎为0，这验证了我们的结论：
当模拟次数足够多并且资产价格预测是服从随机游走时，对该期权进行对冲交易的收益的均值应该在0附近，收益率也几乎为0，也就是说正常情况下，100%的（不留敞口的）delta对冲是无法获得超额收益的。
# 期权端收益
# 基于蒙特卡洛的得到的标的资产价格矩阵计算payoff的现值
def cashflow_perday_observe(price_path, I, T, X, barrier, steps, r):
payoff = np.zeros(I)
for i in range(I):
knocked_out = False # 添加一个变量来跟踪是否发生敲出
for x in range(1, steps + 1):
if price_path[x, i] >= barrier:
knocked_out = True # 如果发生敲出，设置标志为True
payoff[i] = 0
break
if not knocked_out:
# 只有在没有敲出时，观察最后一天的payoff
payoff[i] = max(price_path[steps, i] - X, 0) * np.exp(-r * T)
return payoffreturn_option = -1 * cashflow_perday_observe(spath, I, T, X, barrier, steps, r)
return_option
这一步返回期权端每一次模拟可以带来的收益
option_return.mean() #输出期权端获利的均值
return_hedge = returns_futures + return_option
return_hedge
这一步返回对冲收益：对冲收益=期权收益+期货收益
return_hedge.mean() # 输出对冲的平均收益，我算出来的结果是-92.37433200850103
然而我们前期计算过期权的价格了，期权的权利金我们可以收取92.70639636687538，所以算下来整个对冲可以实现的收益就只有0.33206435837435927元。总结来说，（10000次模拟）对冲后的平均收益为-92.37433200850103，然而期权权利金收入为92.70639636687538
平均起来对冲的收益为0.33约为0，也就是说百分百的纯对冲，收益的期望几乎为0，这是符合我们预期的。
那么最后一步就是，验证商品期货在不同波动率情况下，去做风险对冲实现的收益分别是如何的？
volatility_per_column # 先把每一次模拟路径的真实波动率打印出来
min(volatility_per_column),max(volatility_per_column) #得出最小值和最大值，我计算的分别是(0.09875264990561539, 0.3080380422127937)
# 对于不同波动率按区间划分
bins = np.arange(9,32) / 100
bins# 划分 volatility_per_column 的值到对应的区间，并计算每个区间中 return_hedge 的均值
volatility_bins = np.digitize(volatility_per_column, bins)
mean_return_per_bin = [return_hedge[volatility_bins == i].mean() for i in range(1, len(bins))]# 创建新的数组存储划分后的 volatility_per_column 和对应的均值
volatility_per_bin = [(bins[i-1] + bins[i]) / 2 for i in range(1, len(bins))]
mean_return_per_bin = np.array(mean_return_per_bin)# 绘制二维图形
plt.plot(volatility_per_bin, mean_return_per_bin, marker='o')
plt.xlabel('Volatility')
plt.ylabel('Mean Return')
plt.title('Mean Return vs Volatility in Bins')
plt.grid(True)
plt.show()
根据波动率图像我们可以看到，在不同的波动率上，我们进行delta对冲的收益的均值是不同的。图形显示当波动率小于0.175时，进行对冲的收益的均值是负数，也就是说当我们预测期货市场的波动率不大的时候（至少小于0.175），这时候大量的做对冲交易，最后很可能收益是负数。但当波动率在0.21-0.28之间时，我们大量的进行对冲交易的时候，我们对冲交易的均值就是正值，也就是可以盈利。当波动率大于0.29时，再100%基于delta对冲的收益的均值显著下降，这时候做大量的对冲交易可能会产生较大亏损。
该模型对模拟对冲交易可以提供一些指导性建议，比如当交易员预测波动率小于0.175时，这时的对冲交易的收益均值为负数，交易员可能会想办法减少对delta对冲的依赖，可能会引入gamma等其他的对冲参数，或者当真实波动率在0.21-0.28之间时，交易员可能会更加放心的做对冲，因为平均起来对冲的收益是正的，只要做的够多，那么收益的平均值就越可能大于0。
3825
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
