【50ETF期權】 2. 歷史波動率 · Python 量化交易教程

# 【50ETF期權】 2. 歷史波動率 > 來源：https://uqer.io/community/share/560493a4f9f06c597565ef03 在本文中，我們將通過量化實驗室提供的數據，計算上證50ETF的歷史波動率數據 ```py from CAL.PyCAL import * import numpy as np import pandas as pd import matplotlib.pyplot as plt from matplotlib import rc rc('mathtext', default='regular') import seaborn as sns sns.set_style('white') import math from scipy.stats import mstats ``` 50ETF收盤價 ```py # 華夏上證50ETF secID = '510050.XSHG' begin = Date(2015, 2, 9) end = Date.todaysDate() fields = ['tradeDate', 'closePrice'] etf = DataAPI.MktFunddGet(secID, beginDate=begin.toISO().replace('-', ''), endDate=end.toISO().replace('-', ''), field=fields) etf['tradeDate'] = pd.to_datetime(etf['tradeDate']) etf = etf.set_index('tradeDate') etf.tail(3) ``` | | closePrice | | --- | --- | | tradeDate | | | 2015-09-22 | 2.237 | | 2015-09-23 | 2.180 | | 2015-09-24 | 2.187 | ## 1. EWMA模型計算歷史波動率 EWMA(Exponentially Weighted Moving Average)指數加權移動平均計算歷史波動率： ![](https://box.kancloud.cn/2016-07-30_579cbdba3d1f7.jpg) 其中 ![](https://box.kancloud.cn/2016-07-30_579cbdba51703.jpg) 上式中的 `Si` 為 `i` 天的收盤價，`λ` 為介于0和1之間的常數。也就是說，在第 `n?1` 天估算的第 `n` 天的波動率估計值 `σn` 由第 `n?1` 天的波動率估計值 `σn?1` 和收盤價在最近一天的變化百分比 `un?1` 決定。計算周期為 `N` 天的波動率時， `λ` 可以取為： ![](https://box.kancloud.cn/2016-07-30_579cbdba648b1.jpg) ```py def getHistVolatilityEWMA(secID, beginDate, endDate): cal = Calendar('China.SSE') spotBeginDate = cal.advanceDate(beginDate,'-520B',BizDayConvention.Preceding) spotBeginDate = Date(2006, 1, 1) begin = spotBeginDate.toISO().replace('-', '') end = endDate.toISO().replace('-', '') fields = ['tradeDate', 'preClosePrice', 'closePrice', 'settlePrice', 'preSettlePrice'] security = DataAPI.MktFunddGet(secID, beginDate=begin, endDate=end, field=fields) security['dailyReturn'] = security['closePrice']/security['preClosePrice'] # 日回報率 security['u2'] = (np.log(security['dailyReturn']))**2 # u2為復利形式的日回報率平方 # security['u2'] = (security['dailyReturn'] - 1.0)**2 # u2為日價格變化百分比的平方 security['tradeDate'] = pd.to_datetime(security['tradeDate']) periods = {'hv1W': 5, 'hv2W': 10, 'hv1M': 21, 'hv2M': 41, 'hv3M': 62, 'hv4M': 83, 'hv5M': 104, 'hv6M': 124, 'hv9M': 186, 'hv1Y': 249, 'hv2Y': 497} # 利用pandas中的ewma模型計算波動率 for prd in periods.keys(): # 此處的span實際上就是上面計算波動率公式中lambda表達式中的N security[prd] = np.round(np.sqrt(pd.ewma(security['u2'], span=periods[prd], adjust=False)), 5)*math.sqrt(252.0) security = security[security.tradeDate >= beginDate.toISO()] security = security.set_index('tradeDate') return security ``` ```py secID = '510050.XSHG' start = Date(2015, 2, 9) end = Date.todaysDate() hist_HV = getHistVolatilityEWMA(secID, start, end) hist_HV.tail(2) ``` | | preClosePrice | closePrice | dailyReturn | u2 | hv2M | hv1W | hv1Y | hv3M | hv4M | hv5M | hv2Y | hv1M | hv2W | hv6M | hv9M | | --- | --- | | tradeDate | | | | | | | | | | | | | | | | | 2015-09-23 | 2.237 | 2.180 | 0.974519 | 0.000666 | 0.511318 | 0.304791 | 0.446550 | 0.523224 | 0.519890 | 0.511635 | 0.379718 | 0.449090 | 0.344477 | 0.502269 | 0.472743 | | 2015-09-24 | 2.180 | 2.187 | 1.003211 | 0.000010 | 0.499095 | 0.250658 | 0.444804 | 0.514969 | 0.513699 | 0.506714 | 0.378925 | 0.428453 | 0.312410 | 0.498142 | 0.470203 | ```py secID = '510050.XSHG' start = Date(2007, 1, 1) end = Date.todaysDate() hist_HV = getHistVolatilityEWMA(secID, start, end) ## ----- 50ETF歷史波動率 ----- fig = plt.figure(figsize=(10,12)) ax = fig.add_subplot(211) font.set_size(16) hist_plot = hist_HV[hist_HV.index >= Date(2015,2,9).toISO()] etf_plot = etf[etf.index >= Date(2015,2,9).toISO()] lns1 = ax.plot(hist_plot.index, hist_plot.hv1M, '-', label = u'HV(1M)') lns2 = ax.plot(hist_plot.index, hist_plot.hv2M, '-', label = u'HV(2M)') ax2 = ax.twinx() lns3 = ax2.plot(etf_plot.index, etf_plot.closePrice, '-r', label = '50ETF closePrice') lns = lns1+lns2+lns3 labs = [l.get_label() for l in lns] ax.legend(lns, labs, loc=2) ax.grid() ax.set_xlabel(u"tradeDate") ax.set_ylabel(r"Historical Volatility") ax2.set_ylabel(r"closePrice") ax.set_ylim(0, 0.9) ax2.set_ylim(1.5, 4) plt.title('50ETF Historical EWMA Volatility') ## ----------------------------------- ## ----- 50ETF歷史波動率統計數據 ----- # 注意：該統計數據基于07年以來將近九年的歷史波動率得出 ax3 = fig.add_subplot(212) font.set_size(16) hist_plot = hist_HV[[u'hv2W', u'hv1M', u'hv2M', u'hv3M', u'hv4M', u'hv5M', u'hv6M', u'hv9M', u'hv1Y']] # Calculate the quantiles column wise quantiles = mstats.mquantiles(hist_plot, prob=[0.0, 0.25, 0.5, 0.75, 1.0], axis=0) labels = ['Minimum', '1st quartile', 'Median', '3rd quartile', 'Maximum'] for i, q in enumerate(quantiles): ax3.plot(q, label=labels[i]) # 在統計圖中標出某一天的波動率 date = Date(2015,8,27) last_day_HV = hist_plot.ix[date.toDateTime()].T ax3.plot(last_day_HV.values, 'dr', label=date.toISO()) # 在統計圖中標出最近一天的波動率 last_day_HV = hist_plot.tail(1).T ax3.plot(last_day_HV.values, 'sb', label=Date.fromDateTime(last_day_HV.columns[0]).toISO()) ax3.set_ylabel(r"Volatility") plt.xticks((0,1,2,3,4,5,6,7,8),(0.5,1,2,3,4,5,6,9,12)) plt.xlabel('Periods(Months)') plt.legend() plt.grid() ``` ![](https://box.kancloud.cn/2016-07-30_579cbdba79424.png) 波動率圖中，上圖表示50ETF收盤價格和歷史波動率的走勢關系： + 顯然，短周期波動率對于近期的波動更敏感 + 收盤價的下跌往往伴隨著波動率的上升，兩者的負相關性質明顯波動率圖中，下圖表示50ETF歷史波動率的統計數據，圖中給出了四分位波動率錐： + 8月底時，各個周期歷史波動率均處于歷史高位 + 目前，短周期波動率已經有所回落 ## 2. Close to Close 模型計算歷史波動率 m 天周期的Close to Close波動率： ![](https://box.kancloud.cn/2016-07-30_579cbdbaa113a.jpg) 其中 ![](https://box.kancloud.cn/2016-07-30_579cbdbab4554.jpg) 也就是說，在第 `n?1` 天估算的第 `n` 天的波動率估計值 `σn` 由前面 `m `天的每日收盤價變化百分比 `ui` 的標準差決定。 ```py ## 計算一段時間標的的歷史波動率，返回值包括以下不同周期的波動率： # 一周，半月，一個月，兩個月，三個月，四個月，五個月，半年，九個月，一年，兩年 def getHistVolatilityC2C(secID, beginDate, endDate): cal = Calendar('China.SSE') spotBeginDate = cal.advanceDate(beginDate,'-520B',BizDayConvention.Preceding) spotBeginDate = Date(2006, 1, 1) begin = spotBeginDate.toISO().replace('-', '') end = endDate.toISO().replace('-', '') fields = ['tradeDate', 'preClosePrice', 'closePrice', 'settlePrice', 'preSettlePrice'] security = DataAPI.MktFunddGet(secID, beginDate=begin, endDate=end, field=fields) security['dailyReturn'] = security['closePrice']/security['preClosePrice'] # 日回報率 security['u'] = np.log(security['dailyReturn']) # u2為復利形式的日回報率 security['tradeDate'] = pd.to_datetime(security['tradeDate']) periods = {'hv1W': 5, 'hv2W': 10, 'hv1M': 21, 'hv2M': 41, 'hv3M': 62, 'hv4M': 83, 'hv5M': 104, 'hv6M': 124, 'hv9M': 186, 'hv1Y': 249, 'hv2Y': 497} # 利用方差模型計算波動率 for prd in periods.keys(): security[prd] = np.round(pd.rolling_std(security['u'], window=periods[prd]), 5)*math.sqrt(252.0) security = security[security.tradeDate >= beginDate.toISO()] security = security.set_index('tradeDate') return security ``` ```py secID = '510050.XSHG' start = Date(2007, 1, 1) end = Date.todaysDate() hist_HV = getHistVolatilityC2C(secID, start, end) ## ----- 50ETF歷史波動率 ----- fig = plt.figure(figsize=(10,12)) ax = fig.add_subplot(211) font.set_size(16) hist_plot = hist_HV[hist_HV.index >= Date(2015,2,9).toISO()] etf_plot = etf[etf.index >= Date(2015,2,9).toISO()] lns1 = ax.plot(hist_plot.index, hist_plot.hv1M, '-', label = u'HV(1M)') lns2 = ax.plot(hist_plot.index, hist_plot.hv2M, '-', label = u'HV(2M)') ax2 = ax.twinx() lns3 = ax2.plot(etf_plot.index, etf_plot.closePrice, '-r', label = '50ETF closePrice') lns = lns1+lns2+lns3 labs = [l.get_label() for l in lns] ax.legend(lns, labs, loc=2) ax.grid() ax.set_xlabel(u"tradeDate") ax.set_ylabel(r"Historical Volatility") ax2.set_ylabel(r"closePrice") ax.set_ylim(0, 0.9) ax2.set_ylim(1.5, 4) plt.title('50ETF Historical Close-to-Close Volatility') ## ----------------------------------- ## ----- 50ETF歷史波動率統計數據 ----- # 注意：該統計數據基于07年以來將近九年的歷史波動率得出 ax3 = fig.add_subplot(212) font.set_size(16) hist_plot = hist_HV[[u'hv2W', u'hv1M', u'hv2M', u'hv3M', u'hv4M', u'hv5M', u'hv6M', u'hv9M', u'hv1Y']] # Calculate the quantiles column wise quantiles = mstats.mquantiles(hist_plot, prob=[0.0, 0.25, 0.5, 0.75, 1.0], axis=0) labels = ['Minimum', '1st quartile', 'Median', '3rd quartile', 'Maximum'] for i, q in enumerate(quantiles): ax3.plot(q, label=labels[i]) # 在統計圖中標出某一天的波動率 date = Date(2015,8,27) last_day_HV = hist_plot.ix[date.toDateTime()].T ax3.plot(last_day_HV.values, 'dr', label=date.toISO()) # 在統計圖中標出最近一天的波動率 last_day_HV = hist_plot.tail(1).T ax3.plot(last_day_HV.values, 'sb', label=Date.fromDateTime(last_day_HV.columns[0]).toISO()) ax3.set_ylabel(r"Volatility") plt.xticks((0,1,2,3,4,5,6,7,8),(0.5,1,2,3,4,5,6,9,12)) plt.xlabel('Periods(Months)') plt.legend() plt.grid() ``` ![](https://box.kancloud.cn/2016-07-30_579cbdbac620e.png) 波動率圖中，上圖表示50ETF收盤價格和歷史波動率的走勢關系： + 顯然，短周期波動率對于近期的波動更敏感 + 收盤價的下跌往往伴隨著波動率的上升，兩者的負相關性質明顯波動率圖中，下圖表示50ETF歷史波動率的統計數據，圖中給出了四分位波動率錐： + 8月底時，各個周期歷史波動率均處于歷史高位 + 目前，短周期波動率已經有所回落明顯地，相對于EWMA計算的歷史波動率，Close to Close波動率對于最近價格波動反應比較遲鈍