Estimation of a CIR process with jumps using a closed form approximation likelihood under a strong approximation of order 1

Abstract

We propose here an approach in order to estimate parameters of the CIR model with jumps in the case where the distribution of jump amplitude is estimated non-parametrically. Since the knowledge of the exact distribution of the jump amplitude is a challenge, in this paper we choose not to fix this law in advance but to estimate it on the basis of the available observations. The method of estimation we propose here is based on the approximation of the closed form of transition density. Since the CIR does not have an explicit solution, it is approximated by the second order Milstein scheme in order to have a more accurate approximation. The method of estimation is then applied on real data, which are the Federal Funds rate and 3 Month T-Bill rate. These two sets of data are used to estimated parameters of the CIR model. We then compare our results to those obtained from Vasicek and Brennon–Swartz models with jumps. Results indicate that there is no clear winner of models competitions. Apparently depending on the nature and structural components of the data, there is a winner. The challenge here is that, there is a trade off between the sample size, the number of jumps and the efficiency of estimates. More data involves the likelihood to have more jumps and thereby less efficient are estimates.

Appendices

Appendices

Appendix A: Simulation of a sample from the distribution of the data

The aim of this section is to describe the method used to simulate a sample from the distribution of the data. It is stated as follows: Given that you only have a finite number n of data points, we can turn the empirical c.d.f. into a continuous R.V by using linear interpolation between the $X_{(i)}'s$ by:

