用 pySABR 擬合波動微笑——SABR 模型的 Python 實現
為了模擬一些波動性微笑,我使用了 python 的pySABR包。
我遇到了這樣一種情況,即我有兩個幾乎相同的程式碼,用於兩個不同的波動率微笑,缺少 ATM 報價,而 pySABR 在一種情況下可以正確適應 ATM 波動性,而在另一種情況下則不能。
當一切正常時就是這種情況:
import pysabr import numpy as np from pysabr import Hagan2002LognormalSABR as sabr from pysabr import hagan_2002_lognormal_sabr as hagan2002 testStrikes = np.array([0.04, 0.06, 0.08, 0.10]) testVols = np.array([23.52, 16.24, 20.17, 26.19]) forward_3m_6m = (1/0.25) * (-1 + (1+0.0753*0.5) / (1+0.0747*0.25)) calibration = sabr(f = forward_3m_6m, shift = 0, t = 0.5, beta = 0.5).fit(testStrikes, testVols) smile = [] test_smile = [] for strike in testStrikes: smile.append(sabr(f = forward_3m_6m, shift = 0, t = 0.5, v_atm_n = 136.75/10000.00, beta = 0.5, rho = calibration[1], volvol = calibration[2]).lognormal_vol(strike) * 100.00) test_smile.append(hagan2002.lognormal_vol(strike, forward_3m_6m, 0.5, calibration[0], beta, calibration[1], calibration[2]) * 100.00) print(smile) print(test_smile) print(hagan2002.lognormal_vol(k = 0.0745136, f = forward_3m_6m, t = 0.5, alpha = calibration[0], beta = 0.5, rho = calibration[1], volvol = calibration[2]) * 100.00)
輸出是
[23.52579980722943, 16.22619971530687, 20.186954023608315, 26.176954813043512] [23.52656681608356, 16.227190950406076, 20.188104613955648, 26.178058303454062] 18.369296025878036
第一個和第二個波動率列表之間的區別在於,第一個是使用 ATM 正常波動率報價而不是參數 alpha 建構的。
這是第二種情況的程式碼,當我遇到不同方法的不同結果以及由此產生的 ATM 對數正態波動率值不足的問題時:
import pysabr import numpy as np from pysabr import Hagan2002LognormalSABR as LNsabr from pysabr import hagan_2002_lognormal_sabr as hagan2002LN strikes = np.array([0.05, 0.055, 0.06, 0.0650, 0.07, 0.08, 0.0850, 0.09, 0.095, 0.10]) LogNormalVols = np.array([18.90, 17.30, 16.34, 16.29, 17.19, 20.29, 21.89, 23.42, 24.84, 26.16]) vol_n_ATM = 141.01 / 10000.00 forward_3m_6m = (1 / 0.25) * (- 1 + ( 1+ 0.0756 * 0.5) / (1 + 0.0746 * 0.25)) discount_0m_6m = 0.0753 beta = 0.5 calibration_LN = LNsabr(forward_3m_6m, 0, 0.5, beta).fit(strikes, LogNormalVols) modelVols_LN = [] test_LN = [] for strike in strikes: modelVols_LN.append(LNsabr(forward_3m_6m, 0, 0.5, vol_n_ATM, beta, calibration_LN[1], calibration_LN[2]).lognormal_vol(strike) * 100.00) test_LN.append(hagan2002LN.lognormal_vol(strike, f = forward_3m_6m, t = 0.5, alpha = calibration_LN[0], beta = beta, rho = calibration_LN[1], volvol = calibration_LN[2]) * 100.00) print(modelVols_LN) print(test_LN) print(hagan2002LN.lognormal_vol(k = 0.0753, f = forward_3m_6m, t = 0.5, alpha = calibration_LN[0], beta = beta, rho = calibration_LN[1], volvol = calibration_LN[2]) * 100.00)
輸出是
[20.14451199912703, 18.54257849585578, 17.499371075753768, 17.1951750118971, 17.677189811525658, 20.011518064279755, 21.365538860352675, 22.691679608999, 23.95616514380161, 25.148945356594965] [68.433853990843, 67.98546902874503, 67.8392636042873, 67.91509150013229, 68.15026017672005, 68.91958756908917, 69.39185669606091, 69.89474554951227, 70.41454229587038, 70.94145252003263] 68.64132679555016
人們可以很容易地看到,使用所有四個校準參數來重現波動率微笑會導致對數正態波動率比預期的要高。
我究竟做錯了什麼?第二個程式碼中的錯誤在哪裡?任何幫助將不勝感激。
編輯:您需要在第二個範例中指定關鍵字參數
第二個範例的校準不佳是因為您沒有在
.fit()
函式之前在 LNsabr 對像中定義關鍵字參數。而不是寫:
calibration_LN = LNsabr(forward_3m_6m, 0, 0.5, beta).fit(strikes, LogNormalVols)
您只需要對定義關鍵字參數進行細微的更改:
calibration_LN = LNsabr(f = forward_3m_6m, shift = 0, t = 0.5, beta = beta).fit(strikes, LogNormalVols)
這導致了一個很好的擬合:
我已經擴展了我目前的附錄,以包含上面提供的 Python 程式碼的稍微重寫的版本。你可以自己執行它,看看你是否得到相同的校準參數。我不完全清楚為什麼不指定關鍵字參數時校準會變差。我希望這可以解決您的問題。有空時請提供一些回饋。
我將在下面留下我的舊答案,以供將來的讀者查看。
舊答案:pysabr 包中的優化常式可能有錯誤
錯誤的校準可能來自 Pysabr 包中的優化常式達到局部最小值或類似值。
我使用相同 SABR 對數正態多項式展開的內置函式在 MATLAB 中重新實現了您的第二個範例(導致問題),該函式也在
pysabr
Python 的包中提供(我的程式碼在答案末尾的附錄中提供。)。**可以在此處找到 MATLAB 文件,並在此處記錄 SABR 模型的清晰**校準範例。在後者中,他們介紹了校準 SABR 模型的兩種方式:首先,校準 $ \alpha, \rho, \nu $ 直接最後校準 $ \rho $ 和 $ \nu $ 通過暗示 $ \alpha $ 來自 ATM 的波動性。在 MATLABs 文件中使用您指定的 vols 和strikes 進行第一個校準範例之後,我得到校準參數:$$ \left[\alpha, \rho, \nu\right] = \left[0.0504, : 0.6407, :0.8797\right], $$
在哪裡 $ \nu $ 是 vol 參數的 vol。如果您在**
calibration_LN
**列表對像中插入這些估計值,您將在校準時觀察到非常適合您的 SABR 對數正態模型 $ \alpha $ 自由地:由此,我相信
pysabr
包中的優化程序在校準微笑時可能會出現一些問題。如果我找到包裝產生如此不合適的原因,我將更新我的答案。就目前而言,這是我能提供的唯一幫助。附錄:
Matlab程式碼
%Matlabs inbuilt functions uses settlementdate and exercisedate to %calculate T internally. These are just chosen to get T=0.5: Settle = '01-Jan-2020'; ExerciseDate = '01-Jul-2020'; %Checking if settle and exercise date becomes T=0.5. T = yearfrac(Settle, ExerciseDate, 1); Beta1 = 0.5; %Following the matlab guide: MarketStrikes = [0.05, 0.055, 0.06, 0.0650, 0.07, 0.08, 0.0850, 0.09, 0.095, 0.10]'; MarketVolatilities = [18.90, 17.30, 16.34, 16.29, 17.19, 20.29, 21.89, 23.42, 24.84, 26.16]'/100; CurrentForwardValue = (1 / 0.25) * (- 1 + ( 1+ 0.0756 * 0.5) / (1 + 0.0746 * 0.25)); objFun = @(X) MarketVolatilities - ... blackvolbysabr(X(1), Beta1, X(2), X(3), Settle, ... ExerciseDate, CurrentForwardValue, MarketStrikes); options_SABR = optimoptions(@lsqnonlin, 'Display', 'iter', ... 'MaxFunctionEvaluations', 5000, 'MaxIterations', 10000); %Minimizer: X = lsqnonlin(objFun, [0.5 0 1], [0 -1 0], [Inf 1 Inf], options_SABR); Alpha1 = X(1); Rho1 = X(2); Nu1 = X(3); Sabr_vol = blackvolbysabr(Alpha1, Beta1, Rho1, Nu1, Settle, ExerciseDate, CurrentForwardValue, MarketStrikes); [Alpha1, Rho1, Nu1] plot(MarketStrikes/T, MarketVolatilities) hold on plot(MarketStrikes/T, Sabr_vol) %shifted to see the curvature legend("MarketVolatilities", "SabrVolatilieis")
Python程式碼
import pysabr import numpy as np from pysabr import Hagan2002LognormalSABR as LNsabr from pysabr import hagan_2002_lognormal_sabr as hagan2002LN import matplotlib.pyplot as plt strikes = np.array([0.05, 0.055, 0.06, 0.0650, 0.07, 0.08, 0.0850, 0.09, 0.095, 0.10]) LogNormalVols = np.array([18.90, 17.30, 16.34, 16.29, 17.19, 20.29, 21.89, 23.42, 24.84, 26.16]) vol_n_ATM = 141.01 / 10000.00 forward_3m_6m = (1 / 0.25) * (- 1 + ( 1+ 0.0756 * 0.5) / (1 + 0.0746 * 0.25)) discount_0m_6m = 0.0753 beta = 0.5 calibration_LN = LNsabr(f = forward_3m_6m, shift = 0, t = 0.5, beta = beta).fit(strikes, LogNormalVols) Matlab_LN = [0.0504, 0.6407, 0.8797] #MATLAB calibrated parameters! WORKS! modelVols_LN = [] test_LN = [] for strike in strikes: modelVols_LN.append(LNsabr(forward_3m_6m, 0, 0.5, vol_n_ATM, beta, calibration_LN[1], calibration_LN[2]).lognormal_vol(strike) * 100.00) test_LN.append(hagan2002LN.lognormal_vol(strike, f = forward_3m_6m, t = 0.5, alpha = calibration_LN[0], beta = beta, rho = calibration_LN[1], volvol = calibration_LN[2]) * 100.00) print(np.around(modelVols_LN,3)) print(np.around(test_LN,3)) #print(hagan2002LN.lognormal_vol(k = 0.0753, f = forward_3m_6m, t = 0.5, alpha = calibration_LN[0], beta = beta, rho = calibration_LN[1], volvol = calibration_LN[2]) * 100.00) plt.plot(strikes, modelVols_LN, "r--", strikes, test_LN) print("------------------------------------------------------------------------------------------") print("------------------------------------------------------------------------------------------") print("Matlab calibrated parameters: %s. \nPysabr calibrated parameters: %s" % (Matlab_LN, calibration_LN)) plt.plot(strikes, modelVols_LN, "r--", strikes, test_LN)