State-space models and the Kalman filter

Before this lesson: simple exponential smoothing (lesson 04), AR models and the backshift notation (lesson 05), forecast intervals and ARIMA(0,1,1) (lesson 06). From the probability-and-statistics track: the multivariate Gaussian and its conditional distribution, maximum likelihood, and Bayes' rule with a prior and a posterior. Covariance matrices are positive semidefinite (linear-algebra lesson 07). Time: about 110 minutes, plus the notebook.

---

Why this matters #

• One framework for many models. Local level, trend, seasonal, ARMA and regression models with changing coefficients are all state-space models. One algorithm filters, forecasts and computes the likelihood for all of them. statsmodels computes the exact likelihood of ARIMA and SARIMA models with a Kalman filter (its SARIMAX class is a state-space model).
• Missing data and irregular times need no special code: the filter skips the update when there is no observation.
• Tracking and navigation. The Kalman filter (Kalman, 1960) fuses GPS, inertial sensors and odometry in aircraft, phones and robots. It is the base of SLAM and of data assimilation in weather forecasting (ensemble Kalman filter, see the ml-for-science track).
• Deep learning. The layers of S4 and Mamba are linear state-space recurrences (deep-learning lesson 16). They use the state equation of this lesson without noise and learn $F$ and $H$ by gradient descent.

1. Predict first #

2. The picture #

A hidden state seen through noise #

A state-space model has a state $x_t$ that we do not observe and an observation $y_t$ that we do observe. The state carries everything about the past that matters for the future. The state moves by a linear rule $F$ plus noise; the observation is a linear function $H$ of the state plus noise:

The arrows are the only dependencies: $x_t$ depends on the past only through $x_{t-1}$ (a Markov chain), and $y_t$ depends only on $x_t$. This is the same graph as a hidden Markov model; here the state is a continuous vector and everything is linear and Gaussian.

The example of the figure. A target moves along a line. Its state is $x_t = (\text{position}, \text{velocity})$. In each time step the position grows by the velocity, and a small random acceleration changes both. A sensor measures the position only, with noise of standard deviation 3. The velocity is never measured.

Interactive figure: loading…

The filter repeats two steps.

1. Predict. Move the last estimate forward with the motion rule: new position = old position + velocity. The uncertainty grows, because the random acceleration adds the variance $Q$.
2. Update. A new observation $y_t$ arrives. The difference between it and the predicted position is the innovation $\nu_t$ (the "surprise"). Move the estimate by a fraction of the innovation. The uncertainty shrinks, because we have one more measurement.

The fraction is the Kalman gain $K_t$. It is not a tuning knob: the filter computes it from the two variances. With the true $Q$ and $R$, the filtered position has an RMSE well below the RMSE of the raw observations: the filter averages the noise away and still follows the turns of the target.

One update is Bayes' rule for two Gaussians #

Look at one update in one dimension. Before $y_t$ arrives, the prediction says $x_t \sim \mathcal N(m, P)$: this is the prior. The sensor says $y_t = x_t + v_t$ with $v_t \sim \mathcal N(0, R)$: as a function of $x_t$, this is the likelihood. Bayes' rule multiplies the two. The product of two Gaussian densities in $x$ is again a Gaussian: the posterior.

Interactive figure: loading…

The posterior mean is a weighted average of the prediction and the observation. Each one gets a weight proportional to its precision (one over its variance). So the more precise source wins, and the combined estimate is more precise than either source.

3. Definitions and main results #

3.1 The linear Gaussian state-space model #

Definition. A linear Gaussian state-space model for a series $y_1, \dots, y_n$ (each $y_t$ a vector of size $k$) has a hidden state $x_t$ (a vector of size $d$) and two equations:

The first line is the state equation (or transition equation); the second is the observation equation (or measurement equation).

• $F$ ($d \times d$) is the transition matrix; $H$ ($k \times d$) is the observation matrix.
• $Q$ is the process noise (state noise) covariance; $R$ is the observation noise covariance. Both are positive semidefinite.
• Assumptions: the noises $w_1, w_2, \dots$ and $v_1, v_2, \dots$ are independent of each other, independent over time, and independent of the initial state $x_0 \sim \mathcal N(\hat x_{0\mid 0}, P_{0\mid 0})$.

$F$, $H$, $Q$ and $R$ may depend on $t$ (then we write $F_t$, and so on); every result below still holds. A known input $u_t$ can be added, $x_t = F x_{t-1} + B u_t + w_t$; this only shifts the predicted mean.

Examples. Each model below is one choice of the four matrices.

Model	State $x_t$	$F$	$H$	$Q$	$R$
local level (random walk plus noise)	level $\ell_t$	$1$	$1$	$\sigma_\eta^2$	$\sigma_\varepsilon^2$
local linear trend	$(\ell_t, b_t)$ level and slope	$\begin{pmatrix}1&1\\0&1\end{pmatrix}$	$\begin{pmatrix}1&0\end{pmatrix}$	$\operatorname{diag}(\sigma_\xi^2, \sigma_\zeta^2)$	$\sigma_\varepsilon^2$
constant velocity (tracking)	(position, velocity)	$\begin{pmatrix}1&1\\0&1\end{pmatrix}$	$\begin{pmatrix}1&0\end{pmatrix}$	$q\begin{pmatrix}1/4&1/2\\1/2&1\end{pmatrix}$	$\sigma_v^2$
seasonal, period $m = 4$	$(s_t, s_{t-1}, s_{t-2})$	$\begin{pmatrix}-1&-1&-1\\1&0&0\\0&1&0\end{pmatrix}$	$\begin{pmatrix}1&0&0\end{pmatrix}$	$\operatorname{diag}(\sigma_\omega^2, 0, 0)$	$\sigma_\varepsilon^2$
AR(2), companion form	$(y_t, y_{t-1})$	$\begin{pmatrix}\phi_1&\phi_2\\1&0\end{pmatrix}$	$\begin{pmatrix}1&0\end{pmatrix}$	$\operatorname{diag}(\sigma^2, 0)$	$0$

