4
$\begingroup$

I am actively trying to price an option using bates model on Quantlib.However,when I input my volatility I find the same Black Prices with the basic Heston Model.I wanted to know if my code was right.

import QuantLib as ql
import numpy as np 
import scipy
from scipy.optimize import *

def setup_helpers(engine, expiration_dates, strikes, 
                  data, ref_date, spot, yield_ts, 
                  dividend_ts):
    heston_helpers = []
    grid_data = []
    for i, date in enumerate(expiration_dates):
        for j, s in enumerate(strikes):
            t = (date - ref_date )
            p = ql.Period(t, ql.Days)
            vols = data[i][j]
            helper = ql.HestonModelHelper(
                p, calendar, spot, s, 
                ql.QuoteHandle(ql.SimpleQuote(vols)),
                yield_ts, dividend_ts,ql.BlackCalibrationHelper.ImpliedVolError)
            helper.setPricingEngine(engine)
            heston_helpers.append(helper)
            grid_data.append((date, s))
    return heston_helpers, grid_data

def cost_function_generator1(model, helpers,norm=False):
    def cost_function(params):
        pen = 1
        params_ = ql.Array(list(params))
        model.setParams(params_)
        pen_c = 2 * params[0] * params[1] -params[2]**2
        if pen_c < 0:
            pen = 1e50
        #err =np.sum(spreads*(np.array([(opt.modelValue() - opt.marketValue())**2  for i, opt in enumerate(helpers)])))/np.sum(spreads)#[h.calibrationError() for h in helpers]r h in helpers]
        err =  pen * np.sum(spreads * np.abs([h.calibrationError() for h in helpers])**2)/len(helpers)
        if norm:
            return pen *err #pen * np.sum(spreads *np.abs([h.calibrationError() for h in heston_helpers4])**2/np.array(data).reshape(1,195)[0])/np.sum(spreads)
        else:
            return error
    return cost_function


def cost_function_generator2(model, helpers,norm=False):
    def cost_function(params):
        params_ = ql.Array(list(params))
        model.setParams(params_)
        #err = np.sum(spreads*(np.array([(opt.modelValue() - opt.marketValue())**2  for i, opt in enumerate(helpers)])))/np.sum(spreads)#[h.calibrationError() for h in helpers]
        err = np.sum(spreads * np.abs([h.calibrationError() for h in helpers]))**2/len(helpers)
        if norm:
            return err#len(heston_helpers4)
        else:
            return error
    return cost_function

def calibration_report(helpers, grid_data, detailed=False):
    IVMSE = np.sum(spreads* np.abs([h.calibrationError() for h in helpers])**2)/len(helpers)
    if detailed: print("-"*100)
    summary = "IVMSE  : %1.9f" % (IVMSE)
    print(summary)
    return IVMSE
    
def setup_model(_yield_ts, _dividend_ts, _spot, 
                init_condition=(0.02,0.2,0.5,0.1,0.01,0.1,0.1,0.1)):
    theta, kappa, sigma, rho, v0,lambs,mu,sig_v = init_condition
    process = ql.BatesProcess (_yield_ts, _dividend_ts, 
                           ql.QuoteHandle(ql.SimpleQuote(_spot)), 
                           v0, kappa, theta, sigma, rho,lambs,mu,sig_v)
    model = ql.BatesModel(process)
    engine = ql.BatesEngine(model) 
    return model, engine

def Constraint(x):
    theta, kappa, sigma, rho, v0,lambda_,mu,sig_v = x
    A = 2*kappa*theta - sigma**2 
    return A

const = {'type': 'ineq', 'fun': Cst}
summary= []

This is the initial Setup.

This is the initialising step

model4, engine4 = setup_model(yield_ts, dividend_ts, spot)
heston_helpers4, grid_data4 = setup_helpers(engine4, expiration_dates, strikes, data,calculation_date, spot, yield_ts, dividend_ts)
initial_condition = (0.02, 0.2, 0.5, 0.1, 0.01,0.01,0.01,0.01)
bounds = [(1e-5,10),(1e-5,100), (1e-5,10), (-.99,.99), (1e-5,10),(0,1),(-1,1),(0,1)]
bnds = ((1e-5,10),(1e-5,100), (1e-5,10), (-.99,.99), (1e-5,10),(0,1),(-1,1),(0,1))

And lastly,this is the calibration(it takes 3 mins with 210 data points) :

%%time
cost_function = cost_function_generator1(
    model4, heston_helpers4, norm=True)
cost_function2 = cost_function_generator2(
    model4, heston_helpers4, norm=True)
sol = differential_evolution(cost_function, bounds=bounds,maxiter = 50 , tol = 1e-9,popsize=25)
params_ = ql.Array(list(sol.x))
model4.setParams(params_)
sol = minimize(cost_function2,sol.x, bounds=bnds,tol =1e-9,constraints = const)
params_ = ql.Array(list(sol.x))
model4.setParams(params_)
theta_b, kappa_b, sigma_b, rho_b, v0_b , lambda_ ,mu,sigma_v = model4.params()
error = calibration_report(heston_helpers4, grid_data4)
summary.append(["Scipy DE1", error] + list(model4.params()))

Moreover,Should I use the ql.HestonModelHelper in my case ?

Thank you

$\endgroup$

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.