Analysis, detection and correction of misspecified discrete time state space models

TL;DR

This paper introduces a method to detect and correct parameter misspecifications in discrete-time state space models, improving the accuracy of hidden state inference by leveraging innovation correlation properties.

Contribution

It proposes a novel optimization-based approach that uses innovation correlation to identify and correct parameter errors in weakly nonlinear Gaussian state space models.

Findings

The method effectively detects misspecifications across various models.

It can consistently estimate parameter bias.

Demonstrated success on models of increasing complexity.

Abstract

Misspecifications (i.e. errors on the parameters) of state space models lead to incorrect inference of the hidden states. This paper studies weakly nonlin-ear state space models with additive Gaussian noises and proposes a method for detecting and correcting misspecifications. The latter induce a biased estimator of the hidden state but also happen to induce correlation on innovations and other residues. This property is used to find a well-defined objective function for which an optimisation routine is applied to recover the true parameters of the model. It is argued that this method can consistently estimate the bias on the parameter. We demonstrate the algorithm on various models of increasing complexity.

Equations

$$\left\{\begin{array}[]{ll}x_{t}=b(\theta_{0},x_{t-1})+\beta_{\theta_{0}}\eta_{t},\\ y_{t}=h(\theta_{0},x_{t})+\sigma_{\theta_{0}}\varepsilon_{t}.\end{array}\right.$$

$$\left\{\begin{array}[]{ll}x_{t}=u_{t}(\theta_{0})+A_{\theta_{0}}x_{t-1}+\beta_{\theta_{0}}\eta_{t}\\ y_{t}=d_{t}(\theta_{0})+C_{\theta_{0}}x_{t}+\sigma_{\theta_{0}}\varepsilon_{t}\end{array}\right.$$

$$\left\{\begin{array}[]{ll}u_{t}(\theta_{0})=b(\theta_{0},\hat{x}_{t-1})-A_{\theta_{0}}\hat{x}_{t-1}\\ d_{t}(\theta_{0})=h(\theta_{0},\hat{x}^{-}_{t})-C_{\theta_{0}}\hat{x}^{-}_{t}\end{array}\right.$$

$$\left\{\begin{array}[]{lllll}\hat{x}_{t}^{-}=u_{t}(\theta)+A_{\theta}\hat{x}_{t-1}\\ \hat{x}_{t}=[I-K_{t}C_{\theta}]\hat{x}_{t}^{-}+K_{t}(y_{t}-d_{t}(\theta))\\ P_{t}^{-}=A_{\theta}P_{t-1}A_{\theta}^{\ast}+Q_{\theta}\\ P_{t}=[I-K_{t}C_{\theta}]P_{t}^{-}\\ K_{t}=P_{t}^{-}C_{\theta}^{\ast}[C_{\theta}P_{t}^{-}C_{\theta}^{\ast}+R_{\theta}]^{-1}\end{array}\right.$$

$$\left\{\begin{array}[]{lll}e_{t}:=x_{t}-\hat{x}_{t}\\ \zeta_{t}^{-}:=y_{t}-\hat{y}_{t}^{-}\\ \zeta_{t}:=y_{t}-\hat{y}_{t},\end{array}\right.$$

$$\left\{\begin{array}[]{ll}\hat{y}_{t}^{-}:=d_{t}(\theta)+C_{\theta}\hat{x}_{t}^{-}\\ \hat{y}_{t}:=d_{t}(\theta)+C_{\theta}\hat{x}_{t}\end{array}\right.$$

$$e_t=(I_{n\times n}-K_t{C}_{{\theta }_0})A_{{\theta }_0}{e}_{t-1}-K_t({\sigma }_{{\theta }_0}{\epsilon }_t+C_{{\theta }_0}{\beta }_{{\theta }_0}{\eta }_t)+{\beta }_{{\theta }_0}{\eta }_t$$

$$e_t$$

$${\mathcal{E}}_x^{ϵ}(\theta ,t)=-ϵ\left((I_{n\times n}-K_t{C}_{\theta })\frac{\partial u_t}{\partial \theta }(\theta )-K_t\frac{\partial d_t}{\partial \theta }(\theta )-K_t\frac{\partial C_{\theta }}{\partial \theta }{u}_t(\theta )\right)$$

$${\mathcal{F}}_x^{ϵ}(\theta ,t)=-ϵ\left((I_{n\times n}-K_t{C}_{\theta })\frac{\partial A_{\theta }}{\partial \theta }-K_t\frac{\partial C_{\theta }}{\partial \theta }{A}_{\theta }\right)$$

$${\mathcal{W}}_x^{ϵ}(\theta ,t)=-ϵ\left(\frac{\partial {\beta }_{\theta }}{\partial \theta }{\eta }_t-K_t{C}_{\theta }\frac{\partial {\beta }_{\theta }}{\partial \theta }{\eta }_t-K_t{\beta }_{\theta }\frac{\partial C_{\theta }}{\partial \theta }{\eta }_t-K_t\frac{\partial {\sigma }_{\theta }}{\partial \theta }{\epsilon }_t\right)$$

$${\zeta }_t=C_{\theta }{e}_t+{\sigma }_{\theta }{\epsilon }_t+{\mathcal{E}}_y^{ϵ}(\theta ,t)+{\mathcal{F}}_y^{ϵ}(\theta ,t)x_t+{\mathcal{W}}_y^{ϵ}(\theta ,t)+o(ϵ)$$

$${\mathcal{E}}_y^{ϵ}(\theta ,t)=-ϵ\frac{\partial d_t}{\partial \theta }(\theta ),{\mathcal{F}}_y^{ϵ}(\theta ,t)=-ϵ\frac{\partial C_{\theta }}{\partial \theta },{\mathcal{W}}_y^{ϵ}(\theta ,t)=-ϵ\frac{\partial {\sigma }_{\theta }}{\partial \theta }{\epsilon }_t$$

$$\mathbb{C}\mathrm{ov}({\zeta }_t,{\zeta }_{t-h})$$

$$J(\nu )=\sum _{j=1}^m\sum _{h\ge 0}{\Gamma }_{\nu }^j(h)$$

$$\hat{\Gamma}^{j}_{\nu}(h)=\frac{1}{N-1}\sum_{t=h+1}^{N}\bigg{(}\zeta_{t}^{j}-\overline{\zeta}^{j}\bigg{)}\bigg{(}\zeta_{t-h}^{j}-\overline{\zeta}^{j}\bigg{)}$$

$$\hat{J}(\nu )=\sum _{j=1}^m\sum _{h=1}^{h^{\ast }}{\hat{\Gamma }}_{\nu }^j(h)$$

$$\widehat{ϵ}=\arg \min \hat{J}(\nu ).$$

$$\displaystyle\left\{\begin{array}[]{ll}y_{t}=\alpha x_{t}+\sigma\varepsilon_{t}\\ x_{t}=\gamma x_{t-1}+\beta\eta_{t}\end{array}\right.$$

$$\displaystyle\left\{\begin{array}[]{ll}y_{t}=x_{t}+\sigma\varepsilon_{t}\\ x_{t}=\alpha\sqrt{(x_{t-1}-\gamma)}+\beta\eta_{t}\end{array}\right.$$

$$\displaystyle\left\{\begin{array}[]{ll}\frac{dS_{t}}{S_{t}}=rdt+\sqrt{v_{t}}dW_{t},\qquad S_{0}\geq 0\\ dv_{t}=\kappa(\gamma-v_{t})dt+\beta\sqrt{v_{t}}dW_{t}^{2},\qquad v_{0}\geq 0\end{array}\right.$$

$$p_{{\theta }_0}(v_t|v_{t-1})=2c{\chi }^2(2d+2,2w),$$

$$c=\frac{2\kappa }{{\beta }^2(1-e^{-\kappa \Delta })},w=c{v}_{t-\Delta }{e}^{-\kappa \Delta },d=\frac{2\kappa \gamma }{{\beta }^2}-1.$$

$$\displaystyle\left\{\begin{array}[]{ll}y_{t}=C_{t}(v_{t},S_{t},\theta_{0})+\sigma\varepsilon_{t}\\ v_{t}=\Psi(v_{t-1},\theta_{0},\Delta)+\Phi^{1/2}(v_{t-1},\theta_{0},\Delta)\eta_{t}\end{array}\right.$$

$$\Psi (v_t,{\theta }_0,\Delta )={\mathbb{E}}_{{\theta }_0}[v_{t+1}|v_t]=\gamma (1-e^{-\kappa \Delta })+e^{-\kappa \Delta }{v}_t$$

$$\Phi (v_t,{\theta }_0,\Delta )=\mathbb{V}{\mathrm{ar}}_{{\theta }_0}[v_{t+1}|v_t]=\gamma \frac{{\beta }^2}{2\kappa }(1-e^{-\kappa \Delta }{)}^2+\frac{{\beta }^2}{\kappa }{e}^{-\kappa \Delta }(1-e^{-\kappa \Delta })v_t$$

$$e_t$$

$$y_t-{\hat{y}}_t^{-}$$

$${\hat{x}}_t^{-}=A_{\theta }{\hat{x}}_{t-1}+u_t(\theta )$$

$$y_t-{\hat{y}}_t^{-}$$

Full text

1. Introduction

This paper is concerned with the following family of discrete time state space models with additive Gaussian noises:

[TABLE]

The variables $\eta_{t}\sim\mathcal{N}(0,I_{n\times 1}),\varepsilon_{t}\sim\mathcal{N}(0,I_{m\times 1})$ are assumed to be independent standard normal variables, $t\in{\mathbf{N}}^{\ast}$, $\beta_{\theta_{0}}$ (resp. $\sigma_{\theta_{0}}$) are $n\times n$ (resp. $m\times m$, with $\sigma_{0}\sigma_{0}^{\ast}$ positive definite) matrices, and $\theta_{0}$ stands for the vector of parameters of the model. The functions $b,h,\beta,\alpha$ are assumed to be differentiable. The hidden states (or unobserved signal process) $\left\{x_{t},t\in\mathbb{N}\right\}$ take value in $\mathcal{X}:={\mathbf{R}}^{n}$ and the observations $\left\{y_{t},t\in\mathbf{N}^{\ast}\right\}$ in $\mathcal{Y}:={\mathbf{R}}^{m}$. We also denote the noise covariance matrices $R_{\theta_{0}}:=\sigma_{\theta_{0}}\sigma_{\theta_{0}}^{\ast}$, $Q_{\theta_{0}}:=\beta_{\theta_{0}}\beta_{\theta_{0}}^{\ast}$ where ∗ stands for the transpose.

The aim of filtering is to make inference about the hidden state $x_{t}$ conditionally to the observations $y_{1},\cdots,y_{t}$ denoted $y_{1:t}$ thereafter. In order to do so, there are various ways to estimate the parameters $\theta_{0}$ that, in most situations of interest are unknown and have to be approximated. They may for example be estimated using standard techniques (MLEs…), or be incorporated to the set of random quantities to be estimated. To quote only one example in the recent literature, Particle Gibbs samplers have proven to be a good way to simulate the joint distribution of hidden processes and model parameters in hidden Markov chain models, see e.g. [Andrieu et al., 2010, Chopin and Singh, 2015, del Moral et al., 2016]

Here, we face a different problem: we consider the situation where $\theta_{0}$ has been uncorrectly estimated, for example using a given biased estimator $\hat{\theta}$ such that $\mathbb{E}[\hat{\theta}]=\theta=\theta_{0}+\epsilon$ (the way the estimator has been deviced is of no matter for our purposes). Our interest for these questions originated in the study of random volatility models such as Heston’s, where some parameters are difficult to estimate. We wanted to understand how errors on the model parameters could impact the volatility estimates. The detection of errors method that is the purpose of the present article first arose from statistical phenomena detected in numerical simulations. We realized soon that the phenomena were universal, and related to theoretical properties of mispecified models. Application domains include for example engineering and control where the parameters $\theta$ may be known at inception but may change to a new value $\theta_{0}$, for example due to a mechanical problem, so that $\theta$ becomes a wrong value for the true model parameters. Detecting the change from $\theta$ to $\theta_{0}$ may then be useful not only to improve the inference process, but also to detect the underlying problem.

It is well-known that using such incorrect filter models deteriorates the filter performance and may even cause the filter to diverge. Various results have been obtained in the literature on the impact of $\epsilon$ on the estimator of the hidden state; error covariance matrices have been studied and compared with the covariance matrices of the conditional distribution of $x_{t}$ and $x_{t+1}$ knowing $y_{1:t}$. These results are described in [Jazwinski, 2007], where the reader can also find a survey of the classical literature on the subject.

The aim of the present article is different: we want to take advantage of the theoretical properties of misspecified state space models, not only to understand the impact of $\epsilon$ on the estimation of the hidden states but also, ultimately, to use its statistical properties in order to get a correct set of parameters for the state space model.

The key result underlying our analysis is that misspecifications do certainly induce a biased estimator of the hidden state but also, and most importantly for our purposes, they happen to induce correlation on the innovations and other residues associated to observations. This property is used to find a well-defined objective function for which an optimisation routine is applied to recover the true parameters of the model. It is argued that this method can consistently estimate the bias on the parameter. The method is easy to implement and runs fast. We demonstrate the algorithm on various models of increasing complexity.

Discrete time state space models are notoriously ubiquitous; their use is discussed in most textbooks on filtering from the early [Kalman, 1960, Jazwinski, 2007, Sage and Melsa, 1971, Anderson and Moore, 1979] to the recent literature -we refer e.g. to [Durbin and Koopman, 2012] for a survey. Application domains of our results include, besides finance, control and engineering: ecology, economy, epidemiology, meteorology and neuroscience.

The paper is organized as follows. Section 2 presents the model assumptions and introduces various estimators and processes, including the “interpolation process” (Eq. (4)) that plays a central role in the article. Section 3 states the theoretical results. In section 4, we describe the method and in the following one demonstrate the algorithms on three examples: the first application is largely pedagogical and studies an elementary autoregressive linear model for which our approach can be easily understood. We move then to a nonlinear (square root) model, and, to conclude, apply our approach to a complex and nonlinear model, that is the Heston model, widely used in finance for option pricing and portfolios hedging. The behaviour of this last model when it comes to parameter estimation is notoriously difficult; our method behaves nevertheless quite satisfactorily. We compare finally our method and estimator (based on the interpolation process) with the estimator using the same strategy but based instead on innovations. Some concluding remarks are provided in the last section. The technical proofs are gathered in Appendix A and B.

The theoretical results on misspecified models underlying the constructions in this article were mostly obtained in the first Author’s PhD thesis [El-Kolei, 2012].

Notation: for any continuously differentiable function $g$, $[\partial g/\partial\theta]$ denotes the vector of the partial derivatives of $g$ w.r.t $\theta$.

2. The misspecified (Extended) Kalman Filter

In the linear case, the model (1) reads ($t\in\mathbf{N}^{\ast}$):

[TABLE]

If the vector of parameters $\theta_{0}$ is perfectly known, the optimal filtering $p_{\theta_{0}}(x_{t}|y_{1:t})$ is Gaussian and the Kalman filter gives exactly the two first conditional moments: $\hat{x}_{t}=\mathbb{E}[x_{t}|y_{1:t}]$ and $P_{t}=\mathbb{E}[(x_{t}-\hat{x}_{t})(x_{t}-\hat{x}_{t})^{\ast}|y_{1:t}]$. In particular, the Kalman filter estimator is the BLUE (Best Linear and Unbiased Estimator) among linear estimators.

In most real applications, the linearity assumption of the functions $h$ and $b$ is not satisfied. A linearization by a first order Taylor series expansion can be performed and the Extended Kalman filter (EKF) consists in applying the Kalman filter on this linearized model. Concretely, for the EKF, the matrix $C_{\theta_{0}}$ is the derivative of the function $h$ with respect to (w.r.t.) $x$ computed at the point $(\theta_{0},\hat{x}_{t}^{-})$ where $\hat{x}_{t}^{-}:=\mathbb{E}[x_{t}|y_{1:t-1}]$. The matrix $A_{\theta_{0}}$ is the derivative of the function $b$ w.r.t. $x$ computed at the point $(\theta_{0},\hat{x}_{t-1})$ and the functions $u_{t}(\theta_{0})$ and $d_{t}(\theta_{0})$ are defined as:

[TABLE]

In this paper, we assume that the vector of parameters $\theta_{0}$ is not perfectly known, so that the inference of the hidden state $x_{t}$ conditionally to $y_{1:t}$ is made with a parameter $\theta=\theta_{0}+\epsilon$, where $\epsilon$ stands for the error of specification. This case is frequent in practice since in general the vector of parameters is unknown and need to be estimated by an ordinary method. The resulting estimator can be biased and this bias is propagated on the estimation of the hidden state by the filter.

We now run the Kalman filter (resp. the EKF in the non linear case) with the misspecified model. The filter design reads therefore (take care that we still use the notations $\hat{x}_{t}$, $P_{t}$… for the estimator of $x_{t}$, its variance… but from now on the notation will refer to the estimators build using the biaised parameter $\theta$)

[TABLE]

with initial conditions $\hat{x}_{0}^{-}:=\mathbb{E}[x_{0}]$, $P_{0}^{-}:=Var(x_{0})$ (recall that $R_{\theta}=\sigma_{\theta}\sigma_{\theta}^{\ast}$, $Q_{\theta}=\beta_{\theta}\beta_{\theta}^{\ast}$).

We also introduce the residues, called respectively the filter error, the *innovation *and the *interpolation *processes.

[TABLE]

where

[TABLE]

Notice in particular the introduction of the interpolation process, that we specifically designed for parameters error-tracking purposes.

3. Main result

The empirical and theoretical properties of the interpolation process $(\zeta_{t})_{t\geq 1}$ (precisely, its auto-covariance) are the object of the present section. They will lead to propose a method to detect a misspecified model. Although the detection is useful in practice, we will give also a new method to approximate the bias $\epsilon$ and so to estimate the true parameter $\theta_{0}$.

Let us consider first the linear case. Recall that, by assumption, the functions $C,A,u_{t},d_{t},\beta,\sigma$ are differentiable. Recall also that if the vector of parameters is exactly known, the error a posteriori $e_{t}$ at time $t$ is given by the following formula:

[TABLE]

where $K_{t}$ is the Kalman matrix which minimizes the variance matrix $P_{t}=\mathbb{E}_{\theta_{0}}[(x_{t}-\hat{x}_{t})(x_{t}-\hat{x}_{t})^{\ast}|y_{1:t}]$ of the hidden state ( see [Kalman, 1960]).

The following Theorem gives the propagation of the error a posteriori $e_{t}$ and of $\zeta_{t}$ for the Kalman Filter when $\theta_{0}$ is not exactly known.

Theorem 3.1**.**

Consider the model (2). If $\epsilon<<1$, then:

[TABLE]

with:

[TABLE]

Additionally, the interpolation process $\zeta_{t}$ is equal to:

[TABLE]

with:

[TABLE]

Proof.

See Appendix (A). ∎

We note that the terms depending on $\epsilon$: $\mathcal{E}_{x}^{\epsilon}(\theta,t)$, $\mathcal{F}_{x}^{\epsilon}(\theta,t)$ and $\mathcal{W}_{x}^{\epsilon}(\theta,t)$ (resp. $\mathcal{E}_{y}^{\epsilon}(\theta,t)$, $\mathcal{F}_{y}^{\epsilon}(\theta,t)$ and $\mathcal{W}_{y}^{\epsilon}(\theta,t)$) are the corrective terms coming from the bias of the parameters estimates and they do not appear when the model is well specified.

Besides, we can see in Eq.(6) that at time $t$, the propagation of the state error $e_{t}$ depends on $e_{t-1}$ but also on the state variable $x_{t-1}$. Notice in particular the term $\mathcal{F}_{x}^{\epsilon}(\theta,t)x_{t}$ that contributes non trivially to the auto-correlation of the process $\zeta_{t}$; this term is proportional to $\epsilon$ but, contrary to the other terms contribution to the expansion is not proportional to a filter error term (such as $e_{t}$) or to a noise term (such as $\eta_{t}$).

For linear and gaussian state space models we can express explicitely this auto-covariance for all $t$ and $h>0$:

We can now express the auto-covariance of the interpolation processus $(\zeta_{t})_{t\geq 0}$.

Proposition 3.2**.**

Let $(\zeta_{t})_{t\geq 0}$ defined in (10) and $h>0$, we have, keeping leading contributions

[TABLE]

where the various covariance terms can be computed explicitely.

*The example of the computation of the most complex covariance term ($\mathbb{C}{\rm ov}(e_{t},e_{t-h})$) is detailed in the Appendix B. *

4. Parameter estimation: method

The main idea of the approach consists in minimizing empirically the auto-covariance between $(\zeta_{t})$ in order to reduce as far as possible the corrective terms that appear in the propagation equations (6) and (10). The results obtained on a variety of examples detailed later in the article show the meaningfulness of the approach.

Let us denote $J(\nu)$ the following objective function:

[TABLE]

where $\Gamma^{j}_{\nu}(h)$ denotes the auto-covariance of the $j^{th}$ coordinate $(\zeta^{j}_{t})$ of the vector $(\zeta_{t})$ for the lag $h$ when model parameters are chosen to be $\theta-\nu$ (recall that we know only $\theta$ and want to estimate $\theta_{0}$).

We use as estimator the empirical covariance given by:

Definition 4.1**.**

[TABLE]

where $\overline{\zeta}^{j}$ is the mean of the $\zeta_{t}^{j}$.

We will therefore minimize the following objective function

[TABLE]

As we will see in the numerical application, the choice of the lag range $h^{*}$ has no strong impact on the results.

An estimator of the bias $\epsilon$ is obtained as

[TABLE]

This means that $\nu$ is estimated in function of the tracking error $(\zeta_{t})_{t\geq 1}$.

5. Applications

5.1. Estimation of the linear AR(1) process

Let us consider the following autoregressive process:

[TABLE]

where $\alpha=3$ and $\gamma=0.9$. The noises $\varepsilon_{t}$ and $\eta_{t}$ are supposed i.i.d. with centered and standard Gaussian law. The variances $\sigma^{2}$ and $\beta^{2}$ are equal to $0.2$ and $0.1$ respectively.

We have run a Kalman filter estimation by assuming that the two parameters $\gamma$ and $\alpha$ are biased. We choose $\theta^{(0)}=(\gamma^{(0)},\alpha^{(0)})=(0.8,2.8)$ and $N=500$ to construct the function $\hat{J}(\epsilon)$ and we apply the minimization procedure to estimate the true parameter $\theta_{0}=(\gamma,\alpha)$.

The Mean Squared Error (MSE) was used to measure the quality of the estimation of $\theta_{0}$ with $MC$ (number of Monte Carlo simulations) equal to 100. The result is summarized in Table 5.1.

5.2. Estimation of a weakly Nonlinear model

Let us consider the following nonlinear model

[TABLE]

where $\alpha=5$ and $\gamma=0.008$. The noises $\varepsilon_{t}$ and $\eta_{t}$ are supposed i.i.d. with centered and standard gaussian law. The variances $\sigma^{2}$ and $\beta^{2}$ are equal to $0.2$ and $0.1$ respectively.

Since this model is nonlinear we apply an EKF estimation by assuming that the two parameters $\gamma$ and $\alpha$ are biased. For the initialisation we choose $\theta^{(0)}=(\gamma^{(0)},\alpha^{(0)})=(0.007,5.1)$ and $N=500$ to construct the function $\hat{J}(\epsilon)$ given in (13) and we apply the minimization procedure to estimate the parameter $\theta_{0}$. The results are summarized in Table 5.2.

5.3. Estimation of a strongly Nonlinear model: the Heston model

In 1993, Heston extends the Black-Scholes model by making the volatility parameter stochastic. More precisely, the volatility is modeled by a Cox Ingersoll Ross (CIR) process and the stock price follows the well-known Black-Scholes stochastic differential equation. The Heston stochastic volatility model is widely used in practice for option pricing. The reliability of the calibration of its parameters is important since a possible bias will be repercuted on the volatility estimates and, ultimately, on option prices and hedging strategies.

The model is given by

[TABLE]

where $W=\{W_{t},t\geq 0\}$ and $W^{2}=\{W^{2}_{t},t\geq 0\}$ are two correlated standard Brownian motions such that $\mathbb{C}{\rm ov}(dW_{t},dW_{t}^{2})=\rho dt$, $v_{0}$ is the initial variance, $\kappa$ the mean reversion rate, $\gamma$ the long run variance and $\beta$ the volatility of variance. We set $\theta_{0}:=(\kappa,\gamma,\beta,\rho)$.

The volatility process is always positive and cannot reach zero under the Feller condition $2\kappa\gamma>\beta^{2}$. Furthermore, under this assumption, the process $v_{t}$ has a Gamma invariant distribution $\Gamma(\alpha_{1},\alpha_{2})$ with $\alpha_{1}=\frac{2\kappa\gamma}{\beta^{2}}$ and $\alpha_{2}=\frac{\beta^{2}}{2\kappa}$.

5.3.1. Simulated Data

We sample the trajectory of the variance CIR with a time step $\Delta=1$ day over $t=1,\cdots,N$ days. Conditionally to this trajectory, we sample the trajectory of the logarithm stock price $\log S_{t}$ given by Itô’s formula and discretized by a classical Euler scheme.

For the CIR process we use the discrete time transition equation of a CIR process given by a a non-central chi-square distribution up to a constant:

[TABLE]

where $2d+2$ is the degree of freedom, $2w$ is the parameter of non-centrality and

[TABLE]

We assume that each day $t$, the observation $y_{t}$ corresponds to nine call prices for different strikes $(K_{i},T_{j})_{1\leq i,j\leq 3}$. Here, $K=\left(K_{1},K_{2},K_{3}\right)=\left(90\%,100\%,110\%\right)$ of the stock prices $S_{t}$ and $T=\left(T_{1},T_{2},T_{3}\right)=\left(0.1,0.5,1\right)$. The data length is $N=50$ days.

Then, the discrete time Heston model is given by the following nonlinear state space model with additive noises:

[TABLE]

where the functions $\Psi$ and $\Phi$ (see [Duan and Simonato, 1995]) are given by

[TABLE]

The call prices $C_{t}(v_{t},S_{t},\theta_{0})$ are computed by the Heston formula given in [Heston, 1993]. We assume that these prices are observed with Gaussian measurement error $\varepsilon_{t}$ with zero mean and variance $R=\sigma\sigma^{*}$ independent of $\theta_{0}$. These measurement errors can reflect the presence of different prices (bid-ask prices, closing prices, human errors in data handling) in financial markets.

For the vector of parameters we choose $\theta_{0}=(\kappa,\gamma,\beta,\rho)=(4,0.03,0.4,-0.5)$ which is consistent with empirical applications of daily data (see [Chen, 2007]) and the risk free interest rate $r$ is equal to 0.05.

5.3.2. Empirical detection of misspecified models

Since the Heston model is not linear, we run an EKF estimation of $v_{t}$ by assuming for convenience that only the $i^{{\rm th}}$ coordinate of the estimator of $\theta$ denoted by $\theta_{i}$ is biased, the others $(\theta_{j})_{j=1,\cdots,4}$ are equal to $\theta_{0,j}$ for $j\neq i$.

For each parameter $(\theta_{i})_{i=1,\cdots,4}$, we represent the autocorrelation of the interpolation process $(\zeta^{l}_{t})_{l=1,\ldots,9}$.

For each parameter of the Heston model, we note a presence of correlation of the interpolation process when the model is misspecified (see Figures 5.3.2 up to 5.3.2). We can also remark that this correlation is more important for the mean speed reversion parameter $\kappa$ and for the long run variance $\gamma$.

Furthermore, in order to illustrate the behaviour of the autocorrelation with respect to the bias $\epsilon$ we apply an EKF estimation by considering the three following cases (only for the speed mean reversion parameter, the conclusion is the same for the others parameters): a) $\kappa=4$ (that is the model is well-specified) ; b) $\kappa=4.48$ ; c) $\kappa=4.96$. The comparison is illustrated in Figure 5.3.2. As expected, we observe that no correlation appaers when the model is well-specified (that is in Case a) and in return when a bias is introduced a correlation of the interpolation process appears and most importantly this correlation growths with the bias (see Figure 5.3.2 Case b and c).

5.3.3. Parameter estimation

In Figure 5.3.3, we represented the objective function $\hat{J}(\nu)$ defined in (13) with respect to the parameters of the Heston model. We represent only $\hat{J}(\nu)$ for the long run variance parameter $\gamma$ since for the others parameters the result is the same. We can see that the function $\hat{J}$ is minimal for the true value of $\gamma$, that is $\gamma=0.03$ (see Figure 5.3.3).

5.3.4. Estimation of the long run variance in the Heston model

In a first step, we have run an EKF estimation by assuming that only the long run variance parameter $\gamma$ is biased. We choose $\gamma=0.025$ and we recall that its true value is $0.03$. The number of observations used for the construction of the function $\hat{J}(\nu)$ given in (13) is here $N=100$ and we apply the minimization procedure to recover the parameter $\gamma$.

The MSE was used to measure the quality of the estimation of the parameter $\gamma$ with $MC$ equal to 50. The results are summarized in Table 5.3.4 and Figure 5.3.4.

[TABLE]

Figure 7. Boxplot of the estimation of $\beta$ for $MC=50$ and $h^{*}=2$.

5.3.5. Sensibility w.r.t the lag $h$

In order to see the impact of the lag $h$ on the autocorrelation, we have run our approach for different lags $h$ and compute the MSE (with MC=50). We note (see Table 5.3.5) that the choice of $h$ has not a strong impact on the results. Hence, for the next numerical application we choose $h$ equal to $8$.

In Table 5.3.5, we illustrate the MSE for different number of observations, $N=20$ up to $110$. We see that the estimation is very bad for a small $N$ which is not suprising since in this case the empirical estimator $\hat{J}(\nu)$ of $J(\nu)$ is not consistent. As we expect, the MSE decreases with the number of observations.

5.3.6. Estimation of the Heston model

In this part, we want to estimate all parameters of the Heston model. So, we consider that all parameters are biased with different bias (see $\theta^{(0)}$ on Table 5.3.6 and 5.3.6) and that the true parameter $\theta_{0}$ is given by $\theta_{0}=(4,0.03,0.4,-0.5)$.

In Table 5.3.6 we repeat our procedure of estimation with MC equal to $50$. We note that our approach leads to estimate simultaneously all parameters. We also note that the long run variance parameter $\gamma$ and the speed mean reversion parameter $\kappa$ are easier to estimate than the others parameters. Furthermore, we have seen in Figures 5.3.2 and 5.3.2 of Section 5.3.2 that the correlation was more important for these two parameters.

For $m=1,\cdots,MC$, the choice of the initial condition for $\theta_{m}^{(0)}$ is \theta_{m}^{(0)}\in\bigg{[}3.8+\tilde{\sigma}\mathcal{N}(0,1);0.025+\tilde{\sigma}\mathcal{N}(0,1);0.43+\tilde{\sigma}\mathcal{N}(0,1));-0.55+\tilde{\sigma}\mathcal{N}(0,1)\bigg{]} with $\tilde{\sigma}=10^{-4}$ and where $\mathcal{N}(0,1)$ stands for the centered and standard gaussian law.

6. Comparison with the use of standard innovations

Methods for detection of departures from optimality are usually based on the innovation process $(\zeta^{-}_{t})_{t\geq 1}$. Performance analysis of Kalman filters based on the innovation was introduced in [Wei et al., 1991]. In their papers, the authors propose a test based on the innovations for fault detection and a two-step Kalman filtering procedure to estimate the parameters. Let us mention also [Grewal and Andrews, 2015] where in page 370 the authors give a short discussion on detecting unmodeled state dynamics by Fourier analysis of the filter innovations.

In this part, we compare our minimisation routine (13) with the analogous minimisation routine when one replaces the interpolation process $(\zeta_{t})_{t\geq 1}$ with the innovation process $(\zeta^{-}_{t})_{t\geq 1}$ in order to estimate the parameters.

For this comparison we use the three models defined in the previous section and assume that only one parameter is biased for each model.

We note that the MSE is significantly smaller when one uses the interpolation process instead of the standard innovations and most importantly using the interpolation process to correct the bias is better for complex models with nonlinear effects. The results are summarized in Table 6.

7. Conclusion

In this paper, we propose a new approach to detect and estimate parameters of weakly nonlinear hidden states models. These models are supposed to be misspecified due to the choice of uncorrect parameters. We propose to exploit the autocorrelation of a suitably defined interpolation process based on the estimate of the hidden state with the biased parameters. We vary then the model parameters around the initial misspecified value and apply an optimization procedure to minimize the auto-covariance of this process. We show that this approach leads to detect misspecified models and to estimate the parameters. The computing time is fast and the implementation is easy. Furthermore, we note that the autocorrelation lag parameter $h$ has not a strong impact on the results. All results are illustrated on various models of increasing complexity and in particular on the Heston model widely used in practice for portfolio hedging.

Appendix A Proof of Theorem 3.1:

The proof is essentially based on a first order Taylor expansion of the functions $b$ and $h$ with respect to $\theta$. We have

[TABLE]

where we used $\hat{x}_{t}-\hat{x}_{t}^{-}=-K_{t}C_{\theta}\hat{x}_{t}^{-}-K_{t}(y_{t}-d_{t}(\theta))=-K_{t}(y_{t}-\hat{y}_{t}^{-}).$ Since $u_{t}(\theta_{0})=u_{t}(\theta)-\epsilon\frac{\partial u}{\partial\theta}(\theta)+o(\epsilon)$, and similarly for the other functions of $\theta$, we get

[TABLE]

Furthermore,

[TABLE]

so that:

[TABLE]

Rewriting

[TABLE]

we get:

[TABLE]

Define,

[TABLE]

we obtain:

[TABLE]

By combining Eq.(26) and Eq.(A), we have:

[TABLE]

where,

[TABLE]

One can deduce the Propagation of the interpolations (or residues a posteriori):

[TABLE]

By defining:

[TABLE]

Eq.(10) follows. ∎

Appendix B Covariances in Proposition3.2:

We have, up to $o(\epsilon)$ terms that are neglected

[TABLE]

Denote

[TABLE]

with the convention $\prod_{l=0}^{-1}\equiv 1$, we can rewrite

[TABLE]

and setting

[TABLE]

we get

[TABLE]

In a similar manner, we can rewrite the hidden state $x_{t}$ as

[TABLE]

So,

[TABLE]

By using the fact that

[TABLE]

(an equation that follows by induction from the definition of $x_{n}$), we obtain that Eq.(30) is equal to

[TABLE]

and that Eq.(31) is equal to

[TABLE]

By replacing $x_{l}$ in (32) by (28), we have that Eq.(32) is equal to

[TABLE]

In a similar manner, Eq.(33) is equal to

[TABLE]

And, Eq.(34) is equal to

[TABLE]

In a same way, we obtain $\mathbb{C}{\rm ov}(e_{t},x_{t-h})$ and $\mathbb{C}{\rm ov}(x_{t},e_{t-h})$ ∎.

Acknowledgement: We thank Patricia Reynaud-Bouret and N. Chopin for their suggestions and their interest about this work.