How to read the rows:

• Local level. $\ell_t = \ell_{t-1} + \eta_t$ is a random walk; $y_t = \ell_t + \varepsilon_t$ adds measurement noise. This is the model behind SES (section 3.4).
• Local linear trend. $\ell_t = \ell_{t-1} + b_{t-1} + \xi_t$ and $b_t = b_{t-1} + \zeta_t$: the slope also moves slowly. This is the stochastic version of Holt's linear method (lesson 04).
• Constant velocity. The same $F$. The random acceleration $a_t \sim \mathcal N(0, q)$ adds $a_t/2$ to the position and $a_t$ to the velocity, so $w_t = a_t\,(1/2, 1)^\top$ and $Q = q\,(1/2, 1)^\top(1/2, 1)$.
• Seasonal. The first row says $s_t = -(s_{t-1} + s_{t-2} + s_{t-3}) + \omega_t$: four consecutive seasonal values sum to zero, up to noise. To get "level + season", stack the level and the seasonal state into one vector, put $F$ in block-diagonal form and use $H = (1, 1, 0, 0)$.
• AR(2). The first row of $F x_{t-1}$ is $\phi_1 y_{t-1} + \phi_2 y_{t-2}$; the second row copies $y_{t-1}$ down. The eigenvalues of $F$ are the inverses of the roots of $\phi(z)$ (lesson 05), so the model is stationary exactly when they lie inside the unit circle. Every ARMA($p, q$) model has such a form, with a state of size $\max(p, q + 1)$ (Durbin and Koopman, ch. 3).

Interactive figure: loading…

3.2 The scalar update, step by step #

Start with the local level model. All quantities are numbers. We write $\hat x_{t\mid s}$ for the mean of $x_t$ given $y_1, \dots, y_s$ and $P_{t\mid s}$ for its variance.

Predict. Suppose $x_{t-1} \mid y_{1:t-1} \sim \mathcal N(\hat x_{t-1\mid t-1}, P_{t-1\mid t-1})$. Then $x_t = x_{t-1} + w_t$ with independent $w_t \sim \mathcal N(0, Q)$, so

Update. Write $m = \hat x_{t\mid t-1}$ and $P = P_{t\mid t-1}$. Then:

1. Prior times likelihood, as a function of $x$: $p(x \mid y_{1:t}) \propto \exp\!\big(-\tfrac{(x - m)^2}{2P}\big)\exp\!\big(-\tfrac{(y_t - x)^2}{2R}\big)$.
2. The exponent is $-\tfrac12\big[\tfrac{(x - m)^2}{P} + \tfrac{(x - y_t)^2}{R}\big]$. Expand and keep the terms with $x$: $-\tfrac12\big[x^2\big(\tfrac1P + \tfrac1R\big) - 2x\big(\tfrac mP + \tfrac{y_t}R\big)\big] + \text{const}$.
3. A Gaussian $\mathcal N(\mu, V)$ has exponent $-\tfrac12\big[\tfrac{x^2}{V} - \tfrac{2x\mu}{V}\big] + \text{const}$. Compare the coefficients (this is "completing the square"):

4. So the posterior variance is $V = \dfrac{PR}{P + R}$ and the posterior mean is $\mu = V\big(\tfrac mP + \tfrac{y_t}R\big) = \dfrac{Rm + Py_t}{P + R}$.
5. Write $K = \dfrac{P}{P + R}$. Then $\dfrac{R}{P + R} = 1 - K$, and step 4 becomes

The first form says: a weighted average of prediction and observation. The second form says: the prediction plus the gain times the innovation $\nu_t = y_t - m$. Step 3 says that precisions add: $1/P_{t\mid t} = 1/P_{t\mid t-1} + 1/R$. So $P_{t\mid t}$ is smaller than both $P_{t\mid t-1}$ and $R$.

3.3 The Kalman filter #

For the vector model the same two steps hold with matrices.

Which covariance is which. $P_{t\mid t-1}$ is the uncertainty of the state before we see $y_t$ (the prediction). $P_{t\mid t}$ is the uncertainty after we see it (the filtered estimate). $S_t$ is the uncertainty of the observation $y_t$ before we see it: state uncertainty mapped by $H$, plus sensor noise. With $d = k = 1$ and $F = H = 1$ the box is exactly section 3.2.

Theorem (the Kalman filter is exact). Under the assumptions of section 3.1, $x_t \mid y_{1:t-1} \sim \mathcal N(\hat x_{t\mid t-1}, P_{t\mid t-1})$ and $x_t \mid y_{1:t} \sim \mathcal N(\hat x_{t\mid t}, P_{t\mid t})$ for every $t$, and $y_t \mid y_{1:t-1} \sim \mathcal N(H\hat x_{t\mid t-1}, S_t)$.

Proof by Gaussian conditioning

