The Kalman Filter

"Every instant I am handed two stories: the one my model tells about where I should be, and the one the sensor swears it just measured. I do not pick a winner. I weigh each by how much I trust it, and the weight has a name."

A Kalman Gain Balancing Faith in Model and Measurement

The previous section, Section 7.1, introduced the linear-Gaussian state-space model: a hidden state $\mathbf{z}_t$ that evolves by a linear transition $\mathbf{z}_t = F\mathbf{z}_{t-1} + \mathbf{w}_t$ with process noise $\mathbf{w}_t \sim \mathcal{N}(\mathbf{0}, Q)$, observed only through a linear measurement $\mathbf{x}_t = H\mathbf{z}_t + \mathbf{v}_t$ with measurement noise $\mathbf{v}_t \sim \mathcal{N}(\mathbf{0}, R)$. That section posed the filtering question but did not solve it. The Kalman filter is the solution, and it is the algorithmic heart of Chapter 7. We assume the reader is comfortable with the multivariate Gaussian, with conditional expectation, and with the companion-form linear recursion that Section 6.1 built for the vector autoregression, because the VAR's observed state recursion $\mathbf{z}_t = F\mathbf{z}_{t-1} + \tilde{\boldsymbol{\varepsilon}}_t$ is exactly the transition the Kalman filter now handles when the state is only partially observed. Throughout we use the state-space symbols of the unified notation table in Appendix A: $\mathbf{z}_t$ for the latent state, $\mathbf{x}_t$ for the observation, $F, H, Q, R$ for the four system matrices, and $P_t$ for the state-estimate covariance. We keep state $\mathbf{z}_t$ and observation $\mathbf{x}_t$ throughout Chapter 7, except Section 7.6, which renames the state $\mathbf{h}_t$ to match the recurrent-network literature.

1. The Filtering Problem and the Recursive Bayesian Idea Beginner

Filtering is the task of computing the posterior distribution of the current hidden state given every measurement received up to now, $p(\mathbf{z}_t \mid \mathbf{x}_{1:t})$. The word "filtering" is deliberate: we want the best estimate at the present moment using the present and the past, as opposed to smoothing (estimating a past state using future measurements too, subsection five) or prediction (estimating a future state before its measurement arrives). The naive approach, recomputing the posterior from scratch over the whole history at every new measurement, is hopeless online because the cost grows without bound as the series lengthens. The escape is recursion: maintain a compact summary of everything seen so far and update it incrementally as each measurement lands.

The recursive Bayesian filter makes this precise in two alternating moves that hold for any state-space model, linear or not. Suppose we have the posterior from the previous step, $p(\mathbf{z}_{t-1} \mid \mathbf{x}_{1:t-1})$. The predict step pushes it through the dynamics using the Chapman-Kolmogorov equation, marginalizing out the previous state to obtain the prior for the current state before seeing the new measurement,

The update step then folds in the new measurement $\mathbf{x}_t$ by Bayes' rule, multiplying the predicted prior by the measurement likelihood and renormalizing,

This predict-update cycle is the universal skeleton of recursive estimation. For a general nonlinear or non-Gaussian model the integral and the normalization have no closed form, which is why the particle filters and nonlinear variants of Section 7.3 exist. The miracle of the linear-Gaussian case is that both operations preserve Gaussianity, so the entire belief is forever summarized by a mean and a covariance.

We adopt the standard double-subscript notation that keeps predict and update legible. Write $\hat{\mathbf{z}}_{t \mid t-1}$ and $P_{t \mid t-1}$ for the mean and covariance of the predicted (prior) belief, conditioned on measurements through time $t-1$, and write $\hat{\mathbf{z}}_{t \mid t}$ and $P_{t \mid t}$ for the updated (posterior) belief, conditioned on measurements through time $t$. The filter is then a loop: from the previous posterior $(\hat{\mathbf{z}}_{t-1 \mid t-1}, P_{t-1 \mid t-1})$, predict to get $(\hat{\mathbf{z}}_{t \mid t-1}, P_{t \mid t-1})$, then update with $\mathbf{x}_t$ to get $(\hat{\mathbf{z}}_{t \mid t}, P_{t \mid t})$, and repeat. Figure 7.2.1 shows the cycle.

2. The Predict Step: Propagating the Prior Intermediate

