Maximum-Likelihood Calibration of a CIR Process for Heston-Type Rate Models

# Calibration of Heston version of CIR


# Calibration of Heston version of CIR












I'd like to calibrate a variant of Heston model for interest rates which is describe by this couple of SDE

\begin{aligned}dr_t&=a(b-r_t)+\sqrt{r_t}\sigma_t dW_t^1 \\ d\sigma_t&=k(\theta-\sigma_t)+\sqrt{\sigma_t}\zeta dW_t^2 \end{aligned}

with real markets data.

First of all I discretize the equations getting

\begin{aligned}r_t-r_{t-1}&=a(b-r_t)+\sqrt{r_{t-1}}\sigma_t \varepsilon_t \\ \sigma_t-\sigma_{t-1}&=k(\theta-\sigma_t)+\sqrt{\sigma_{t-1}}\zeta \epsilon_t \end{aligned}

And now I don't exactly how to proceed. I have some ideas:

- Try to apply a VAR (Vector AutoRegressive) model: rearraing the discretized equations I should get

\begin{aligned} \frac{r_t}{\sqrt{r_{t-1}}}&=ab+(a-1)\sqrt{r_{t-1}}+\sigma_t \varepsilon_t \\ \frac{\sigma_t}{\sqrt{\sigma_{t-1}}}&=k\theta+(k-1)\sqrt{\sigma_{t-1}}+\zeta \epsilon_t \end{aligned} or \begin{aligned} \frac{r_t-r_{t-1}}{\sqrt{r_{t-1}}}&=\frac{a(b-r_t)}{\sqrt{r_{t-1}}}+\sigma_t \varepsilon_t \\ \frac{\sigma_t-\sigma_{t-1}}{\sqrt{\sigma_{t-1}}}&=\frac{k(\theta-\sigma_t)}{\sqrt{\sigma_{t-1}}}+\zeta\epsilon_t \end{aligned} but I don't know whether this choice is suitable and which of these two proposals is the best one because of the square root and because the left sides of the equations are a trasformation of the right ones

- Apply OLS to the $\sigma_t$ equation, plug the results into the $r_t$ equation and the run again OLS or ML. The first part should not be a problem as the volatility is constant, $\zeta$, and OLS should be good rearraning in this way

$$\frac{\sigma_t-\sigma_{t-1}}{\sqrt{\sigma_t}}=\frac{k(\theta-\sigma_t)}{\sqrt{\sigma_t}}+\zeta\varepsilon_t $$ In the second equation ($r_t$ with $\sigma_t$ plugged in) I don't think OLS are a good choice because the model is heteroskedastic and $\varepsilon_t$ and $\epsilon_t$ may be correlated a priori.

Do you have any suggestion for the choice of the calibration technique? Moreover, how can I practical implement this calibration? I can use MatLab, Stata and Python (which would be my optimal choice if possible).

## Answer by AlexAbrahams (score 1)

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

I have attached below a Python implementation of the CIR calibration step by MLL. After this has been run you can continue to estimate $\mu$ (drift) and $\rho$ (correlation between Brownian motions).

```
# CIR Calibration Script
import numpy as np
from numpy import abs, r_, var, ones, sqrt, exp,log, real, mean
from numpy.linalg import solve
from scipy.optimize import fmin
from scipy.special import ive

# Parameters
t_ = len(data) # Length time series
p = (1 / t_) * ones((1, t_))  # Probabilities
t_obs = 252 * [...]  # Years
delta_t = [...]  # Time-step

def CIR_LOGLIK(parameters, dt, data, p):
    kappa = parameters[0]
    s2 = parameters[1]
    eta = parameters[2]
    c = (2 * kappa) / ((eta ** 2) * (1 - exp(-kappa * dt)))
    q = ((2 * kappa * s2) / (eta ** 2)) - 1
    u = c * exp(-kappa * dt) * data[:-1]
    v = c * data[1:]

    maxll = np.sum(- p[1:] * log(c) + p[1:] * (u + v) - p[1:] * log(v / u) * q / 2 - p[1:] * log(
        ive(q, 2 * sqrt(u * v)))
                 - p[1:] * abs(real(2 * sqrt(u * v))))
    return maxll

def FitCIR(data, dt, parameters=None, p=None):
    if parameters is None:
        x = r_[ones((1, t_ - 1)), data[np.newaxis, :t_ - 1]]
        ols = solve(x @ x.T, x @ data[1:t_].T)
        m = mean(data)
        v = var(data, ddof=1)
        parameters = r_[-log(ols[1]) / dt, m, sqrt(2 * ols[1] * v / m)]

    CIRparams, fevals = fmin(CIR_LOGLIK, parameters, args=(dt, data, p[0]), maxfun=1000, xtol=1e-8, ftol=1e-8,
                                disp=False, retall=True)
    return CIRparams

par_CIR = FitCIR(y[0, -t_obs:], delta_t, None, p)
kappa = par_CIR[0]
y_ = par_CIR[1]
eta = par_CIR[2]
```