We use one fact about jointly Gaussian vectors $(a, b)$: the conditional distribution of $a$ given $b$ is Gaussian with mean $\mu_a + \Sigma_{ab}\Sigma_{bb}^{-1}(b - \mu_b)$ and covariance $\Sigma_{aa} - \Sigma_{ab}\Sigma_{bb}^{-1}\Sigma_{ba}$ (the multivariate Gaussian lesson of the probability-and-statistics track). The proof is by induction over $t$.

1. Start. $x_0 \sim \mathcal N(\hat x_{0\mid 0}, P_{0\mid 0})$ by assumption.
2. Predict. Assume $x_{t-1} \mid y_{1:t-1} \sim \mathcal N(\hat x_{t-1\mid t-1}, P_{t-1\mid t-1})$. The noise $w_t$ is independent of $x_{t-1}$ and of $y_{1:t-1}$. A linear map of a Gaussian plus an independent Gaussian is Gaussian, with mean $F\hat x_{t-1\mid t-1}$ and covariance $FP_{t-1\mid t-1}F^\top + Q$.
3. Joint of state and observation. Given $y_{1:t-1}$, the pair $(x_t, y_t)$ with $y_t = Hx_t + v_t$ is jointly Gaussian ($v_t$ is independent of everything else), with

$$\mu = \begin{pmatrix}\hat x_{t\mid t-1}\\ H\hat x_{t\mid t-1}\end{pmatrix}, \qquad \Sigma = \begin{pmatrix} P_{t\mid t-1} & P_{t\mid t-1}H^\top \\ HP_{t\mid t-1} & HP_{t\mid t-1}H^\top + R\end{pmatrix}.$$

The lower-right block is $S_t$. This already proves the claim about $y_t \mid y_{1:t-1}$. 4. Condition on $y_t$. With $a = x_t$ and $b = y_t$, the fact above gives the mean $\hat x_{t\mid t-1} + P_{t\mid t-1}H^\top S_t^{-1}(y_t - H\hat x_{t\mid t-1}) = \hat x_{t\mid t-1} + K_t\nu_t$ and the covariance $P_{t\mid t-1} - P_{t\mid t-1}H^\top S_t^{-1}HP_{t\mid t-1} = (I - K_tH)P_{t\mid t-1}$. This is the update step, and the induction continues. $\blacksquare$

Two remarks. First, the covariances $P_{t\mid t-1}$, $P_{t\mid t}$, $S_t$ and the gain $K_t$ do not depend on the data $y_t$: they can be computed before the data arrive. Second, without the Gaussian assumption (noises only uncorrelated, with the same means and covariances), the same recursion still gives the best linear estimate of $x_t$ in mean squared error.

Why the velocity changes when only the position is observed. In the constant-velocity model, $H = (1, 0)$, so $K_t = P_{t\mid t-1}H^\top/S_t$ is the first column of $P_{t\mid t-1}$ divided by $S_t$. The second entry of that column is the covariance of position and velocity. The predict step makes it nonzero (a fast target moves far), so a surprise in the position also corrects the velocity. The worked example in section 4 shows this with numbers.

Numerical form and cost

The form $P_{t\mid t} = (I - K_tH)P_{t\mid t-1}$ can lose symmetry and positive semidefiniteness in floating point. The Joseph form $P_{t\mid t} = (I - K_tH)P_{t\mid t-1}(I - K_tH)^\top + K_tRK_t^\top$ is a sum of two positive semidefinite matrices, so it stays positive semidefinite, and it is correct for any gain $K_t$, not only the optimal one. Square-root filters propagate a Cholesky factor of $P$ instead (linear-algebra lesson 07). Solve with $S_t$ (Cholesky) instead of forming $S_t^{-1}$. One step costs $O(d^3 + k^3)$ operations; this is why a state of size $10^8$ (a weather model) needs the ensemble methods of the ml-for-science track.

3.4 Steady state: the local level filter is exponential smoothing #

For the local level model the covariance recursion does not involve the data. Combine predict and update into one map for $P_t = P_{t\mid t-1}$:

Theorem (steady-state gain). Let $Q > 0$, $R > 0$ and $q = Q/R$ (the signal-to-noise ratio). For every start $P_{0\mid 0} \ge 0$, $P_{t\mid t-1}$ converges to $P = \tfrac12\big(Q + \sqrt{Q^2 + 4QR}\big)$, and the gain converges to

Proof.

1. Fixed point. $P = g(P)$ means $P(P + R) = PR + Q(P + R)$, that is $P^2 - QP - QR = 0$. The only positive root is $P = \tfrac12\big(Q + \sqrt{Q^2 + 4QR}\big)$.
2. Gain. Divide by $R$: $p = P/R$ solves $p^2 = qp + q$, so $p = \tfrac12\big(q + \sqrt{q^2 + 4q}\big)$. The claim $K_\infty = p/(p + 1) = p - q$ is equivalent to $p = (p - q)(p + 1) = p^2 + p - qp - q$, which is $p^2 = qp + q$. So $K_\infty = p - q = \tfrac12\big(\sqrt{q^2 + 4q} - q\big)$.
3. Convergence. $P_1 = P_{0\mid 0} + Q \ge Q$ and $g(P) \ge Q$ for every $P \ge 0$, so every $P_t \ge Q$, and the fixed point is $\ge Q$ too. The derivative is $g'(P) = R^2/(P + R)^2$, so on $[Q, \infty)$ it is at most $R^2/(Q + R)^2 = 1/(1 + q)^2 < 1$. By the mean value theorem,

so the error shrinks geometrically. $\blacksquare$

The link to SES. With a constant gain $K$, the update is