Start from the previous posterior, a Gaussian with mean $\hat{\mathbf{z}}_{t-1 \mid t-1}$ and covariance $P_{t-1 \mid t-1}$. The dynamics say $\mathbf{z}_t = F\mathbf{z}_{t-1} + \mathbf{w}_t$ with $\mathbf{w}_t \sim \mathcal{N}(\mathbf{0}, Q)$ independent of the state. We push the Gaussian through this linear-plus-noise map and read off the new mean and covariance. The mean of a linear transform is the linear transform of the mean, and the additive noise has mean zero, so the predicted state mean is simply the deterministic part of the dynamics evaluated at the previous estimate,

(If the model carries a control input $\mathbf{u}_t$ entering as $B\mathbf{u}_t$, it is added here; we suppress it for clarity, since the autonomous case carries every idea.) The covariance requires the standard transformation rule: a linear map $F$ acting on a random vector with covariance $P$ produces a vector with covariance $FPF^\top$, and because the process noise is independent of the previous state its covariance $Q$ adds on. Hence the predicted covariance is

Read this equation as the precise statement that prediction increases uncertainty. The term $FP_{t-1\mid t-1}F^\top$ propagates the old uncertainty through the dynamics (it may grow or shrink depending on whether $F$ stretches or contracts), and then $Q$ adds the irreducible new uncertainty injected because the model is an imperfect description of reality. With no measurement to anchor it, a pure predict loop would let $P$ grow until the estimate is worthless, which is exactly the situation in dead reckoning between GPS fixes: a vehicle that integrates its wheel-speed and heading without external correction watches its position covariance balloon until a fresh fix collapses it again. The update step is what supplies that collapse.

3. The Update Step: Innovation, Gain, and the Trust Trade-off Intermediate

Now the measurement $\mathbf{x}_t$ arrives, and we fold it in by Bayes' rule. Because both the predicted prior and the measurement likelihood are Gaussian, the posterior is Gaussian and its mean and covariance follow in closed form. The derivation organizes around three quantities. First, what did the model expect to measure, and how wrong was it? The predicted measurement is $H\hat{\mathbf{z}}_{t \mid t-1}$, and the innovation (or measurement residual) is the surprise, the gap between what arrived and what was expected,

Second, how surprised should we have been? The innovation has its own covariance, called the innovation covariance, which combines the predicted state uncertainty mapped into measurement space with the measurement noise,

The two pieces of $S_t$ have clean meanings: $HP_{t\mid t-1}H^\top$ is how uncertain the model is about what it should measure (state uncertainty seen through $H$), and $R$ is how noisy the sensor itself is. Third, and centrally, the Kalman gain $K_t$ converts a measurement-space surprise into a state-space correction,

With these three objects the update writes itself. The posterior mean is the prior mean nudged along the gain in proportion to the innovation, and the posterior covariance is the prior covariance reduced by the information the measurement supplied,

These three update lines are not assumed; they fall out of completing the square in the product of the two Gaussians, or equivalently from the standard formula for the conditional distribution of one jointly-Gaussian block given another. The joint distribution of the state and the next measurement, both before $\mathbf{x}_t$ is seen, is Gaussian with cross-covariance $\operatorname{Cov}(\mathbf{z}_t, \mathbf{x}_t) = P_{t\mid t-1}H^\top$ and measurement marginal covariance $S_t = HP_{t\mid t-1}H^\top + R$. The Gaussian conditioning formula then states that the mean of $\mathbf{z}_t$ given $\mathbf{x}_t$ is the prior mean plus (cross-covariance)(marginal-covariance)$^{-1}$ times the residual, which is exactly $\hat{\mathbf{z}}_{t\mid t-1} + P_{t\mid t-1}H^\top S_t^{-1}\tilde{\mathbf{x}}_t$, and the conditional covariance shrinks by the same factor, giving $P_{t\mid t-1} - P_{t\mid t-1}H^\top S_t^{-1}H P_{t\mid t-1} = (I - K_t H)P_{t\mid t-1}$. So the Kalman gain $K_t = P_{t\mid t-1}H^\top S_t^{-1}$ is precisely the regression coefficient of the hidden state on the measurement residual; the update is a single step of Gaussian conditioning, no more and no less.

The structure of the gain is the entire epigraph rendered as algebra. Read $K_t = P_{t\mid t-1}H^\top S_t^{-1}$ as a ratio: the numerator $P_{t\mid t-1}H^\top$ is the model's uncertainty (mapped into measurement space), and the denominator $S_t$ is the total uncertainty (model plus measurement). When the model is very uncertain relative to the sensor, the gain is large and the update leans hard on the measurement, snapping the estimate toward $\mathbf{x}_t$. When the sensor is very noisy relative to the model, the gain is small and the update barely moves, keeping faith with the prediction. The gain is the formal, optimal answer to "how much should I trust this new measurement?", and it is recomputed every step from the current balance of uncertainties. This is why the persona of the epigraph refuses to pick a winner: the correct weight is never zero or one but a precise blend that the gain delivers.

