This setion of the documentation illustrates how to calculate implied volatilities and how to calibrate a model to VSTOXX volatility index call option quotes. The example implements the calibration for a total of one month worth of data.
In [1]:
import numpy as np
import pandas as pd
We start by loading VSTOXX data from a pandas HDFStore into DataFrame objects (source: Eurex, cf. http://www.eurexchange.com/advanced-services/).
In [2]:
h5 = pd.HDFStore('./data/vstoxx_march_2014.h5', 'r')
vstoxx_index = h5['vstoxx_index']
vstoxx_futures = h5['vstoxx_futures']
vstoxx_options = h5['vstoxx_options']
h5.close()
VSTOXX index for the first quarter of 2014.
In [3]:
%matplotlib inline
vstoxx_index['V2TX'].plot()
Out[3]:
The VSTOXX futures data (8 futures maturities/quotes per day).
In [4]:
vstoxx_futures.info()
In [5]:
vstoxx_futures.tail()
Out[5]:
The VSTOXX options data. This data set is quite large due to the large number of European put and call options on the VSTOXX.
In [6]:
vstoxx_options.info()
In [7]:
vstoxx_options.tail()
Out[7]:
As a helper function we need a function to calculate all relevant third Fridays for all relevant maturity months of the data sets.
In [8]:
import datetime as dt
import calendar
def third_friday(date):
day = 21 - (calendar.weekday(date.year, date.month, 1) + 2) % 7
return dt.datetime(date.year, date.month, day)
In [9]:
third_fridays = {}
for month in set(vstoxx_futures['EXP_MONTH']):
third_fridays[month] = third_friday(dt.datetime(2014, month, 1))
In [10]:
third_fridays
Out[10]:
Often calibration efforts are undertaken to replicate the market implied volatilities or the so-called volatility surface as good as possible. With DX Analytics and the BSM_european_option class, you can efficiently calculate (i.e. numerically estimate) implied volatilities. For the example, we use the VSTOXX futures and call options data from 31. March 2014.
In [11]:
from dx import *
Some definitions, the pre-selection of option data and the pre-definition of the market environment needed.
In [12]:
V0 = 17.6639 # VSTOXX level on 31.03.2014
futures_data = vstoxx_futures[vstoxx_futures.DATE == '2014/3/31'].copy()
options_data = vstoxx_options[(vstoxx_options.DATE == '2014/3/31')
& (vstoxx_options.TYPE == 'C')].copy()
me = market_environment('me', dt.datetime(2014, 3, 31))
me.add_constant('initial_value', 17.6639) # index on 31.03.2014
me.add_constant('volatility', 2.0) # for initialization
me.add_curve('discount_curve', constant_short_rate('r', 0.01)) # assumption
options_data['IMP_VOL'] = 0.0 # initialization new iv column
The following loop now calculates the implied volatilities for all those options whose strike lies within the defined tolerance level.
In [13]:
%%time
tol = 0.3 # tolerance level for moneyness
for option in options_data.index:
# iterating over all option quotes
forward = futures_data[futures_data['MATURITY'] == \
options_data.loc[option]['MATURITY']]['PRICE'].values
# picking the right futures value
if (forward * (1 - tol) < options_data.loc[option]['STRIKE']
< forward * (1 + tol)):
# only for options with moneyness within tolerance
call = options_data.loc[option]
me.add_constant('strike', call['STRIKE'])
me.add_constant('maturity', call['MATURITY'])
call_option = BSM_european_option('call', me)
options_data.loc[option, 'IMP_VOL'] = \
call_option.imp_vol(call['PRICE'], 'call', volatility_est=0.6)
A selection of the results.
In [14]:
options_data[60:70]
Out[14]:
And the complete results visualized.
In [15]:
import matplotlib.pyplot as plt
%matplotlib inline
plot_data = options_data[options_data.IMP_VOL > 0]
plt.figure(figsize=(8, 6))
for maturity in sorted(set(options_data['MATURITY'])):
data = plot_data.isin({'MATURITY': [maturity,]})
data = plot_data[plot_data.MATURITY == maturity]
# select data for this maturity
plt.plot(data['STRIKE'], data['IMP_VOL'],
label=maturity.date(), lw=1.5)
plt.plot(data['STRIKE'], data['IMP_VOL'], 'r.')
plt.grid()
plt.xlabel('strike')
plt.ylabel('implied volatility of volatility')
plt.legend()
plt.show()
This sub-section now implements the model calibration based on selected options data. In particular, we choose, for a given pricing date, the following options data:
The following following returns the relevant market data per calibration date:
In [16]:
tol = 0.2
def get_option_selection(pricing_date, maturity, tol=tol):
''' Function selects relevant options data. '''
forward = vstoxx_futures[(vstoxx_futures.DATE == pricing_date)
& (vstoxx_futures.MATURITY == maturity)]['PRICE'].values[0]
option_selection = \
vstoxx_options[(vstoxx_options.DATE == pricing_date)
& (vstoxx_options.MATURITY == maturity)
& (vstoxx_options.TYPE == 'C')
& (vstoxx_options.STRIKE > (1 - tol) * forward)
& (vstoxx_options.STRIKE < (1 + tol) * forward)]
return option_selection, forward
Given the options and their respective quotes to which to calibrate the model, the function get_option_models returns the DX Analytics option models for all relevant options. As risk factor model we use the square_root_diffusion class.
In [17]:
def get_option_models(pricing_date, maturity, option_selection):
''' Models and returns traded options for given option_selection object. '''
me_vstoxx = market_environment('me_vstoxx', pricing_date)
initial_value = vstoxx_index['V2TX'][pricing_date]
me_vstoxx.add_constant('initial_value', initial_value)
me_vstoxx.add_constant('final_date', maturity)
me_vstoxx.add_constant('currency', 'EUR')
me_vstoxx.add_constant('frequency', 'W')
me_vstoxx.add_constant('paths', 10000)
csr = constant_short_rate('csr', 0.01)
# somewhat arbitrarily chosen here
me_vstoxx.add_curve('discount_curve', csr)
# parameters to be calibrated later
me_vstoxx.add_constant('kappa', 1.0)
me_vstoxx.add_constant('theta', 1.2 * initial_value)
me_vstoxx.add_constant('volatility', 1.0)
vstoxx_model = square_root_diffusion('vstoxx_model', me_vstoxx)
# square-root diffusion for volatility modeling
# mean-reverting, positive process
# option parameters and payoff
me_vstoxx.add_constant('maturity', maturity)
payoff_func = 'np.maximum(maturity_value - strike, 0)'
option_models = {}
for option in option_selection.index:
strike = option_selection['STRIKE'].ix[option]
me_vstoxx.add_constant('strike', strike)
option_models[option] = \
valuation_mcs_european_single(
'eur_call_%d' % strike,
vstoxx_model,
me_vstoxx,
payoff_func)
return vstoxx_model, option_models
The function calculate_model_values estimates and returns model value estimates for all relevant options given a parameter set for the square_root_diffusion risk factor model.
In [18]:
def calculate_model_values(p0):
''' Returns all relevant option values.
Parameters
===========
p0 : tuple/list
tuple of kappa, theta, volatility
Returns
=======
model_values : dict
dictionary with model values
'''
kappa, theta, volatility = p0
vstoxx_model.update(kappa=kappa,
theta=theta,
volatility=volatility)
model_values = {}
for option in option_models:
model_values[option] = \
option_models[option].present_value(fixed_seed=True)
return model_values
The calibration of the pricing model is based on the minimization of the mean-squared error (MSE) of the model values vs. the market quotes. The MSE calculation is implemented by the function mean_squared_error which also penalizes economically implausible parameter values.
In [19]:
i = 0
def mean_squared_error(p0):
''' Returns the mean-squared error given
the model and market values.
Parameters
===========
p0 : tuple/list
tuple of kappa, theta, volatility
Returns
=======
MSE : float
mean-squared error
'''
if p0[0] < 0 or p0[1] < 5. or p0[2] < 0 or p0[2] > 10.:
return 100
global i, option_selection, vstoxx_model, option_models, first, last
pd = dt.datetime.strftime(
option_selection['DATE'].iloc[0].to_pydatetime(),
'%d-%b-%Y')
mat = dt.datetime.strftime(
option_selection['MATURITY'].iloc[0].to_pydatetime(),
'%d-%b-%Y')
model_values = calculate_model_values(p0)
option_diffs = {}
for option in model_values:
option_diffs[option] = (model_values[option]
- option_selection['PRICE'].loc[option])
MSE = np.sum(np.array(option_diffs.values()) ** 2) / len(option_diffs)
if i % 150 == 0:
# output every 0th and 100th iteration
if i == 0:
print '%12s %13s %4s %6s %6s %6s --> %6s' % \
('pricing_date', 'maturity_date', 'i', 'kappa',
'theta', 'vola', 'MSE')
print '%12s %13s %4d %6.3f %6.3f %6.3f --> %6.3f' % \
(pd, mat, i, p0[0], p0[1], p0[2], MSE)
i += 1
return MSE
The function get_parameter_series calibrates the model to the market data for every date contained in the pricing_date_list object for all maturities contained in the maturity_list object.
In [20]:
import scipy.optimize as spo
def get_parameter_series(pricing_date_list, maturity_list):
global i, option_selection, vstoxx_model, option_models, first, last
# collects optimization results for later use (eg. visualization)
parameters = pd.DataFrame()
for maturity in maturity_list:
first = True
for pricing_date in pricing_date_list:
option_selection, forward = \
get_option_selection(pricing_date, maturity)
vstoxx_model, option_models = \
get_option_models(pricing_date, maturity, option_selection)
if first is True:
# use brute force for the first run
i = 0
opt = spo.brute(mean_squared_error,
((0.5, 2.51, 1.), # range for kappa
(10., 20.1, 5.), # range for theta
(0.5, 10.51, 5.0)), # range for volatility
finish=None)
i = 0
opt = spo.fmin(mean_squared_error, opt,
maxiter=200, maxfun=350, xtol=0.0000001, ftol=0.0000001)
parameters = parameters.append(
pd.DataFrame(
{'pricing_date' : pricing_date,
'maturity' : maturity,
'initial_value' : vstoxx_model.initial_value,
'kappa' : opt[0],
'theta' : opt[1],
'sigma' : opt[2],
'MSE' : mean_squared_error(opt)}, index=[0,]),
ignore_index=True)
first = False
last = opt
return parameters
This completes the set of necessary function to implement such a larger calibration effort. The following code defines the dates for which a calibration shall be conducted and for which maturities the calibration is required.
In [21]:
%%time
pricing_date_list = pd.date_range('2014/3/1', '2014/3/31', freq='B')
maturity_list = [third_fridays[7]]
parameters = get_parameter_series(pricing_date_list, maturity_list)
The results are now stored in the pandas DataFrame called parameters. We set a new index and inspect the last results. Throughout the MSE is pretty low indicated a good fit of the model to the market quotes.
In [22]:
paramet = parameters.set_index('pricing_date')
paramet.tail()
Out[22]:
This is also illustrated by the visualization of the time series data for the calibrated/optimal parameter values. The MSE is below 0.01 throughout.
In [23]:
%matplotlib inline
paramet[['kappa', 'theta', 'sigma', 'MSE']].plot(subplots=True, color='b', figsize=(10, 12))
plt.tight_layout()
The following generates a plot of the calibration results for the last calibration day. The absolute price differences are below 0.10 EUR for all options.
In [24]:
index = paramet.index[-1]
opt = np.array(paramet[['kappa', 'theta', 'sigma']].loc[index])
option_selection = get_option_selection(index, maturity_list[0], tol=tol)[0]
model_values = np.sort(np.array(calculate_model_values(opt).values()))[::-1]
import matplotlib.pyplot as plt
%matplotlib inline
fix, (ax1, ax2) = plt.subplots(2, sharex=True, figsize=(8, 8))
strikes = option_selection['STRIKE'].values
ax1.plot(strikes, option_selection['PRICE'], label='market quotes')
ax1.plot(strikes, model_values, 'ro', label='model values')
ax1.set_ylabel('option values')
ax1.grid(True)
ax1.legend(loc=0)
wi = 0.25
ax2.bar(strikes - wi / 2., model_values - option_selection['PRICE'],
label='market quotes', width=wi)
ax2.grid(True)
ax2.set_ylabel('differences')
ax2.set_xlabel('strikes')
Out[24]:
Copyright, License & Disclaimer
© Dr. Yves J. Hilpisch | The Python Quants GmbH
DX Analytics (the "dx library") is licensed under the GNU Affero General Public License version 3 or later (see http://www.gnu.org/licenses/).
DX Analytics comes with no representations or warranties, to the extent permitted by applicable law.
http://www.pythonquants.com | analytics@pythonquants.com | http://twitter.com/dyjh
Python Quant Platform | http://quant-platform.com
Derivatives Analytics with Python (Wiley Finance) | http://derivatives-analytics-with-python.com
Python for Finance (O'Reilly) | http://shop.oreilly.com/product/0636920032441.do