This is the SES recursion $\ell_t = \alpha y_t + (1 - \alpha)\ell_{t-1}$ of lesson 04 with $\alpha = K_\infty$. The forecasts agree too: both are flat, $\hat y_{t+h\mid t} = \hat x_{t\mid t}$. So SES is the steady-state Kalman filter of the local level model, and $\alpha$ is decided by the signal-to-noise ratio $q$:

$q = Q/R$	0.01	0.1	0.5	1	10
$K_\infty = \alpha$	0.095	0.270	0.500	0.618	0.916

Small $q$ (the level moves little, the noise is large) gives a small $\alpha$: average over a long past. Large $q$ gives $\alpha$ near 1: trust the last observation. As $q \to 0$, $K_\infty \approx \sqrt q$; as $q \to \infty$, $K_\infty \to 1$.

Interactive figure: loading…

SES intervals and ARIMA(0,1,1)

In steady state the innovations $\nu_t$ have constant variance $\sigma^2 = S = P + R$ and are white noise. From $y_t = \hat x_{t-1\mid t-1} + \nu_t$ and $\hat x_{t\mid t} = \hat x_{t-1\mid t-1} + \alpha\nu_t$ we get $y_t - y_{t-1} = \nu_t - (1 - \alpha)\nu_{t-1}$: an ARIMA(0,1,1) with $\theta = \alpha - 1$, as in lesson 06. The $h$-step forecast variance of the state-space model is $P_{t\mid t} + hQ + R$. It equals the ARIMA formula $\sigma^2\big(1 + (h - 1)\alpha^2\big)$: both are $\sigma^2$ for $h = 1$, and each extra step adds $Q$ in one form and $\sigma^2\alpha^2 = P^2/(P + R)$ in the other. These are equal, because the fixed-point equation says $Q = P - PR/(P + R) = P^2/(P + R)$.

3.5 Forecasting and missing data #

Forecasting. To forecast from origin $t$, run the predict step $h$ times without updates:

so $\hat x_{t+h\mid t} = F^h\hat x_{t\mid t}$. The forecast of the observation is $\hat y_{t+h\mid t} = H\hat x_{t+h\mid t}$ with covariance $HP_{t+h\mid t}H^\top + R$. Under the model, a 95 % prediction interval for a scalar $y_{t+h}$ is $\hat y_{t+h\mid t} \pm 1.96\sqrt{HP_{t+h\mid t}H^\top + R}$. For the local level model, $P_{t+h\mid t} = P_{t\mid t} + hQ$: the variance grows linearly in $h$, and the fan widens like $\sqrt h$. These intervals assume that $F, H, Q, R$ are known; using estimated values makes them a little too narrow (lesson 10 checks the coverage of intervals).

Missing data. If $y_t$ is missing, skip the update: $\hat x_{t\mid t} = \hat x_{t\mid t-1}$ and $P_{t\mid t} = P_{t\mid t-1}$. This is exact, because "no observation" means "condition on nothing new". For the local level model the variance grows by $Q$ in each missing step. At the first observation after the gap, $P_{t\mid t-1}$ is large, so the gain is near 1 and the estimate jumps to the data. The likelihood (next section) simply has no term for a missing time. If only some entries of a vector $y_t$ are missing, keep the rows of $H$ and $R$ that belong to the observed entries.

Irregular time steps. If the time between observations is $\Delta_t$, use $F$ and $Q$ for that step. For the local level model, a random walk observed after $\Delta_t$ units has $Q_t = \Delta_t\,\sigma_\eta^2$.

3.6 The likelihood: the prediction-error decomposition #

To estimate unknown parameters $\theta$ (for example $Q$ and $R$, or the $\phi_i$ of an AR model), we need the likelihood $p(y_1, \dots, y_n \mid \theta)$. The observations are dependent, so it is not a product of marginal densities. The chain rule of probability gives

and the theorem of section 3.3 says $y_t \mid y_{1:t-1} \sim \mathcal N(H\hat x_{t\mid t-1}, S_t)$. So

This is the prediction-error decomposition: the filter turns dependent observations into innovations $\nu_t$ that are independent with known covariances $S_t$. One pass of the filter gives the exact Gaussian log-likelihood in $O(n)$ steps. Maximise it with a numerical optimiser; write variances as $Q = e^{\theta_1}$, $R = e^{\theta_2}$ so that they stay positive.

Checking the model. If the model is correct, the standardised innovations $\nu_t/\sqrt{S_t}$ are independent $\mathcal N(0, 1)$. Check them like residuals: plot them, look at their ACF and run the Ljung–Box test (lesson 02).

The Nile example (Durbin and Koopman, ch. 2). The annual flow of the Nile at Aswan, 1871–1970, in $10^8\,\text{m}^3$, fitted with a local level model. Maximising the likelihood gives $\hat R = 15\,099$ and $\hat Q = 1\,469$, so $\hat q = 0.097$ and $K_\infty = 0.267$: the level moves slowly compared with the year-to-year noise. The maximum log-likelihood is $-632.55$ (with the diffuse start of section 3.7).

Interactive figure: loading…

3.7 Initialisation #

The filter needs $\hat x_{0\mid 0}$ and $P_{0\mid 0}$. Often we know nothing about the start. A diffuse prior uses a very large $P_{0\mid 0}$. For the local level model, $P_{1\mid 0}$ is then huge, $K_1 \approx 1$, and so $\hat x_{1\mid 1} \approx y_1$ and $P_{1\mid 1} \approx R$: the first estimate is the first observation, with the uncertainty of one observation. The first innovation has an almost infinite variance and carries no information, so the diffuse log-likelihood leaves out the first term.