A worked single cycle makes the abstraction concrete. The numeric example below walks one full predict-update on a scalar local-level model, the simplest non-trivial state-space model, where the state is an unobserved level that drifts and the measurement is that level plus noise.

4. Optimality, the Steady-State Filter, and the Log-Likelihood Advanced

The Kalman filter is not merely a reasonable estimator; for the linear-Gaussian model it is the optimal one, in a strong and precise sense. Among all estimators, linear or nonlinear, that use the measurements $\mathbf{x}_{1:t}$, the Kalman filter's posterior mean $\hat{\mathbf{z}}_{t\mid t}$ is the minimum-mean-squared-error estimate of the state, because it is exactly the conditional expectation $\mathbb{E}[\mathbf{z}_t \mid \mathbf{x}_{1:t}]$ and the conditional expectation is the MMSE estimator for any distribution. It is also unbiased, $\mathbb{E}[\hat{\mathbf{z}}_{t\mid t} - \mathbf{z}_t] = \mathbf{0}$, and $P_{t\mid t}$ is its true error covariance. If we drop the Gaussian assumption but keep the means and covariances, a second, weaker optimality survives: the Kalman filter is still the best linear unbiased estimator (the minimum-variance estimator among all linear functions of the data), so it remains the right answer even when the noise is merely white with known covariance rather than truly Gaussian. This dual guarantee, exactly optimal under Gaussianity and best-linear otherwise, is why the filter is the default estimator across engineering and statistics.

A striking practical fact is that for a time-invariant model (constant $F, H, Q, R$) the covariance recursion forgets its starting point. The covariance update $P_{t\mid t}$ does not depend on the measurements at all, only on the system matrices, so it can be run offline. Iterating predict-then-update on the covariance alone, $P^- = FPF^\top + Q$ followed by $P = (I - KH)P^-$, drives $P$ to a fixed point $P_\infty$ that satisfies the discrete-time algebraic Riccati equation. At that fixed point the gain is constant, $K_\infty = P_\infty^- H^\top S_\infty^{-1}$, and the filter becomes a fixed linear recursion $\hat{\mathbf{z}}_{t\mid t} = (I - K_\infty H)F\hat{\mathbf{z}}_{t-1\mid t-1} + K_\infty\mathbf{x}_t$ with no per-step matrix inverse. This steady-state Kalman filter is what many embedded systems ship: precompute $K_\infty$ once, then run a cheap constant-gain update forever. The convergence of $P$ to $P_\infty$ is what the worked example will display as "the gain settling."

The filter quietly hands us a second prize that turns it into an estimation engine, not just a tracker. Running the filter produces, at each step, the innovation $\tilde{\mathbf{x}}_t$ and its covariance $S_t$, and under the model the innovations are independent zero-mean Gaussians with covariance $S_t$. The log-likelihood of the entire observed series therefore factors into a sum of one-step innovation terms, the prediction-error decomposition,

where $d$ is the measurement dimension and $\boldsymbol{\theta}$ collects the unknown entries of $F, H, Q, R$. This is the key that unlocks parameter estimation: to fit the unknown system matrices (most commonly the noise variances $Q$ and $R$), run the filter to get the innovations, evaluate this log-likelihood, and maximize it numerically. The same information criteria of Section 5.7 then apply directly to compare competing state-space specifications: AIC and BIC built on this exact log-likelihood let us choose, for instance, between a local-level and a local-linear-trend model, just as they chose ARIMA orders. The filter is thus both the inference algorithm and the likelihood evaluator, which is why a single state-space fit yields tracking, forecasting, and model selection in one pass.

5. Smoothing: Refining the Past with the Future Advanced

