写在前面，Tushare的使用链接：https://tushare.pro/register?reg=438650
一、微笑波动率曲线成因分析
根据6月30日的市场数据，以沪深300ETF期权为例，作出如下微笑波动率曲面。
在BS期权定价模型假设下，期权波动率通常为常数，然而这与实际情况有较大差别。如上表1，通过传统BS期权定价模型计算出来的隐含波动率呈现出一种被称为“波动率微笑”的现象。具有相同到期日和标的资产而执行价格不同的期权，其执行价格偏离标的资产现货价格越远，隐含波动率越大。具有相同执行价格和标的资产而到期日不同的期权，其执行价到期日越近，隐含波动率越大。原因如下：
资产价格非正态分布
标准BS模型假定标的资产价格服从对数正态分布，收益率服从正态分布。但是大量实证检验发现，现实市场中，金融资产的收益率分布更加显示出尖峰肥尾的特征。这种分布下，收益率出现极端值的概率高于正态分布，而在公式中采用收益率正态分布的前提假设，会低估了到期时期权价值变为实值与虚值出现的概率，相应也低估了深度实值和深度虚值期权的价格。资产价格跳跃
BS模型采用的是风险中性定价，并假设资产价格服从带漂移项的布朗过程，忽略了现实市场上资产价格在一定冲击下发生跳跃的可能。例如价格在期权临近到期前发生跳跃，且空方根据变化后的价格调整标的资产头寸并持有到期，到期时复制组合与期权价值将可能出现较大偏差，使得期权空方面临额外风险。这种风险无法分散化，空方必须要求相应补偿，造成期权市场价格对理论价格的溢价。Gamma风险和Vega风险
期权空方的Delta套期保值中还面临Camma风险和Vega风险。Gamma是期权价格对标定资产价格的二阶导数，Gamma越大，Delta对标的资产价格的变动越敏感，复制组合价值越容易偏离理论期权价值。期权空方因交易成本无法连续调整标定资产头寸时，Gamma风险也就不可避免。Vega是期权价格对标的资产波动率的一阶导数。标的资产的实际波动率并非BS模型假设的常数，而是一个随机变量。当市场其他条件不变时，Vega越大，期权理论价值越容易发生变化，复制组合的价值也越容易与其发生偏离，导致Vega风险。Gamma风险和Vega风险均可通过引入同种标的资产、相同执行价格但期限较短的期权来进行套期保值，因组合中增加了新的期权，还需要相应重新调整标的资产头寸。这都增加了期权空方的Delta套期保值成本。标的资产价格预期
隐含波动率的“微笑”现象与市场对标的资产未来价格走势的预期有关。假定某标的资产的当前价格为S0，市场预期将很快下跌至S1。此时任何执行价格下的看跌期权价格均上升，看涨期权价格均下降。其中执行价格在S1到S0之间的看跌期权因将从虚值转为实值，期权卖方面临的风险将更大，Delta套期保值的成本更高，价格上升最多。同理，执行价格在S1到S0之间的看涨期权因将从实值转为虚值，价格下跌最多。可见，在当前资产价格尚未发生变动的条件下，虚值看跌期权的隐含波动率上升幅度大于实值看跌期权，实值看涨期权的隐含波动率下降幅度大于虚值看涨期权，体现在波动率“微笑”曲线上，均表现为曲线的左半部分高于右半部分。期权市场溢价
从理论上来看，期权从平值状态变为实值状态和虚值状态的概率应该基本相同，并且在平值状态时其时间价值最大。深度实值期权的Delta接近1，在投资中的杠杆作用最大，相应市场需求量很大。但是除非投资者预期标定资产的价格会有一个根本性的变动，一般不会出售深度实值期权，因此，供给量较小。溢价期权、折价期权分别处于实值和虚值状态，其带给投资者未来较大收益的可能性比平值期权要小，其时间价值也会比平值期权小。深度实值期权的溢价较高，其隐含波动率也较高。对相同执行价格的看涨期权和看跌期权，当一个处于深度实值状态时，另一个必然处于深度虚值状态。根据看涨看跌平价关系，这两个期权的波动率应当大致相同。可见，实值看涨期权的溢价也会造成虚值看跌期权的溢价，造成隐含波动率的“微笑”。标的资产和期权交易成本
标定资产的交易成本是期权空方Delta套期保值额外成本的重要来源之一。在保值成本增加相同的条件下，深度实值和深度虚值期权的隐含波动率增加更多，呈现出隐含波动率“微笑”曲线。期权本身也存在交易成本。深度实值和深度虚值期权的流动性较差，交易成本也较大，这个效应通过期权的Gamma风险保值，可引发波动率“微笑”。平价期权的Gamma风险最大，如只用标的资产保值，其头寸调整最为频繁，引致的额外成本最大。但是另外两个效应减轻了这个影响：第一，平价期权的Gamma随时间单调衰减的速度非常快，即Gamma风险下降的速度很快；第二，平价期权可利用短期平价期权保值，后者的交易成本相对较低。相比之下，短期深度实值和深度虚值期权的交易成本较高，深度实值和深度虚值期权更倾向于频繁调整标定资产头寸，受标的资产交易成本的影响更大。可见，交易成本对较长期深度实值和深度虚值期权的影响是双重的，大于对平价期权的影响。交易成本不对称
做市商机制下有可能出现买卖价差的不对称。市场上买方力量较强时，做市商在卖出上收取较多佣金，在买入时收取较少佣金，真实的资产价格将更接近买入价而非买入卖出的中间价。卖方力量较强时正好相反。作为一种保值工具，期权市场需求通常大于市场供给，在做市商的作用下，更易出现买卖价差的不对称，且其真实价格更接近买入价。实证研究中多采用中间价，高估了期权价格，也高估了隐含波动率。对深度实值和深度虚值期权，因其买卖差价更大，中间价对实际价格的高估更多。报价机制和价格误差
资产的实际价格是连续变动的，但在市场交易中标的资产、期权和无风险利率都存在报价的最小间距，造成实际价格与报价之间存在误差。资产的真实价格应在当前报价水平上加上或减去一个微小的差价。此外，标的资产与期权报价间的不同步、当前指数价格或者成分股当前价格的不同步等也会影响期权价格的准确性。
二、期权定价模拟：利用蒙特卡洛模拟估计欧式看涨期权的价值
同样，以沪深300ETF期权为标的，利用蒙特卡洛方法为其进行模拟定价，分别选取1000，5000，10000频次，主要参考Black-Scholes-Merton(BSM)模型，在模型中到期的指数水平是一个随机变量，通过到期指数公式计算：
ST=S0exp((r−12σ2)T+σT‾‾√z)ST=S0exp((r−12σ2)T+σTz)
其中z是一个服从标准正态分布的随机变量．
蒙特卡洛算法的描述：
１．生成 I 个服从标准正态分布的随机数 i={1,2,3,I}, I为随机模拟的次数
２．通过上面的到期指数计算公式，计算出所有模拟结果的到期指数ST(i)ST(i)
３．计算到期期权的每一个模拟可能的内在价值hT(i)hT(i)
hT(i)=max(ST(i)−K,0)hT(i)=max(ST(i)−K,0)
４．通过蒙特卡洛估算公式计算出期权现值
C0≈e−rT1T∑IihT(i)
代码附录：
1.微笑波动率曲线成因分析
1）参数选择：
参数 tushare函数 tushare参数
执行价格（K） opt_basic exercise_price
到期时间（t） opt_basic delist_date
标的资产价格（S） fund_daily close
利率（rf） shibor 1w
期权最新价格（c） opt_daily settle
2）代码：
#设置plotly用户名和API-Keyimport chart_studiochart_studio.tools.set_credentials_file(username='EasyLake', api_key='GrIorDtlVNRQjFlclvdH')#导入库import pandas as pdimport datetimeimport tushare as tsimport numpy as npfrom math import log,sqrt,expfrom scipy import statsimport plotly.graph_objects as goimport plotlyimport plotly.express as pximport matplotlib.pyplot as plt#调用apipro = ts.pro_api()plotly.offline.init_notebook_mode(connected=True)#提取数据并清洗def extra_data(date): # 提取数据# 提取沪深300ETF期权合约基础信息df_basic = pro.opt_basic(exchange='SSE', fields='ts_code,name,call_put,exercise_price,list_date,delist_date')df_basic = df_basic.loc[df_basic['name'].str.contains('沪深300ETF')]df_basic = df_basic[(df_basic.list_date<=date)&(df_basic.delist_date>date)] # 提取当天市场上交易的期权合约df_basic = df_basic.drop(['name','list_date'],axis=1)df_basic['date'] = date# 提取日线行情数据df_cal = pro.trade_cal(exchange='SSE', cal_date=date, fields = 'cal_date,is_open,pretrade_date')if df_cal.iloc[0, 1] == 0:date = df_cal.iloc[0, 2] # 判断当天是否为交易日，若否则选择前一个交易日opt_list = df_basic['ts_code'].tolist() # 获取沪深300ETF期权期权合约列表df_daily = pro.opt_daily(trade_date=date,exchange = 'SSE',fields='ts_code,trade_date,settle')df_daily = df_daily[df_daily['ts_code'].isin(opt_list)]# 提取沪深300ETF指数数据df_300etf = pro.fund_daily(ts_code='510300.SH', trade_date = date,fields = 'close')s = df_300etf.iloc[0, 0] # 提取无风险利率数据（用一周shibor利率表示）df_shibor = pro.shibor(date = date,fields = '1w')rf = df_shibor.iloc[0,0]/100# 数据合并df = pd.merge(df_basic,df_daily,how='left',on=['ts_code'])df['s'] = sdf['r'] = rfdf = df.rename(columns={'exercise_price':'k', 'settle':'c'})#print(df)return dfdef data_clear(df): # 数据清洗def days(df): # 计算期权到期时间start_date = datetime.datetime.strptime(df.date,"%Y%m%d")end_date = datetime.datetime.strptime(df.delist_date,"%Y%m%d")delta = end_date - start_datereturn int(delta.days)/365def iq(df): # 计算隐含分红率#q = -log((df.settle+df.exercise_price*exp(-df.interest*df.delta)-df.settle_p)/(df.s0))/df.deltaq = -log((df.c+df.k*exp(-df.r*df.t)-df.c_p)/(df.s))/df.treturn qdf['t'] = df.apply(days,axis = 1)df = df.drop(['delist_date','date'],axis = 1)# 计算隐含分红率df_c = df[df['call_put']=='C']df_p = df[df['call_put']=='P']df_p = df_p.rename(columns={'c':'c_p','ts_code':'ts_code_p','call_put':'call_put_p'})df = pd.merge(df_c,df_p,how='left',on=['trade_date','k','t','r','s'])df['q'] = df.apply(iq,axis = 1)c_list = [x for x in range(8)]+[11]df_c = df.iloc[:,c_list]df_p = df[['ts_code_p','trade_date','c_p','call_put_p','k','t','r','s','q']]df_p = df_p.rename(columns={'c_p':'c','ts_code_p':'ts_code','call_put_p':'call_put'})df_c = df_c.append(df_p)#print(df_c)return df_c #计算隐含波动率#根据BS公式计算期权价值def bsm_value(s,k,t,r,sigma,q,option_type):d1 = ( log( s/k ) + ( r -q + 0.5*sigma**2 )*t )/( sigma*sqrt(t) )d2 = ( log( s/k ) + ( r -q - 0.5*sigma**2 )*t )/( sigma*sqrt(t) )if option_type.lower() == 'c':value = (s*exp(-q*t)*stats.norm.cdf( d1) - k*exp( -r*t )*stats.norm.cdf( d2))else:value = k * exp(-r * t) * stats.norm.cdf(-d2) - s*exp(-q*t) * stats.norm.cdf(-d1)return value#二分法求隐含波动率def bsm_imp_vol(s,k,t,r,c,q,option_type):c_est = 0 # 期权价格估计值top = 1 #波动率上限floor = 0 #波动率下限sigma = ( floor + top )/2 #波动率初始值count = 0 # 计数器while abs( c - c_est ) > 0.00001:c_est = bsm_value(s,k,t,r,sigma,q,option_type) #根据价格判断波动率是被低估还是高估，并对波动率做修正count += 1 if count > 50: # 时间价值为0的期权是算不出隐含波动率的，因此迭代到一定次数就不再迭代了sigma = 0breakif c - c_est > 0: #f(x)>0floor = sigmasigma = ( sigma + top )/2else:top = sigmasigma = ( sigma + floor )/2return sigma def cal_iv(df): # 计算主程序option_list = df.ts_code.tolist()df = df.set_index('ts_code')alist = []for option in option_list:s = df.loc[option,'s']k = df.loc[option,'k']t = df.loc[option,'t']r = df.loc[option,'r']c = df.loc[option,'c']q = df.loc[option,'q']option_type = df.loc[option,'call_put']sigma = bsm_imp_vol(s,k,t,r,c,q,option_type)alist.append(sigma)df['iv'] = alistreturn df def data_pivot(df): # 数据透视df = df.reset_index()option_type = 'C' # 具有相同执行价格、相同剩余到期时间的看涨看跌期权隐含波动率相等，因此算一个就够了df = df[df['call_put']==option_type]df = df.drop(['ts_code','trade_date','c','s','r','call_put','q'],axis=1)df['t'] = df['t']*365df['t'] = df['t'].astype(int)df = df.pivot_table(index=["k"],columns=["t"],values=["iv"])df.columns = df.columns.droplevel(0)df.index.name = Nonedf = df.reset_index()df = df.rename(columns={'index':'k'})return dfdef fitting(df): # 多项式拟合 col_list = df.columnsfor i in range(df.shape[1]-1):x_col = col_list[0]y_col = col_list[i+1]df1 = df.dropna(subset=[y_col])x = df1.iloc[:,0]y = df1.iloc[:,i+1]degree = 2weights = np.polyfit(x, y, degree)model = np.poly1d(weights)predict = np.poly1d(model)x_given_list = df[pd.isnull(df[y_col]) == True][x_col].tolist()# 所有空值对应的k组成列表for x_given in x_given_list:y_predict = predict(x_given)df.loc[df[x_col]==x_given, y_col] = y_predictreturn dfdef im_surface(df): # 波动率曲面作图# df = plot_df()df = fitting(df) #df.to_excel('iv_fitting.xlsx')df = df.set_index('k')y = np.array(df.index)x = np.array(df.columns)fig = go.Figure(data=[go.Surface(z=df.values, x=x, y=y)])fig.update_layout(scene = dict(xaxis_title='剩余期限',yaxis_title='执行价格',zaxis_title='隐含波动率'),width=1400,margin=dict(r=200, b=100, l=100, t=100))#fig.write_image("fig1.jpg")plotly.offline.plot(fig)def smile_plot(df): # 波动率微笑作图# df = plot_df()df = df.set_index('k')df = df.stack().reset_index()df.columns = ['k', 'days', 'iv']fig = px.line(df, x="k", y="iv", color="days",line_shape="spline")plotly.offline.plot(fig)def main():date = '0630'# plot_df()df = extra_data(date) # 提取数据df = data_clear(df) # 数据清洗df = cal_iv(df) # 计算隐含波动率df = data_pivot(df) # 数据透视表smile_plot(df) # 波动率微笑im_surface(df) # 波动率曲面if __name__ == '__main__':main()2.期权定价模拟1）参数选择：参数tushare函数tushare参数执行价格（K）opt_basicexercise_price到期时间（t）opt_basicdelist_date标的资产价格（S）fund_dailyclose利率（rf）shibor1w期权最新价格（c）opt_dailysettle3）代码：#设置plotly用户名和API-Keyimport chart_studiochart_studio.tools.set_credentials_file(username='EasyLake', api_key='GrIorDtlVNRQjFlclvdH')#导入库import pandas as pdimport datetimeimport tushare as tsimport numpy as npfrom math import log,sqrt,expfrom scipy import statsimport plotly.graph_objects as goimport plotlyimport plotly.express as pximport matplotlib.pyplot as plt#调用apipro = ts.pro_api()plotly.offline.init_notebook_mode(connected=True)# 提取数据并清洗def extra_data(date): # 提取数据# 提取沪深300ETF期权合约基础信息df_basic = pro.opt_basic(exchange='SSE', fields='ts_code,name,call_put,exercise_price,list_date,delist_date')df_basic = df_basic.loc[df_basic['name'].str.contains('沪深300ETF')]df_basic = df_basic[(df_basic.list_date<=date)&(df_basic.delist_date>date)] # 提取当天市场上交易的期权合约df_basic = df_basic.drop(['name','list_date'],axis=1)df_basic['date'] = date# 提取日线行情数据df_cal = pro.trade_cal(exchange='SSE', cal_date=date, fields = 'cal_date,is_open,pretrade_date')if df_cal.iloc[0, 1] == 0:date = df_cal.iloc[0, 2] # 判断当天是否为交易日，若否则选择前一个交易日opt_list = df_basic['ts_code'].tolist() # 获取沪深300ETF期权期权合约列表df_daily = pro.opt_daily(trade_date=date,exchange = 'SSE',fields='ts_code,trade_date,settle')df_daily = df_daily[df_daily['ts_code'].isin(opt_list)]# 提取沪深300ETF指数数据df_300etf = pro.fund_daily(ts_code='510300.SH', trade_date = date,fields = 'close')s = df_300etf.iloc[0, 0] # 提取无风险利率数据（用一周shibor利率表示）df_shibor = pro.shibor(date = date,fields = '1w')rf = df_shibor.iloc[0,0]/100# 数据合并df = pd.merge(df_basic,df_daily,how='left',on=['ts_code'])df['s'] = sdf['r'] = rfdf = df.rename(columns={'exercise_price':'k', 'settle':'c'})#print(df)return dfdef data_clear(df): # 数据清洗def days(df): # 计算期权到期时间start_date = datetime.datetime.strptime(df.date,"%Y%m%d")end_date = datetime.datetime.strptime(df.delist_date,"%Y%m%d")delta = end_date - start_datereturn int(delta.days)/365def iq(df): # 计算隐含分红率#q = -log((df.settle+df.exercise_price*exp(-df.interest*df.delta)-df.settle_p)/(df.s0))/df.deltaq = -log((df.c+df.k*exp(-df.r*df.t)-df.c_p)/(df.s))/df.treturn qdf['t'] = df.apply(days,axis = 1)df = df.drop(['delist_date','date'],axis = 1)# 计算隐含分红率df_c = df[df['call_put']=='C']df_p = df[df['call_put']=='P']df_p = df_p.rename(columns={'c':'c_p','ts_code':'ts_code_p','call_put':'call_put_p'})df = pd.merge(df_c,df_p,how='left',on=['trade_date','k','t','r','s'])df['q'] = df.apply(iq,axis = 1)c_list = [x for x in range(8)]+[11]df_c = df.iloc[:,c_list]df_p = df[['ts_code_p','trade_date','c_p','call_put_p','k','t','r','s','q']]df_p = df_p.rename(columns={'c_p':'c','ts_code_p':'ts_code','call_put_p':'call_put'})df_c = df_c.append(df_p)#print(df_c)return df_c #计算隐含波动率#根据BS公式计算期权价值def bsm_value(s,k,t,r,sigma,q,option_type):d1 = ( log( s/k ) + ( r -q + 0.5*sigma**2 )*t )/( sigma*sqrt(t) )d2 = ( log( s/k ) + ( r -q - 0.5*sigma**2 )*t )/( sigma*sqrt(t) )if option_type.lower() == 'c':value = (s*exp(-q*t)*stats.norm.cdf( d1) - k*exp( -r*t )*stats.norm.cdf( d2))else:value = k * exp(-r * t) * stats.norm.cdf(-d2) - s*exp(-q*t) * stats.norm.cdf(-d1)return value#二分法求隐含波动率def bsm_imp_vol(s,k,t,r,c,q,option_type):c_est = 0 # 期权价格估计值top = 1 #波动率上限floor = 0 #波动率下限sigma = ( floor + top )/2 #波动率初始值count = 0 # 计数器while abs( c - c_est ) > 0.00001:c_est = bsm_value(s,k,t,r,sigma,q,option_type) #根据价格判断波动率是被低估还是高估，并对波动率做修正count += 1 if count > 50: # 时间价值为0的期权是算不出隐含波动率的，因此迭代到一定次数就不再迭代了sigma = 0breakif c - c_est > 0: #f(x)>0floor = sigmasigma = ( sigma + top )/2else:top = sigmasigma = ( sigma + floor )/2return sigma def MC(s,k,t,r,sigma):I = 50# 随机数的个数，模拟计算到期指数的次数z = np.random.standard_normal(I)# 股票到期指数水平S_T = s * np.exp((r - 0.5 * sigma ** 2)* t + sigma * np.sqrt(t) * z)# 期权的内在价值h_T = np.maximum(S_T - k, 0)# 期权现值C = np.exp(-r * t)* np.sum(h_T)/I# 期权现值print ('value of the European call option %.6f'%C)print('行权价=%.6f'%k)def draw_pic(position,y,label,plt):plt.subplot(position)plt.plot(y, label=label)plt.grid(True)plt.legend(loc=0)plt.ylabel(label) %matplotlib inlineplt.figure(figsize=(12,7))draw_pic(211, S_T, 'S_T', plt)draw_pic(212, h_T, 'h_T', plt)return Cdef cal_mc(df): # 计算主程序option_list = df.ts_code.tolist()df = df.set_index('ts_code')alist = []for option in option_list:s = df.loc[option,'s']k = df.loc[option,'k']t = df.loc[option,'t']r = df.loc[option,'r']c = df.loc[option,'c']q = df.loc[option,'q']option_type = df.loc[option,'call_put']sigma = bsm_imp_vol(s,k,t,r,c,q,option_type)C=MC(s,k,t,r,sigma)alist.append(C)df['C'] = alistreturn df def data_pivot(df): # 数据透视df = df.reset_index()option_type = 'C' # 具有相同执行价格、相同剩余到期时间的看涨看跌期权隐含波动率相等，因此算一个就够了df = df[df['call_put']==option_type]df = df.drop(['ts_code','trade_date','c','s','r','call_put','q'],axis=1)df['t'] = df['t']*365df['t'] = df['t'].astype(int)df = df.pivot_table(index=["k"],columns=["t"],values=["C"])df.columns = df.columns.droplevel(0)df.index.name = Nonedf = df.reset_index()df = df.rename(columns={'index':'k'})return dfdef fitting(df): # 多项式拟合 col_list = df.columnsfor i in range(df.shape[1]-1):x_col = col_list[0]y_col = col_list[i+1]df1 = df.dropna(subset=[y_col])x = df1.iloc[:,0]y = df1.iloc[:,i+1]degree = 2weights = np.polyfit(x, y, degree)model = np.poly1d(weights)predict = np.poly1d(model)x_given_list = df[pd.isnull(df[y_col]) == True][x_col].tolist()# 所有空值对应的k组成列表for x_given in x_given_list:y_predict = predict(x_given)df.loc[df[x_col]==x_given, y_col] = y_predictreturn dfdef im_surface(df): # 波动率曲面作图# df = plot_df()df = fitting(df) #df.to_excel('iv_fitting.xlsx')df = df.set_index('k')y = np.array(df.index)x = np.array(df.columns)fig = go.Figure(data=[go.Surface(z=df.values, x=x, y=y)])fig.update_layout(scene = dict(xaxis_title='剩余期限',yaxis_title='执行价格',zaxis_title='期权价格'),width=1400,margin=dict(r=200, b=100, l=100, t=100))#fig.write_image("fig1.jpg")plotly.offline.plot(fig)def main():date = '0630'# plot_df()df = extra_data(date) # 提取数据df = data_clear(df) # 数据清洗df = cal_mc(df) # 计算期权定价df = data_pivot(df) # 数据透视表im_surface(df) # 波动率曲面if __name__ == '__main__':main()
【Python+Tushare·期权】微笑波动率曲线成因分析 并利用蒙特卡洛模拟估计欧式看涨期权的价值