In practice: (1) the exact diffuse start for the local level model is $\hat x_{1\mid 1} = y_1$, $P_{1\mid 1} = R$, and the likelihood starts at $t = 2$ (this lesson and Durbin and Koopman use it); or (2) set $P_{0\mid 0}$ to a very large number, for example $10^7$ for the Nile data, and leave out the first likelihood term (an approximate diffuse start). Different starts change the estimates a little. For the Nile data, statsmodels 0.15 UnobservedComponents gives $\hat R = 15\,078$, $\hat Q = 1\,479$. The first filtered values have wide bands for the same reason: early on, there is little data. For a stationary state (such as an AR part) use instead its stationary distribution as the start.

3.8 Smoothing: using the whole series #

The filter estimates $x_t$ from $y_1, \dots, y_t$. After the whole series is observed, we can estimate each $x_t$ from all of $y_1, \dots, y_n$: this is smoothing. The Rauch–Tung–Striebel (RTS) smoother runs backward after the filter, for $t = n - 1, \dots, 1$:

It starts at $t = n$, where smoother and filter agree. In words: the smoother corrects the filtered estimate by what the future revealed about the prediction error at $t + 1$.

Smoothed variances are smaller. $P_{t\mid n} \preceq P_{t\mid t}$ (in the order of positive semidefinite matrices): conditioning a Gaussian on more data never increases its covariance. On the Nile data, the smoothed variance in the middle of the series is about 2 330, against 4 030 for the filter.

Filter or smoother? The smoother uses future data. Use it to describe the past: decomposition into trend and season, filling gaps, finding the date of a level shift, or anomaly analysis (lesson 11). Never use it to make or evaluate forecasts: in a backtest (lesson 03), a smoothed state at the forecast origin has seen the test period. That is leakage.

Beyond linear and Gaussian

• Extended Kalman filter (EKF). For $x_t = f(x_{t-1}) + w_t$, $y_t = h(x_t) + v_t$, linearise $f$ and $h$ around the current estimate (Jacobians) and run the Kalman steps. Works when the nonlinearity is mild over the size of the uncertainty.
• Unscented Kalman filter (UKF). Push $2d + 1$ chosen "sigma points" through $f$ and $h$ and match the mean and covariance. No Jacobians; often more accurate than the EKF.
• Particle filter. Represent the distribution by many weighted samples; predict by simulating the model, update by weighting with the likelihood, and resample. Works for any model, but needs many particles when the state has many dimensions.
• Ensemble Kalman filter and 4D-Var. For states with millions of entries (weather), estimate $P$ from an ensemble of model runs, or solve a variational problem over a time window. See the ml-for-science track.
• Deep state-space layers. S4 and Mamba use $x_t = \bar Ax_{t-1} + \bar Bu_t$, $y_t = Cx_t + Du_t$ (deep-learning lesson 16): $\bar A$ is our $F$, $C$ is our $H$, $u_t$ is a known input. There is no noise and no update step; the matrices are learned by gradient descent.

Answers to section 1. (1) About 2.1 m. Precisions add: $1/9 + 1/9 = 2/9$, so the variance is $4.5$ and the standard deviation is $\sqrt{4.5} \approx 2.12$ m, which is $3/\sqrt 2$. (2) During the gap only the predict step runs, so the variance grows by $Q$ (local level) or faster (models with a velocity) at every step. At the first observation after the gap the gain is large, so the estimate jumps toward the data and the variance drops at once. (3) Yes: the local level model (random walk plus noise). SES with $\alpha = K_\infty$ is its steady-state Kalman filter, and $\alpha$ is decided by the signal-to-noise ratio $q = Q/R$.

4. Worked example #

(a) Three steps of the local level filter. $Q = 1$, $R = 2$, start $\hat x_{0\mid 0} = 10$, $P_{0\mid 0} = 3$, observations $y_1 = 13$, $y_2 = 11$, $y_3 = 14$.

Step $t = 1$.

1. Predict: $\hat x_{1\mid 0} = 10$, $P_{1\mid 0} = 3 + 1 = 4$.
2. Innovation: $\nu_1 = 13 - 10 = 3$, with variance $S_1 = 4 + 2 = 6$.
3. Gain: $K_1 = 4/6 = 0.667$.
4. Update: $\hat x_{1\mid 1} = 10 + 0.667 \cdot 3 = 12$, $P_{1\mid 1} = (1 - 0.667) \cdot 4 = 1.333$.

Step $t = 2$.

1. Predict: $\hat x_{2\mid 1} = 12$, $P_{2\mid 1} = 1.333 + 1 = 2.333$.
2. Innovation: $\nu_2 = 11 - 12 = -1$, $S_2 = 2.333 + 2 = 4.333$.
3. Gain: $K_2 = 2.333/4.333 = 0.538$.
4. Update: $\hat x_{2\mid 2} = 12 - 0.538 = 11.462$, $P_{2\mid 2} = 0.462 \cdot 2.333 = 1.077$.

Step $t = 3$.

1. Predict: $\hat x_{3\mid 2} = 11.462$, $P_{3\mid 2} = 2.077$.
2. Innovation: $\nu_3 = 14 - 11.462 = 2.538$, $S_3 = 4.077$.
3. Gain: $K_3 = 2.077/4.077 = 0.509$.
4. Update: $\hat x_{3\mid 3} = 11.462 + 0.509 \cdot 2.538 = 12.755$, $P_{3\mid 3} = 1.019$.