The filter estimates $\mathbf{z}_t$ using measurements up to time $t$ only. Often we are analyzing a recorded series after the fact and want the best estimate of every past state given the entire dataset, $p(\mathbf{z}_t \mid \mathbf{x}_{1:T})$ with $T > t$, which is the smoothing problem. The Rauch-Tung-Striebel (RTS) smoother solves it with an elegant backward pass that reuses the filter's stored quantities. Run the forward Kalman filter once, storing every filtered $(\hat{\mathbf{z}}_{t\mid t}, P_{t\mid t})$ and predicted $(\hat{\mathbf{z}}_{t+1\mid t}, P_{t+1\mid t})$, then sweep backward from $t = T-1$ to $1$, correcting each filtered estimate using the already-smoothed future state through a smoother gain $C_t = P_{t\mid t}F^\top P_{t+1\mid t}^{-1}$ and the recursion $\hat{\mathbf{z}}_{t\mid T} = \hat{\mathbf{z}}_{t\mid t} + C_t(\hat{\mathbf{z}}_{t+1\mid T} - \hat{\mathbf{z}}_{t+1\mid t})$, with a matching covariance update that always satisfies $P_{t\mid T} \preceq P_{t\mid t}$. Intuitively, future measurements carry information about where the state was, so the smoothed estimate is never less certain than the filtered one and usually much sharper, especially at the start of the series where the filter had little history. Smoothing is the right tool for offline analysis, for the expectation step when fitting state-space parameters by EM, and for imputing the irregularly-sampled clinical series that thread through this book's healthcare examples (Section 7.1); for real-time tracking, where the future is not yet available, the filter is all we have.

6. Worked Example: A Kalman Filter From Scratch, Then the Library Pair Advanced

We now build the entire filter from scratch and run it on a concrete tracking problem, then reproduce it in a few lines with filterpy and statsmodels. The system is a one-dimensional constant-velocity tracker, the classic introductory Kalman target: an object moves along a line at a roughly constant velocity that drifts slowly (a small random acceleration is the process noise), and a sensor reports its noisy position once per time step. The hidden state is two-dimensional, position and velocity, $\mathbf{z}_t = (p_t, v_t)^\top$, even though we measure only position, so the filter must infer the unobserved velocity from the position stream alone, which is the entire point of having a latent state. Code 7.2.1 sets up the model matrices and simulates the ground truth and the noisy measurements.

import numpy as np

rng = np.random.default_rng(7)

# State z = (position, velocity). Constant-velocity dynamics with dt = 1.
dt = 1.0
F = np.array([[1.0, dt],      # position += velocity * dt
              [0.0, 1.0]])    # velocity persists
H = np.array([[1.0, 0.0]])    # we measure position only, not velocity

# Process noise: a small random acceleration perturbs the velocity (and,
# through it, the position). q scales the continuous-acceleration variance.
q = 0.004
Q = q * np.array([[dt**3 / 3, dt**2 / 2],
                  [dt**2 / 2, dt]])
R = np.array([[0.25]])        # sensor position-noise variance (std = 0.5)

# Simulate T steps of ground truth and noisy position measurements.
T = 80
z_true = np.zeros((T, 2)); z_true[0] = [0.0, 1.0]   # start at origin, speed 1
y = np.zeros((T, 1))
chol_Q = np.linalg.cholesky(Q)
for t in range(1, T):
    z_true[t] = F @ z_true[t-1] + chol_Q @ rng.standard_normal(2)
y[:, 0] = z_true[:, 0] + np.sqrt(R[0, 0]) * rng.standard_normal(T)
print("true final state (pos, vel):", z_true[-1].round(3))
print("first 4 noisy measurements :", y[:4, 0].round(3))

true final state (pos, vel): [79.642  0.94 ]
first 4 noisy measurements : [ 0.371  1.245  2.4    2.879]

Code 7.2.2 implements the filter itself: a direct transcription of the predict equations of subsection two and the update equations of subsection three, looping over the measurement stream and recording the gain at each step so we can watch it settle toward the steady-state value of subsection four.

def kalman_filter(y, F, H, Q, R, z0, P0):
    """Run a linear-Gaussian Kalman filter over a (T, d_y) measurement array.
    Returns filtered state means, covariances, and the per-step Kalman gain."""
    T = y.shape[0]
    n = F.shape[0]
    z = z0.copy(); P = P0.copy()
    I = np.eye(n)
    z_filt = np.zeros((T, n))
    P_filt = np.zeros((T, n, n))
    gains  = np.zeros((T, n, H.shape[0]))
    for t in range(T):
        # --- Predict: push the belief through the dynamics, grow uncertainty.
        z_pred = F @ z                       # z_{t|t-1} = F z_{t-1|t-1}
        P_pred = F @ P @ F.T + Q             # P_{t|t-1} = F P F^T + Q
        # --- Update: fold in the measurement, shrink uncertainty.
        innov = y[t] - H @ z_pred            # innovation (the surprise)
        S = H @ P_pred @ H.T + R             # innovation covariance
        K = P_pred @ H.T @ np.linalg.inv(S)  # Kalman gain = trust weight
        z = z_pred + K @ innov               # posterior mean
        P = (I - K @ H) @ P_pred             # posterior covariance
        z_filt[t] = z; P_filt[t] = P; gains[t] = K
    return z_filt, P_filt, gains

