0
$\begingroup$

I am trying to compute the PDF for the Heston model using the Breeden Litzenberger formula.

I have calculated the the Heston implied volatilities for a strike range (which i have interpolated using cubic spline interpolation) using python:

enter image description here

In order to get the PDF, I am using the Breeden and Litzenberger forumla: enter image description here

For the derivative of the implied volatility w.r.t the strike price, I have used numerical differentiation. Doing this together with the Breeden and Litzenberger forumla I get a PDF looking like this:

enter image description here

This does not look right. I have used the same approach but for the SABR model which results in a good looking PDF. I am pretty sure my code is correct so I was wondering if there is something I am missing about the PDF of the Heston model??

For the Heston model I have used the following inputs (S=1 and a strike range K$\in$(0.8,1.3):

enter image description here

I am actually also not sure what parameter to use for the Black Scholes sigma, but it should not make or break the above PDF?

Improve this question
$\endgroup$
9
  • $\begingroup$ What is the forward at expiry? Usually I would guess it was close to the bottom of the vol smile, but given your large negative rho it could be some way to the left hand side. In that case, you would need to extend your vol smiles down to lower strikes and you might see your pdf peak and then fall. $\endgroup$
    – StackG
    Sep 6, 2020 at 12:14
  • $\begingroup$ @StackG I have used the Spot S=1 such that I can think of the inputs K as the moneyness. I also just tried extrapolating such that the strike range was from 0.2-1.3 and changed the rho to -0.170152. This just made the left tail rise even higher! $\endgroup$
    – Klein
    Sep 6, 2020 at 12:23
  • 1
    $\begingroup$ If you have the right parameters (see @JohnDoe's comment), you can evaluate the characteristic function and invert it numerically to backout the density, without any partial derivatives $\endgroup$
    – Alex
    Sep 6, 2020 at 13:12
  • 1
    $\begingroup$ @Klein I think JohnDoe meant is that $\theta$ is the long-term of the variance. A reasonable long-term mean of volatility is something like 20%. This however means that $\theta=0.2^2=0.04$, i.e. much smaller than your value. Your number suggests that the long-term volatility is $\sqrt{\theta}\approx60\%$. That's a bit high perhaps? $\endgroup$
    – Alex
    Sep 6, 2020 at 13:33
  • 1
    $\begingroup$ @Klein I meant that the Heston closed form solution is $C=SP_1-Ke^{-rT}P_2$. This $P_2$ is the risk-neutral probability of $S_T$ being greater or equal than $K$ and an integral in terms of the characteristic function (Fourier transform of the RN density). You can obtain the RN density this way. $\endgroup$
    – Alex
    Sep 6, 2020 at 13:35

2 Answers 2

Reset to default
3
$\begingroup$

I can't quite even re-create your vol smile... when I plug in the parameters you've provided (at $\tau = 0.12$) I get a downward sloping vol smile that doesn't have a minimum at the strikes I looked at

Vols, Option Prices and PDF

I then backed out the options prices at each of a close-up grid of strikes and calculated the curvature of the prices, which is very close to the rn pdf (just need to correct by a factor of the dcf, which is close to 1 for such short times), and it looks roughly as expected

I've attached my code below, it should be very easy for you to play with the parameters and try to work out what is going wrong in your script (if you share the code, we might be able to help out more)

import QuantLib as ql
import numpy as np
import pandas as pd
from matplotlib import pyplot as plt

# Your parameters
tau = 0.1219; r = 0.0457; sigma = 0.4433; rho = -0.6175; nu = 0.474; theta = 0.3737; kappa = 1.042

today = ql.Date(1, 9, 2020)
expiry_date = today + ql.Period(int(365*tau), ql.Days)

# Setting up discount curve
risk_free_curve = ql.FlatForward(today, r, ql.Actual365Fixed())
flat_curve = ql.FlatForward(today, 0.0, ql.Actual365Fixed())
riskfree_ts = ql.YieldTermStructureHandle(risk_free_curve)
dividend_ts = ql.YieldTermStructureHandle(flat_curve)

# Setting up a Heston model
spot = 1

# I guess this is the correct mapping?
v0, sigma = nu, sigma

heston_process = ql.HestonProcess(riskfree_ts, dividend_ts, ql.QuoteHandle(ql.SimpleQuote(spot)), v0, kappa, theta, sigma, rho)
heston_model = ql.HestonModel(heston_process)
heston_handle = ql.HestonModelHandle(heston_model)
heston_vol_surface = ql.HestonBlackVolSurface(heston_handle)

# Now doing some pricing and curvature calculations
strikes = np.arange(0.5, 1.6, 0.01)
vols = [heston_vol_surface.blackVol(tau, x) for x in strikes]

option_prices = []

for strike in strikes:
    option = ql.EuropeanOption( ql.PlainVanillaPayoff(ql.Option.Call, strike), ql.EuropeanExercise(expiry_date))

    heston_engine = ql.AnalyticHestonEngine(heston_model)
    option.setPricingEngine(heston_engine)

    option_prices.append(option.NPV())

prices = pd.DataFrame([strikes, option_prices]).transpose()
prices.columns = ['strike', 'price']
prices['curvature'] = (-2 * prices['price'] + prices['price'].shift(1) + prices['price'].shift(-1)) / 0.01**2

# And plotting...
fig = plt.figure()
ax = fig.add_subplot(111)
ax2 = ax.twinx()

ax.plot(strikes, vols, label='Black Vols')
ax2.plot(strikes, option_prices, label='Option Prices', color='orange')
ax2.plot(prices['strike'], prices['curvature'], label='dC/dK (~pdf)', color='purple')

ax.legend(loc="lower left")
ax2.legend(loc="upper right")
ax.grid()
Improve this answer
$\endgroup$
2
  • 1
    $\begingroup$ Nice anser (+1): Would you agree that an alternative to the second partial derivative approach is to simply Fourier invert the characteristic function, which we know in closed-form for the Heston model? $\endgroup$
    – Alex
    Sep 6, 2020 at 13:43
  • $\begingroup$ I agree with the statement... all expect for the word "simply"! $\endgroup$
    – StackG
    Sep 6, 2020 at 13:49
0
$\begingroup$

First, couple of corrections (I am not sure, just guessing):

  1. $X_T$ - is it strike or price of forward underlying?

Let it be strike, $X$ and the underlying is $S_t$ with forward: $Fwd_T=S_t/D_T$.

  1. Breeden and Litzenberger formula:

No, B&L formula is this: $PDF(S_T)=D_T\cdot\frac{d^2C(X)}{dX^2}$, where $D_T$ is discount factor.

  1. Finally, my recipe to get RND has two steps:

a. Use a standard B-S function to price European Call with $\sigma(X)$, i.e $Call(S_t,X) = D\cdot (Fwd \cdot N(d_1(\sigma(X))) - K \cdot N(d_2(\sigma(X))))$. Or use Put pricing.

b. Assume smiled volatility and apply B&L formula to B-S Call (or Put) to calculate RND.

Very likely, to get your formula Taylor (2007) twice differentiated B-S by strike and made the same just in one step. Then my recipe can be considered as a cross-check of that derivation (maybe there is typo in formula?).

Is it Heston model or SABR or something else it does not matter much because you start with parameterized volatility shape (smile).

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.