I am doing an Extended Kalman Filter for a vehicle traking problem, I have read some articles about this and I have a doubt related to how the error covariance matrix (Q) is initialized.
Let suppose that the dynamic model models a vehicle moving at constant velocity, where the state is $x_{k} = [lon,lat,vx,vy]$ and $P_{k}$ is the covariance matrix. In this case the error covariance matrix, Q, is initialized as follow:
\begin{equation}
Q = q_{0}\Delta t
\begin{bmatrix}
\Delta t^{2}/3 & \Delta t/2 \\
\Delta t/2 & 1
\end{bmatrix}
\end{equation}
In the article The Interacting Multiple Model Algorithm for Accurate State Estimation of Maneuvering Targets, the scalar $q_{0}$ is called as filter plant noise spectral density.
I don't know how to proof this initialization and how to use it in another kind of dynamic model.
In summary my doubt is how to correctly initialize the error covariance matrix in an Extended Kalman Filter.
Best Answer
There are different approaches to initialize the Q matrix.
First off all let's make it clear for the initialization example in your question.
\begin{equation} Q = q_{0}\Delta t \begin{bmatrix} \Delta t^{2}/3 & \Delta t/2 \\ \Delta t/2 & 1 \end{bmatrix} \end{equation}
This approach is well explained in Kalman and Bayesian Filters in Python (Chapter 7.3.1 Continuous White Noise Model).
It works for systems with state vectors containing derivatives such as position, velocity, acceleration, jerk... The last element in such a sequence is meant to stay constant during the prediction step. The system noise is modelled as a continuous white noise $Q_c$ with a spectral density $q_0$ and should correct this modelling issue (that's why the spectral density is applied only to the right bottom element). For a state vector containing position and velocity it would look like
\begin{equation} Q_c = \begin{bmatrix} 0 & 0\\ 0 & 1\\ \end{bmatrix} q_{0} \end{equation}
To define the amount of noise for a discrete time slot $\Delta t$ the continuous noise should be integrated over the interval $[0 .. \Delta t]$ :
\begin{equation} Q=\int_{0}^{\Delta t} F Q_c F^{T} dt = q_0 \begin{bmatrix} \Delta t^{3}/3 & \Delta t^{2}/2\\ \Delta t^{2}/2 & \Delta t\\ \end{bmatrix} \end{equation}
F is the transition matrix
\begin{equation} F = \begin{bmatrix} 1 & \Delta t\\ 0 & 1\\ \end{bmatrix} \end{equation}
To use this approach you need to make sure your state elements are in a right order.
If your state vector is
\begin{equation} x = \begin{bmatrix} x \\ vx \\ y \\ vy \end{bmatrix} \end{equation}
then Q would be
\begin{equation} Q = q_{0} \begin{bmatrix} \Delta t^{3}/3 & \Delta t^{2}/2 & 0 & 0\\ \Delta t^{2}/2 & \Delta t & 0 & 0\\ 0 & 0 & \Delta t^{3}/3 & \Delta t^{2}/2\\ 0 & 0 & \Delta t^{2}/2 & \Delta t \end{bmatrix} \end{equation}
So you need to sort the elements of Q corresponding to your state vector.
Regarding the spectral density: the exact value is usually unknown and is used to tune the system, til the performance fits to your expectations.
Another approach
Another easier intuition when dealing with Q is to interpret its elements as uncertainty being added to your state elements during the prediction step. For each element of the state vector ask yourself what is the biggest prediction error when using the transition matrix F. Assuming that the error is normally distributed you can use the 3-sigma rule. Sigma squared gives you the corresponding diagonal element of the Q matrix.
\begin{equation} Q_{ii} = (max(\varepsilon_{i}) /3)^{2} \end{equation}
The better the model, the smaller the prediction error, and the smaller the Q-element.