z0 = np.array([0.0, 0.0])          # ignorant prior: guess zero position/velocity
P0 = np.eye(2) * 10.0              # large prior variance: "I really don't know"
z_filt, P_filt, gains = kalman_filter(y, F, H, Q, R, z0, P0)

pos_rmse_raw  = np.sqrt(np.mean((y[:, 0]      - z_true[:, 0])**2))
pos_rmse_filt = np.sqrt(np.mean((z_filt[:, 0] - z_true[:, 0])**2))
print(f"position RMSE, raw measurements : {pos_rmse_raw:.3f}")
print(f"position RMSE, Kalman-filtered  : {pos_rmse_filt:.3f}")
print("recovered final velocity (never measured):", round(z_filt[-1, 1], 3))
print("position-gain at steps 0, 5, 40 :",
      gains[[0, 5, 40], 0, 0].round(3))

position RMSE, raw measurements : 0.503
position RMSE, Kalman-filtered  : 0.235
recovered final velocity (never measured): 0.967
position-gain at steps 0, 5, 40 : [0.989 0.46  0.356]

The numbers tell the story but a picture makes it vivid. Figure 7.2.2 plots, on simulated-style coordinates, the noisy measurements scattering around the true trajectory and the filtered estimate threading smoothly between them, alongside the position gain decaying from one to its steady-state value.

That was the transparent machinery, every matrix operation visible. In production you would not hand-roll the predict-update loop, the covariance bookkeeping, and the numerical-stability safeguards; mature libraries package all of it. Code 7.2.3 is the library pair, the identical filter through filterpy and, for the state-space-model framing with built-in maximum-likelihood parameter estimation, statsmodels.

# --- Shortcut A: filterpy, a direct object wrapper of the same recursion ---
from filterpy.kalman import KalmanFilter

kf = KalmanFilter(dim_x=2, dim_z=1)
kf.F, kf.H, kf.Q, kf.R = F, H, Q, R          # same matrices as Code 7.2.1
kf.x = np.array([0.0, 0.0]); kf.P = np.eye(2) * 10.0
xs = []
for t in range(T):
    kf.predict()                              # one call = our predict block
    kf.update(y[t])                           # one call = our update block
    xs.append(kf.x.copy())
xs = np.array(xs)
print("filterpy final state:", xs[-1].round(3))   # matches z_filt[-1]

# --- Shortcut B: statsmodels, which ALSO estimates Q and R by MLE ---
import statsmodels.api as sm
# UnobservedComponents fits a local-linear-trend (position + velocity) model
# and maximizes the prediction-error log-likelihood of subsection 4 for us.
mod = sm.tsa.UnobservedComponents(y[:, 0], level="local linear trend")
res = mod.fit(disp=False)                      # MLE of the noise variances
print("statsmodels smoothed final level:", round(res.states.smoothed[-1, 0], 3))
print("log-likelihood (prediction-error decomposition):", round(res.llf, 2))

filterpy final state: [79.51  0.967]
statsmodels smoothed final level: 79.583
log-likelihood (prediction-error decomposition): -71.34

Step back and see what the section assembled. We started from the recursive Bayesian filter, the universal predict-update skeleton, and used the linear-Gaussian conjugacy to collapse it to the propagation of a mean and a covariance. The predict step pushed the belief through the dynamics and grew its uncertainty by $Q$; the update step folded in the measurement through the Kalman gain, the precise precision-weighted blend of model and measurement that the epigraph's persona embodies. We established the filter's minimum-variance optimality, the steady-state constant-gain regime, and the prediction-error log-likelihood that turns the filter into a parameter-estimation and model-selection engine, then sketched the RTS smoother that sharpens the past with the future. The from-scratch implementation tracked a moving target through noise, recovered a velocity it never measured, and watched its gain settle, before filterpy and statsmodels reproduced it in a handful of lines and added the maximum-likelihood parameter fitting we had assumed away. The filter's one structural limitation is its linearity assumption: the predict and update steps both required linear $F$ and $H$ so that Gaussians stay Gaussian. The next section, Section 7.3, removes that assumption, building the extended and unscented Kalman filters that handle nonlinear dynamics and measurements while preserving the predict-update structure we have just made second nature.