$$\begin{aligned} F_X(x)=\left\{ \begin{array}{ll} 0 &{}\text { If } x< X_{(1)}\\ &{}\\ \dfrac{i-1}{n-1}+\dfrac{x-X_{(i)}}{(n-1)(X_{(i+1)}-X_{(i)})}&{}\text { If } X_{(i)} \le x<X_{(i+1)},\forall i\\ &{}\\ 1&{}\text { If } x\ge X_{(n)}\\ \end{array} \right. \end{aligned}$$

where $X_{(1)}\le \ldots \le X_{(n)}$ represents the order statistics.

Now to find the inverse of F, we proceed as follows: For a given x and y such that $F(x)=y$, it follows that if $X_{(i)} \le x<X_{(i+1)}$ then $F(X_{(i)}) \le y <F(X_{(i+1)})$ and its follows that

$$\begin{aligned} \dfrac{i-1}{n-1}\le y< \dfrac{i}{n-1}\text { and this implies that } i =[y(n-1)+1]. \end{aligned}$$

Since it is well known that if $U\sim \mathcal {U}([0,1])$ then $F_X^{-1}(U)$ has the same distribution with the r.v X. One random value is then simulated from the distribution of the data as follows:

• 1. We simulate u from uniform distribution U([0, 1])

• 2. We then take the integer part (I) of $p=u(n-1)+1$, that is $I=[p]$

• 3. The expected value is then obtained by $x =x_{(I)}+(p-I+1)(x_{(I+1)}-x_{(I)})$

Appendix B: OLS approach

The aim of this section is to find the estimates of $\kappa$, $\mu$ and $\sigma$ that serve as initial parameters of $k,\mu$ and $\sigma$ for our method. To illustrate the principle of the method, we consider the following equation

$$\begin{aligned} dr_t=\kappa (\mu -r_t)dt +\sigma r_t^\gamma dW_t,\text { with }r_0\text { which known } \end{aligned}$$

(35)

In this paper, we consider the where $\gamma =0$, $\gamma =1/2$ and $\gamma =1$ as initial values in the Vasiseck, CIR and CKLS model. Since the method is used just to get an idea of the initial values of parameters $k,\mu$ and $\sigma$, we consider her the Euler scheme, also called the Euler–Maruyama approximation of (35) which is written as

$$\begin{aligned} \hat{r}(t_{i+1})= \hat{r}(t_i)+\kappa (\mu - \hat{r}(t_i))\Delta t+\sigma \mid \hat{r}(t_i)\mid ^\gamma \sqrt{\Delta t}\Delta Z_i, \end{aligned}$$

(36)

where $\hat{r}_0=r_0$ is given and $Z_1,\;Z_2,\ldots ,Z_{n-1}$ are standard independent $\mathcal {N}(0,1)$ Guassian variables. We also have $i=1,\ldots ,n-1$, $\Delta Z_i=Z_{i+1}-Z_i$ and we observe that here $r(t_1 ),\ldots , r(t_n )$ are empirical data while, $\hat{r}(t_1 ),\ldots , \hat{r}(t_n )$ are simulated data. After some transformations, we can write that

$$\begin{aligned} \frac{r(t_{i+1})-r(t_i)}{\mid r(t_i)\mid ^\gamma }=\kappa \mu \frac{\Delta t}{ \mid r(t_i)\mid ^\gamma } -\kappa \frac{ r(t_i)}{ \mid r(t_i)\mid ^\gamma }\Delta t +\sigma \sqrt{\Delta }Z_i. \end{aligned}$$

(37)

By setting $y_i=\frac{r(t_{i+1})-r(t_i)}{\mid r(t_i)\mid ^\gamma }$, $\beta _1=\kappa \mu$, $\beta _2=-\kappa$, $z_{1i}=\frac{\Delta t}{ \mid r(t_i)\mid ^\gamma }$, $z_{2i}=\frac{ r(t_i)}{ \mid r(t_i)\mid ^\gamma }\Delta t$, (37) can be written as $y_i=\beta _1z_{1i}+\beta _2z_{2i}+\epsilon _i,$ which is equivalent to

$$\begin{aligned}Y=\beta Z+\epsilon ,\text { with }\end{aligned}$$

$$\begin{aligned} Y=\left[ \begin{array}{c} y_1\\ y_2\\ . \\ .\\ .\\ y_{n-1}\\ \end{array} \right] , Z=\left[ \begin{array}{lr} \frac{\Delta t}{ \mid r(t_1)\mid ^\gamma }&{}\frac{ r(t_1}{ \mid r(t_1)\mid ^\gamma }\Delta t\\ \frac{\Delta t}{ \mid r(t_2)\mid ^\gamma }&{}\frac{ r(t_2)}{ \mid r(t_2)\mid ^\gamma }\Delta t\\ . &{}. \\ .&{}.\\ .&{}.\\ \frac{\Delta t}{ \mid r(t_{n-1})\mid ^\gamma }&{}\frac{ r(t_{n-1})}{ \mid r(t_{n-1})\mid ^\gamma }\Delta t\\ \end{array} \right] , \beta =\left[ \begin{array}{c} \beta _1\\ \beta _2\\ \end{array} \right] , \epsilon =\sigma \sqrt{\Delta t}\left[ \begin{array}{c} \mathcal {N}_{1}(0,1)\\ \mathcal {N}_{2}(0,1)\\ . \\ .\\ .\\ \mathcal {N}_{n-1}(0,1)\\ \end{array} \right] , \end{aligned}$$

where $\mathcal {N}_{1}(0,1),\;\mathcal {N}_{2}(0,1),\ldots , \mathcal {N}_{n-1}(0,1)$ are independent $\mathcal {N}(0,1).$ The OLSE of $\beta$ is then

$$\begin{aligned} \hat{\beta }=\arg \min \limits _{\beta }\mid \mid Y-Z\beta \mid \mid , \end{aligned}$$

where $\mid \mid .\mid \mid$ denotes Euclidean distance. It is well known that

$$\begin{aligned} \hat{\beta }=\left[ \begin{array}{c} \hat{\beta _1}\\ \hat{\beta _2}\\ \end{array} \right] =(Z^TZ)^{-1}Z^TY, \end{aligned}$$

The expression of $\hat{\beta }$ allows to deduce that

$$\begin{aligned}\hat{\kappa }=-\hat{\beta _2},\text { }\hat{\mu }=\frac{\hat{\beta _1}}{\hat{\kappa }},\text { }\hat{\sigma }=\frac{1}{\sqrt{n\Delta t}}\mid \mid Y-Z\hat{\beta }\end{aligned}$$

Appendix C: Simulation process

In this section, we present the way that the discrete process given by Eq. (38) has been. Let’s recall here that, the discrete form of the process is given by

$$\begin{aligned} \widehat{X}_{t_{i+1}}= \widehat{X}_{t^-_{i+1}}+ \sum \limits _{N_{t_i}+1}^{N_{t_{i+1}}}J_{t_i} \end{aligned}$$

(38)

With respect to equation (38), the main difficulty in simulating $\widehat{X}_{t_{i}}$ lies in the presence of the compound Poison process $\sum \nolimits _{N_{t_i}+1}^{N_{t_{i+1}}}J_{t_i}$. In fact, in the absence of this term (which is known as the jump term), the discrete process is easily simulated using Eq. (6). To Simulate the full process including the jump term, we use a rejection sampling method which can be described as follows:

• Step 1 Simulate a $N_T$-sample $\{J_i\}_{i=1}^{N_T}$ using the method described in Appendix.

• Step 2 Simulate $N_T$ uniform random variables $\{U_i\}_{i=1}^{N_T}$ on [0, 1]

• Step 3 define the interval $[(1-\ell )/2,(1-\ell )/2]$ where $\ell =\lambda h$, $\lambda$ represents the frequency of jumps and h the step interval ($h=1/12$ for monthly observations, $h=1/52$ for weekly observations and $h=1/365$ for daily observations.)

• Step 4 For each each observation date $t_j$,

  • either $U_i\in [(1-\ell )/2,(1+\ell )/2]$ and we have $\widehat{X}_{t_{i+1}}= \widehat{X}_{t^-_{i+1}}+ J_i$

  • or $U_i\notin [(1-\ell )/2,(1-\ell )/2]$ and we have $\widehat{X}_{t_{i+1}}= \widehat{X}_{t^-_{i+1}}$