Steady state. The gains $0.667, 0.538, 0.509$ approach $K_\infty$. Here $q = 1/2$, so

with $P = \tfrac12(1 + \sqrt{1 + 8}) = 2$. The filter becomes SES with $\alpha = 0.5$.

Forecast from $t = 3$. $\hat y_{3+h\mid 3} = 12.755$ for every $h$. The variance of $y_{3+h}$ is $P_{3\mid 3} + hQ + R = 1.019 + h + 2$: $4.019$, $5.019$, $6.019$ for $h = 1, 2, 3$. The 95 % interval for $h = 1$ is $12.755 \pm 1.96\sqrt{4.019} = 12.755 \pm 3.93$.

Log-likelihood. Each term is $-\tfrac12\big(\log(2\pi S_t) + \nu_t^2/S_t\big)$: $-2.565$, $-1.768$, $-2.412$, so $\log L = -6.744$.

Missing $y_2$. Skip the update at $t = 2$: $\hat x_{2\mid 2} = 12$, $P_{2\mid 2} = P_{2\mid 1} = 2.333$. Then $P_{3\mid 2} = 3.333$, $K_3 = 3.333/5.333 = 0.625$ (larger than before) and $\hat x_{3\mid 3} = 12 + 0.625 \cdot (14 - 12) = 13.25$.

(b) One step of constant-velocity tracking. $F = \begin{pmatrix}1&1\\0&1\end{pmatrix}$, $H = (1, 0)$, $Q = \operatorname{diag}(0, 1)$, $R = 1$, start $\hat x_{0\mid 0} = (0, 1)$ (position 0, velocity 1), $P_{0\mid 0} = I$. We observe $y_1 = 2.5$.

1. Predict the mean: $\hat x_{1\mid 0} = F\hat x_{0\mid 0} = (0 + 1, 1) = (1, 1)$.
2. Predict the covariance: $FP_{0\mid 0}F^\top = FF^\top = \begin{pmatrix}2&1\\1&1\end{pmatrix}$, and adding $Q$ gives $P_{1\mid 0} = \begin{pmatrix}2&1\\1&2\end{pmatrix}$. The off-diagonal 1 is new: position and velocity are now correlated.
3. Innovation: $\nu_1 = 2.5 - 1 = 1.5$; $S_1 = HP_{1\mid 0}H^\top + R = 2 + 1 = 3$.
4. Gain: $K_1 = P_{1\mid 0}H^\top/S_1 = (2, 1)^\top/3$, that is $(0.667, 0.333)^\top$.
5. Update: $\hat x_{1\mid 1} = (1, 1) + 1.5\,(0.667, 0.333) = (2, 1.5)$.
6. Covariance:

The target is further ahead than predicted, so the filter raises the velocity from 1 to 1.5, although the velocity was not measured. The position variance fell from 2 to $2/3$; the velocity variance fell only from 2 to $5/3$.

5. In code #

One function covers every model of this lesson. It treats a row of NaN as a missing observation.

import numpy as np
from scipy.optimize import minimize

def kalman_filter(y, F, H, Q, R, x0, P0):
    """Filtered means/covariances and log-likelihood. y: (n, k) array; NaN rows are missing."""
    x, P = np.asarray(x0, float), np.asarray(P0, float)
    xs, Ps, loglik = [], [], 0.0
    for yt in y:
        x, P = F @ x, F @ P @ F.T + Q                   # predict
        if not np.any(np.isnan(yt)):                     # update (skip if missing)
            nu = yt - H @ x                              # innovation
            S = H @ P @ H.T + R                          # its covariance
            K = P @ H.T @ np.linalg.inv(S)               # Kalman gain
            x, P = x + K @ nu, (np.eye(len(x)) - K @ H) @ P
            loglik -= 0.5 * (np.log(np.linalg.det(2 * np.pi * S)) + nu @ np.linalg.solve(S, nu))
        xs.append(x); Ps.append(P)
    return np.array(xs), np.array(Ps), loglik

# worked example (a): local level, Q = 1, R = 2, prior N(10, 3)
I = np.eye(1)
y = np.array([[13.], [11.], [14.]])
xs, Ps, ll = kalman_filter(y, I, I, 1 * I, 2 * I, [10.], 3 * I)
print(xs.ravel().round(4), Ps.ravel().round(4), round(ll, 4))
# [12.     11.4615 12.7547] [1.3333 1.0769 1.0189] -6.7442
y_miss = np.array([[13.], [np.nan], [14.]])
print(kalman_filter(y_miss, I, I, I, 2 * I, [10.], 3 * I)[0].ravel().round(4))
# [12.   12.   13.25]

# worked example (b): one step of constant-velocity tracking
F = np.array([[1., 1.], [0., 1.]]); H = np.array([[1., 0.]])
xs, Ps, _ = kalman_filter(np.array([[2.5]]), F, H, np.diag([0., 1.]), np.eye(1), [0., 1.], np.eye(2))
print(xs[0], Ps[0].round(4).tolist())
# [2.  1.5] [[0.6667, 0.3333], [0.3333, 1.6667]]

# Nile: maximum likelihood for (Q, R); NILE holds the 100 values (they are in the notebook)
def neg_loglik(theta):
    Q, R = np.exp(theta)                       # log scale keeps the variances positive
    # exact diffuse start: x_{1|1} = y_1, P_{1|1} = R; the likelihood starts at t = 2
    return -kalman_filter(NILE[1:, None], I, I, Q * I, R * I, [NILE[0]], R * I)[2]
