SABR Smile Calibration: Use Explicit Keywords to Avoid Parameter Errors

# Fitting a volatility smile with pySABR -- Python implementation of SABR model


# Fitting a volatility smile with pySABR -- Python implementation of SABR model












In order to model some volatility smiles I'm using the python's pySABR package.

I ran into a situation when I have two almost identical pieces of code for two different volatility smiles missing the ATM quotes and the pySABR can properly fit the ATM volatility in one case and can't in another.

This is the case when everything is working just fine:

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

The output is

```
[23.52579980722943, 16.22619971530687, 20.186954023608315, 26.176954813043512]

[23.52656681608356, 16.227190950406076, 20.188104613955648, 26.178058303454062]

18.369296025878036
```

The difference between the first and the second volatility lists is that the first was built using the ATM normal volatility quote instead of the parameter alpha.

This is the code for the second case when I ran into a problem with different results for different methods and inadequate value of the resulting ATM lognormal volatility:

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

The output is

```
[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
```

One can easily see that reproducing the volatility smile using all four of the calibrated parameters leads to a way higher lognormal volatilities than expected.

What am I doing wrong? Where is the error in the second code? Any help will be appreciated.

## Answer by Pleb (score 5, accepted)

https://quant.stackexchange.com/a/68195

### Edit: You need to specify the keyword arguments in your second example

The poor calibration of your second example comes from the fact, that you didn't define the keyword arguments in the LNsabr-object before the `.fit()` function.

Instead of writing:

```
calibration_LN = LNsabr(forward_3m_6m, 0, 0.5, beta).fit(strikes, LogNormalVols)
```

You only need to make the slight change of defining the keyword arguments:

```
calibration_LN = LNsabr(f = forward_3m_6m, shift = 0, t = 0.5, beta = beta).fit(strikes, LogNormalVols)
```

This results in a good fit:

I have extended my current appendix to contain a slightly rewritten version of your Python code provided above. You can run this yourself and see whether you get the same calibrated parameters. I am not completely apparent why the calibration becomes poor when not specifying the keyword arguments. I hope this resolves your issue. Please provide some feedback when you have time.

I will leave my old answer below for future readers to see.

### Old answer: There might be an error with the optimization routine in the pysabr package

The bad calibration might come from the optimization routine in the Pysabr package hitting a local minimum or similar.

I reimplemented your second example (which caused issues) in MATLAB using their inbuilt function of the same SABR lognormal polynomial expansion which is also provided in the `pysabr` package in Python (My code is provided in an appendix at the end of my answer.). The MATLAB documentation can be found here, and a clear-cut calibration example of the SABR model is documented here. In the latter, they walk through both ways of calibrating the SABR model: First, calibrating $\alpha, \rho, \nu$ directly and lastly calibrating $\rho$ and $\nu$ by implying $\alpha$ from the ATM volatility. Following the first calibration example in MATLABs documentation with your specified vols and strikes, I get the calibrated parameters to:

$$\left[\alpha, \rho, \nu\right] = \left[0.0504, \: 0.6407, \:0.8797\right],$$

where $\nu$ is the vol of vol parameter. If you insert these estimates in the `calibration_LN` list object, you will observe a very good fit for your SABR log-normal model, when calibrating $\alpha$ freely:

From this, I believe the optimization routine in the `pysabr` package might have some issues when calibrating the smile. I will update my answer if I find the reason why the package produces such a bad fit. As for now, this is the only help I can provide.

### Appendix:

#### Matlab code

```
%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 code

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