Kalman filter specification in Python

Question

I'm learning about state space models and I've written a simple code to apply a Kalman filter. I want to apply it to the following time-dependent model: \begin{align*} x_{t} &=\Phi_tx_{t-1}+ \Upsilon_tu_t + w_t\\ y_t &=Ax_{t} + v_t\\ w_t &\sim \mathcal{N}(0,Q)\\ v_t &\sim \mathcal{N}(0,R) \end{align*} The number of hidden variables in my model is 3, while the number of observed variables is 2. $\Phi_t$ is a $3\times 3$ matrix, $\Upsilon_t$ is a $3\times 1$ matrix, $A$ is a $2\times 3$ matrix, $Q$ is a $3\times 3$ covariance matrix and $R$ is a $2\times 2$ covariance matrix. $u_t$ is an external variable.

The code that I've written is located here. From the random data that I generated, it seems as though the confidence interval makes no sense: I plotted what I think should be the 95% confidence interval, but it is clear that (way) too many points are located outside. It seems that I am doing something wrong, but I cannot seem to find the mistake. Any help would be much appreciated.

Answer 1

I managed to get your code working! You had the vast majority of it correct, I just noticed a couple little errors.

The big problem is this line in your dynamics update

P_pred[t] = Phi[t] * P_filt[t - 1] * Phi[t].T + Q

Remember that matrix multiplication is @ not * in python, so you want (I make this mistake all the time)

P_pred[t] = Phi[t] @ P_filt[t - 1] @ Phi[t].T + Q

At the start of your filter you'll want to have the below, since P_filt[0] is what the filter uses initially.

P_filt[0] = P0 # this, 
P_pred[0] = P0 # and this.

For your Q matrix, you'll also want to have something that's positive definite. The way that you have it right now there's a zero eigenvalue. Even something like the below would be fine

Q = np.array([
    [1, 0, -0.5],
    [0, 1, -0.5],
    [-0.5, -0.5, 2],
])