res = minimize(neg_loglik, x0=np.log([1000., 10000.]), method="Nelder-Mead")
Q_hat, R_hat = np.exp(res.x)
q = Q_hat / R_hat
print(round(Q_hat), round(R_hat), round(-res.fun, 2))         # 1469 15099 -632.55
print(round(q, 4), round((np.sqrt(q**2 + 4*q) - q) / 2, 3))   # 0.0973 0.267

For real work, statsmodels has the same model as sm.tsa.UnobservedComponents(y, "local level") and general models with sm.tsa.statespace.MLEModel. The notebook compares the results in an optional cell.

6. Where it shows up #

• Classical forecasting software. statsmodels fits ARIMA, SARIMAX, unobserved-components and dynamic factor models through the Kalman filter; this is how it handles missing values and computes the exact likelihood. The R package KFAS and the structural models of Harvey (1989) work the same way.
• Structural time-series models. Write a series as level + trend + season + regression, each a block of the state, and read the smoothed components. Time-varying regression coefficients ($y_t = z_t^\top\beta_t + v_t$ with $\beta_t$ a random walk) are a state-space model with $H_t = z_t^\top$; economists use them for changing relations and nowcasting.
• Navigation and robotics. GPS + inertial sensor fusion, target tracking with radar, and the EKF in SLAM.
• Weather and climate. Data assimilation combines a physical model forecast with millions of observations; the ensemble Kalman filter and 4D-Var are the large-scale versions of the update step (see the ml-for-science track).
• Machine learning. Probabilistic state-space models with a neural network inside (deep Kalman filters, DeepState) use the same predict–update structure. S4 and Mamba layers use the linear state equation without noise (deep-learning lesson 16). A Gaussian process with a Matérn kernel in time can be written as a state-space model, so a Kalman filter computes it in $O(n)$ instead of $O(n^3)$ (Hartikainen and Särkkä, 2010).

7. Common mistakes #

• Mixing up the state and the observation. $x_t$ is hidden; $y_t$ is measured. The filtered level $\hat x_{t\mid t}$ is not the data, and the interval for $y_{t+h}$ must add $R$ to the state variance.
• Mixing up the covariances. $P_{t\mid t-1}$ (before $y_t$), $P_{t\mid t}$ (after $y_t$) and $S_t$ (of the observation) are three different things. The gain uses $P_{t\mid t-1}$, not $P_{t\mid t}$.
• Thinking the gain is a free parameter. $K_t$ follows from $F, H, Q, R$ and $P_{0\mid 0}$. To get a smoother or a more reactive filter, change $Q$ or $R$ (or better, estimate them).
• Too small a $Q$. With $Q$ near 0 the filter "trusts the model" and becomes slow; it can lag behind a real change for a long time. The innovations then have long runs of one sign: check them.
• Using the smoother for forecasts or in a backtest. It uses future data. Only the filter is allowed at a forecast origin.
• Forgetting the start. A small $P_{0\mid 0}$ with a wrong $\hat x_{0\mid 0}$ makes the first estimates wrong for a long time, and the first likelihood terms can dominate the fit. Use a diffuse start when you know nothing.

8. Check yourself #

1. Write the model $z_t = \phi z_{t-1} + \varepsilon_t$, $y_t = z_t + v_t$ (an AR(1) measured with noise) in state-space form.
   Answer

   State $x_t = z_t$, $F = \phi$, $H = 1$, $Q = \sigma_\varepsilon^2$, $R = \sigma_v^2$. With $|\phi| < 1$ a natural start is the stationary distribution: $\hat x_{0\mid 0} = 0$, $P_{0\mid 0} = \sigma_\varepsilon^2/(1 - \phi^2)$ (lesson 02).
2. In the local level model set $Q = 0$ and use a diffuse start. What does the filter compute?
   Answer

   The running mean. With $Q = 0$ the level is a constant $\mu$. After the diffuse start, $\hat x_{1\mid 1} = y_1$ and $P_{1\mid 1} = R$. By induction, if $P_{t-1\mid t-1} = R/(t-1)$, then $K_t = \frac{R/(t-1)}{R/(t-1) + R} = 1/t$ and $\hat x_{t\mid t} = \hat x_{t-1\mid t-1} + \frac1t(y_t - \hat x_{t-1\mid t-1})$, the update of the mean of $y_1, \dots, y_t$, with $P_{t\mid t} = R/t$. The gain goes to 0, which matches $K_\infty = 0$ for $q = 0$.
3. In worked example (b), why did the velocity estimate change although $H = (1, 0)$ measures only the position?
   Answer

   The predict step made position and velocity correlated: $P_{1\mid 0}$ has off-diagonal entry 1. The gain is $P_{1\mid 0}H^\top/S_1$, the first column of $P_{1\mid 0}$ over $S_1$, so its second entry is $1/3 \ne 0$. A position surprise of $+1.5$ therefore raises the velocity by $1.5/3 = 0.5$.
4. Why does the likelihood use the innovations $\nu_t$ and not the observations $y_t$ themselves?
   Answer

   The $y_t$ are dependent, so their joint density is not a product of marginals. The chain rule writes it as a product of one-step conditional densities $p(y_t \mid y_{1:t-1})$, and each one is the Gaussian density of the innovation $\nu_t$ with covariance $S_t$. Under the model the innovations are independent, so the log-likelihood is a simple sum.
