計算 Heston 模型的機率密度函式 (PDF)

我正在嘗試使用 Breeden Litzenberger 公式計算 Heston 模型的 PDF。

我已經使用 python 計算了罷工範圍的赫斯頓隱含波動率（我使用三次樣條插值插值）：

為了獲得 PDF，我正在使用 Breeden 和 Litzenberger 論壇：

對於執行價格的隱含波動率的導數，我使用了數值微分。與 Breeden 和 Litzenberger 論壇一起執行此操作，我得到了一個如下所示的 PDF：

這看起來不對。我使用了相同的方法，但 SABR 模型可以生成好看的 PDF。我很確定我的程式碼是正確的，所以我想知道 Heston 模型的 PDF 是否遺漏了什麼？

對於 Heston 模型，我使用了以下輸入（S=1 和打擊範圍 K $ \in $ (0.8,1.3):

我實際上也不確定 Black Scholes sigma 使用什麼參數，但它不應該製作或破壞上述 PDF？

我什至無法重新創建你的 vol 微笑……當我插入你提供的參數時（在 $ \tau = 0.12 $ ）我得到一個向下傾斜的微笑，在我看到的罷工中沒有最低限度

然後，我在每個近距離行權網格中取消期權價格併計算價格的曲率，它非常接近 rn pdf（只需要修正一個接近 1 的 dcf 因子在這麼短的時間內），它看起來大致如預期的那樣

我在下面附上了我的程式碼，您應該很容易使用參數並嘗試找出腳本中出了什麼問題（如果您共享程式碼，我們可能會提供更多幫助）

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()

首先，一些更正（我不確定，只是猜測）：

1. $ X_T $ - 是遠期標的物的行使價還是價格？

讓它罷工， $ X $ 基礎是 $ S_t $ 與前鋒： $ Fwd_T=S_t/D_T $ .

2. >

Breeden 和 Litzenberger 公式：

不，B&L 公式是這樣的： $ PDF(S_T)=D_T\cdot\frac{d^2C(X)}{dX^2} $ ， 在哪裡 $ D_T $ 是折扣係數。

3. 最後，我獲得 RND 的方法有兩個步驟：

一個。使用標準 BS 函式對歐洲呼叫定價 $ \sigma(X) $ ， IE $ Call(S_t,X) = D\cdot (Fwd \cdot N(d_1(\sigma(X))) - K \cdot N(d_2(\sigma(X)))) $ . 或使用看跌期權定價。

灣。假設微笑波動率並將 B&L 公式應用於 BS 看漲（或看跌）來計算 RND。

很可能，讓您的公式 Taylor (2007) 通過罷工兩次區分 BS 並一步完成。那麼我的食譜可以被認為是對該推導的交叉檢查（也許公式中有錯字？）。

是 Heston 模型還是 SABR 或其他什麼都無關緊要，因為您從參數化波動率形狀（微笑）開始。

引用自：https://quant.stackexchange.com/questions/57827