5. For the Nile data the filter at 1970 gives $\hat x_{1970\mid 1970} = 798.4$ and $P_{1970\mid 1970} = 4\,032$, with $\hat Q = 1\,469$ and $\hat R = 15\,099$. Give the 95 % prediction interval for the flow in 1971.
   Answer

   The variance is $P_{1970\mid 1970} + Q + R = 4\,032 + 1\,469 + 15\,099 = 20\,600$, so the standard deviation is $143.5$. The interval is $798 \pm 1.96 \cdot 143.5 = 798 \pm 281$, about $[517, 1080]$. Most of the width comes from $R$, the year-to-year noise, not from the uncertainty about the level.
6. At which time do the filter and the smoother give the same estimate, and why?
   Answer

   At the last time $t = n$. There, both condition on the same data $y_1, \dots, y_n$, so $\hat x_{n\mid n}$ is both the filtered and the smoothed estimate. The RTS recursion starts from it and goes backward.

Recall cardsAnswer in your head first, then reveal. Add the cards to your review deck to see them again just before you would forget.

1. Write the two equations of a linear Gaussian state-space model.

   Show answer

2. What does the predict step of the Kalman filter compute?

   Show answer

3. Give the Kalman gain in the vector case.

   Show answer

4. What does the update step compute?

   Show answer

5. Prior $\mathcal N(10, 4)$, observation $y = 14$ with $R = 4$. What is the posterior?

   Show answer

6. What is the steady-state gain of the local level model?

   Show answer

7. How are SES and the Kalman filter related?

   Show answer

8. How does the Kalman filter handle a missing observation?

   Show answer

9. State the prediction-error decomposition of the log-likelihood.

   Show answer

10. Why may only the filter, not the smoother, be used at a forecast origin?

    Show answer

Add these cards to my reviewOpen the review →

9. Interview angle #

• "Explain the Kalman filter in one minute." — A hidden state follows a linear model with noise and is measured with noise. Predict: push the estimate and its covariance through the model. Update: correct with the innovation times the gain, which weighs prediction uncertainty against sensor noise. For linear Gaussian models it is the exact posterior.
• "What is the intuition for the Kalman gain?" — It is the relative trust in the measurement: $P/(P + R)$ in one dimension. Precisions add. Large $R$ → small gain; large prediction uncertainty → gain near 1.
• "How do you choose $Q$ and $R$?" — Estimate them by maximum likelihood with the prediction-error decomposition, then check that the standardised innovations look like white noise with variance 1. Hand tuning is a last resort.
• "What is the link between SES and the Kalman filter?" — SES is the steady-state filter of the local level model, with $\alpha$ fixed by $q = Q/R$. ARIMA(0,1,1) with $\theta = \alpha - 1$ gives the same forecasts.
• "Filter vs smoother vs forecast?" — $\hat x_{t\mid t}$, $\hat x_{t\mid n}$ and $\hat x_{t+h\mid t}$. Only the filter and the forecast are allowed in real-time use and in backtests.
• "What do you do for a nonlinear model?" — EKF (linearise), UKF (sigma points) or a particle filter (samples); for very large states the ensemble Kalman filter.

10. Practice #

The notebook for this lesson runs in Colab. It has five exercises, each with a check cell and a folded solution:

1. Write kalman_filter(y, F, H, Q, R, x0, P0) that returns filtered means, covariances and the log-likelihood and treats np.nan as missing; check it on worked example (a) and against Bayes' rule computed on a grid.
2. Track a constant-velocity target from position-only observations and compare the RMSE of the filtered position with the RMSE of the raw observations.
3. Compute the steady-state gain of the local level model and show that SES with $\alpha = K_\infty$ gives the same values as the filter after the start-up.
4. Fit the local level model to the Nile data: maximise the log-likelihood over $(\log Q, \log R)$ with scipy.optimize.minimize and compare with $\hat R \approx 15\,099$, $\hat Q \approx 1\,469$.
5. Write rts_smoother and check that the smoothed variances are never larger than the filtered ones.

11. Watch and read #

Visually Explained: Kalman FiltersVisually Explained · 11 min

The predict and update steps as pictures of Gaussians; matches section 2.

Special Topics - The Kalman Filter (1 of 55) What is a Kalman Filter?Michel van Biezen · 6 min

The first video of a long series with small numerical examples; useful if you want more worked steps.

02417 Lecture 11 part A: Introduction to state space modelsLasse Engbo Christiansen · 8 min

State-space models from the time-series view (DTU course 02417), with the same two equations.

The Kalman Filter [Control Bootcamp]Steve Brunton · 6 min

The control-engineering view: the filter as an optimal observer; uses continuous time.

• Durbin and Koopman, Time Series Analysis by State Space Methods, 2nd ed., Oxford University Press, 2012 (publisher page): ch. 2 is the local level model with the Nile data; ch. 3 puts ARIMA and structural models in state-space form.
• Simo Särkkä, Bayesian Filtering and Smoothing, Cambridge University Press 2013, free PDF from the author: ch. 4 (Kalman filter as Bayesian filtering) and ch. 8 (RTS smoother).
• Hyndman and Athanasopoulos, Forecasting: Principles and Practice (3rd ed.), §8.5 Innovations state space models for exponential smoothing.
• Greg Welch and Gary Bishop, An Introduction to the Kalman Filter, UNC technical report TR 95-041: a short, practical introduction with the discrete filter equations.
• Shumway and Stoffer, Time Series Analysis and Its Applications, 4th ed., Springer 2017, ch. 6 (state-space models), with R code on the authors' GitHub page.
• statsmodels documentation: UnobservedComponents and state-space methods.