Skip to content

Latest commit

 

History

History
1174 lines (907 loc) · 44.5 KB

File metadata and controls

1174 lines (907 loc) · 44.5 KB

Derivation of the Kalman Filter

This document derives the discrete-time linear Kalman filter from first principles. The derivation establishes what the filter is, why it works, and what it assumes. The reference implementation in kalman_filter.py is a direct expression of the mathematics here.

1. Problem Statement

1.1 State estimation as a problem class

Many systems of engineering interest share a common structure. An underlying state evolves over time according to known dynamics. That state is observed through measurements that are noisy, indirect, or both. The state itself is hidden. What is available are the measurements and a model of how the state and measurements relate.

The state estimation problem is to recover the best estimate of the hidden state from the sequence of measurements observed so far. "Best" must be defined precisely. The answer depends on what assumptions hold about the system and the noise.

Examples that fit this structure:

  • Tracking a moving object whose position and velocity must be inferred from noisy position measurements
  • Estimating the orientation of a spacecraft from gyroscope and star tracker readings
  • Recovering a signal from a noisy sensor whose dynamics are known
  • Estimating the parameters of a slowly time-varying system from measurements of its outputs

In each case, the state is what we want to know. The measurements are what we have. The dynamics and measurement model are what we assume.

1.2 The discrete-time linear-Gaussian setting

The Kalman filter solves this problem under three specific assumptions.

Time is discrete. The state evolves in steps: $k = 0, 1, 2, \ldots$ Measurements arrive at the same discrete times. Continuous-time formulations exist. This document does not cover them.

The dynamics and measurement model are linear. The state at time $k$ is a linear function of the state at time $k-1$ plus noise. The measurement at time $k$ is a linear function of the state at time $k$ plus noise. Nonlinear extensions exist, including the Extended Kalman Filter and the Unscented Kalman Filter. This document does not cover them.

The noise is Gaussian. Both the process noise and the measurement noise are Gaussian with known covariance. The process noise represents uncertainty in the dynamics. The measurement noise represents uncertainty in observations.

Under these three assumptions, the optimal state estimator has a closed-form recursive structure. That structure is the Kalman filter.

1.3 What we are solving and what we are assuming

The problem this document derives a solution to:

Given a sequence of measurements $z_1, z_2, \ldots, z_k$ and a model of how the state evolves and how measurements relate to the state, compute the conditional mean and covariance of the state at time $k$ given all measurements through time $k$.

The conditional mean is the best estimate in a sense made precise in Section 6. The conditional covariance quantifies the uncertainty remaining in that estimate.

The assumptions this derivation depends on:

  • The state evolves according to a linear process model with additive Gaussian process noise (Section 2.1)
  • Measurements are a linear function of the state with additive Gaussian measurement noise (Section 2.2)
  • The process noise and measurement noise are mutually independent, white, and have known covariance matrices (Section 2.3)
  • The initial state is Gaussian with known mean and covariance

These assumptions are restrictive. Many real systems violate them. The Kalman filter remains useful in those cases as an approximation and as the foundation for more general filters. Those extensions are out of scope here. The derivation that follows holds exactly under the assumptions above.

2. The State-Space Model

2.1 The process model

The state of the system at time $k$ is represented by a vector $\mathbf{x}_k \in \mathbb{R}^n$. The state evolves from one time step to the next according to:

$$ \mathbf{x}_k = \mathbf{F}_k \mathbf{x}_{k-1} + \mathbf{w}_k \tag{2.1} $$

where $\mathbf{F}_k \in \mathbb{R}^{n \times n}$ is the state transition matrix and $\mathbf{w}_k \in \mathbb{R}^n$ is the process noise.

The state transition matrix $\mathbf{F}_k$ encodes the dynamics of the system. It specifies how the state at time $k-1$ propagates to time $k$ in the absence of noise. The subscript $k$ allows $\mathbf{F}$ to vary with time. In many systems $\mathbf{F}$ is constant, in which case the subscript is dropped.

The process noise $\mathbf{w}_k$ accounts for everything the model does not capture exactly. Unmodeled forces, discretization error, and stochastic variation in the dynamics are all absorbed into $\mathbf{w}_k$. The size of $\mathbf{w}_k$ encodes how much the model is trusted at each time step.

2.2 The measurement model

The state $\mathbf{x}_k$ is not observed directly. What is observed at time $k$ is a measurement $\mathbf{z}_k \in \mathbb{R}^m$ related to the state by:

$$ \mathbf{z}_k = \mathbf{H}_k \mathbf{x}_k + \mathbf{v}_k \tag{2.2} $$

where $\mathbf{H}_k \in \mathbb{R}^{m \times n}$ is the measurement matrix and $\mathbf{v}_k \in \mathbb{R}^m$ is the measurement noise.

The measurement matrix $\mathbf{H}_k$ specifies what aspects of the state are observed and how. A position-only sensor observing a state of position and velocity has $\mathbf{H} = [1, 0]$. A sensor that observes a linear combination of state components has rows of $\mathbf{H}$ corresponding to that combination. The measurement need not be a one-to-one function of the state; $m < n$ is common and expected.

The measurement noise $\mathbf{v}_k$ accounts for sensor error, quantization, and any other source of corruption between the true state and the recorded measurement. Like the process noise, its size encodes how much the measurement is trusted.

2.3 The Gaussian noise assumptions

The process noise and measurement noise are assumed to be Gaussian with zero mean and known covariance:

$$ \mathbf{w}_k \sim \mathcal{N}(\mathbf{0}, \mathbf{Q}_k) \tag{2.3} $$

$$ \mathbf{v}_k \sim \mathcal{N}(\mathbf{0}, \mathbf{R}_k) \tag{2.4} $$

where $\mathbf{Q}_k \in \mathbb{R}^{n \times n}$ is the process noise covariance and $\mathbf{R}_k \in \mathbb{R}^{m \times m}$ is the measurement noise covariance. Both matrices are symmetric and positive semi-definite. $\mathbf{R}_k$ is typically positive definite in practice; a measurement with zero noise in some direction is rarely physical.

Zero mean. The noise has no systematic bias. Any bias in the dynamics or measurements is assumed to be modeled in $\mathbf{F}_k$ or $\mathbf{H}_k$ rather than absorbed into the noise. A noise process with non-zero mean is equivalent to a deterministic forcing term and should be treated as such.

Known covariance. The covariance matrices $\mathbf{Q}_k$ and $\mathbf{R}_k$ are specified by the modeler. They are not estimated by the filter during operation. Misspecified covariances degrade filter performance and are detectable through innovation sequence diagnostics. Section 7 develops the statistical properties of the innovation that make this detection possible.

Whiteness. The noise is uncorrelated across time:

$$ \mathbb{E}[\mathbf{w}_k \mathbf{w}_j^T] = \mathbf{Q}_k \delta_{kj} \tag{2.5} $$

$$ \mathbb{E}[\mathbf{v}_k \mathbf{v}_j^T] = \mathbf{R}_k \delta_{kj} \tag{2.6} $$

where $\delta_{kj}$ is the Kronecker delta. Each noise realization is independent of every other.

Mutual independence. The process noise and measurement noise are independent of each other at all times:

$$ \mathbb{E}[\mathbf{w}_k \mathbf{v}_j^T] = \mathbf{0} \quad \text{for all } k, j \tag{2.7} $$

Initial state. The initial state $\mathbf{x}_0$ is Gaussian with known mean and covariance:

$$ \mathbf{x}_0 \sim \mathcal{N}(\hat{\mathbf{x}}_0, \mathbf{P}_0) \tag{2.8} $$

The initial state is independent of all process and measurement noise.

These assumptions are restrictive. Real systems violate them in specific, characterizable ways. Section 7 develops the diagnostics that detect when the assumptions hold and when they fail.

2.4 Notation table

The notation used throughout this derivation:

Symbol Meaning Dimension
$\mathbf{x}_k$ True state at time $k$ $n$
$\hat{\mathbf{x}}_k$ Estimated state at time $k$ $n$
$\hat{\mathbf{x}}_k^-$ Predicted state at time $k$ (before measurement) $n$
$\mathbf{P}_k$ State estimate covariance at time $k$ $n \times n$
$\mathbf{P}_k^-$ Predicted state covariance at time $k$ $n \times n$
$\mathbf{F}_k$ State transition matrix $n \times n$
$\mathbf{Q}_k$ Process noise covariance $n \times n$
$\mathbf{z}_k$ Measurement at time $k$ $m$
$\mathbf{H}_k$ Measurement matrix $m \times n$
$\mathbf{R}_k$ Measurement noise covariance $m \times m$
$\mathbf{w}_k$ Process noise $n$
$\mathbf{v}_k$ Measurement noise $m$
$\boldsymbol{\nu}_k$ Innovation (defined in Section 5.1) $m$
$\mathbf{S}_k$ Innovation covariance (defined in Section 5.1) $m \times m$
$\mathbf{K}_k$ Kalman gain (defined in Section 5.2) $n \times m$
$n$ State dimension scalar
$m$ Measurement dimension scalar

State estimates are distinguished from true state by the hat: $\hat{\mathbf{x}}$ denotes an estimate of $\mathbf{x}$. Predicted quantities (before incorporating a measurement) carry a superscript minus: $\hat{\mathbf{x}}_k^-$ is the state estimate at time $k$ before the measurement $\mathbf{z}_k$ has been processed. Updated quantities (after incorporating the measurement) carry no superscript: $\hat{\mathbf{x}}_k$ is the state estimate after $\mathbf{z}_k$ has been processed.

The transition between predicted and updated quantities is the heart of the filter. The next sections derive that transition.

3. The Bayesian Foundation

3.1 Recursive Bayesian estimation

The state estimation problem can be stated in Bayesian terms. At time $k$, we have observed the measurements $\mathbf{z}_{1:k} = {\mathbf{z}_1, \mathbf{z}_2, \ldots, \mathbf{z}_k}$. We want the posterior distribution of the state given those measurements:

$$ p(\mathbf{x}_k \mid \mathbf{z}_{1:k}) \tag{3.1} $$

This distribution captures everything that can be known about the state at time $k$ from the available data. Its mean is the best estimate. Its covariance is the uncertainty in that estimate.

The recursive structure comes from a simple observation. The posterior at time $k$ can be computed from the posterior at time $k-1$ without revisiting all past measurements. This is what makes the filter a filter rather than a batch estimator.

The recursion proceeds in two steps:

$$ p(\mathbf{x}_{k-1} \mid \mathbf{z}_{1:k-1}) ;\xrightarrow{\text{predict}}; p(\mathbf{x}_k \mid \mathbf{z}_{1:k-1}) ;\xrightarrow{\text{update}}; p(\mathbf{x}_k \mid \mathbf{z}_{1:k}) \tag{3.2} $$

The prediction step takes the posterior from the previous time step forward through the dynamics, producing a prior at the current time step. The update step incorporates the new measurement, producing the new posterior. Each step has a precise mathematical form derived in the sections that follow.

3.2 The prediction-update structure

The prediction step is an application of the law of total probability. Given the posterior at time $k-1$ and the process model from Section 2.1, the predicted distribution at time $k$ is:

$$ p(\mathbf{x}_k \mid \mathbf{z}_{1:k-1}) = \int p(\mathbf{x}_k \mid \mathbf{x}_{k-1}) , p(\mathbf{x}_{k-1} \mid \mathbf{z}_{1:k-1}) , d\mathbf{x}_{k-1} \tag{3.3} $$

The integral marginalizes over the previous state. The transition density $p(\mathbf{x}k \mid \mathbf{x}{k-1})$ encodes the process model. The result is the distribution of $\mathbf{x}_k$ conditioned on all measurements up to but not including time $k$.

The update step is an application of Bayes' rule. Given the predicted distribution and the measurement model from Section 2.2, the posterior after observing $\mathbf{z}_k$ is:

$$ p(\mathbf{x}_k \mid \mathbf{z}_{1:k}) = \frac{p(\mathbf{z}_k \mid \mathbf{x}_k) , p(\mathbf{x}_k \mid \mathbf{z}_{1:k-1})}{p(\mathbf{z}_k \mid \mathbf{z}_{1:k-1})} \tag{3.4} $$

The likelihood $p(\mathbf{z}_k \mid \mathbf{x}_k)$ encodes the measurement model. The denominator is a normalizing constant. The numerator is the product of the prior at time $k$ and the likelihood of the observed measurement.

Equations (3.3) and (3.4) are general. They hold for any process and measurement model with well-defined transition and likelihood densities. The Kalman filter is what these equations become under the linear-Gaussian assumptions of Section 2.

3.3 Why Gaussian and linear produce closed form

The integrals and ratios in equations (3.3) and (3.4) are generally intractable. For arbitrary process and measurement models, the posterior at time $k$ may have no analytical form. Particle filters and other Monte Carlo methods exist precisely because most state estimation problems do not admit closed-form solutions.

The linear-Gaussian case is the exception. Under the assumptions of Section 2:

  • The transition density $p(\mathbf{x}k \mid \mathbf{x}{k-1})$ is Gaussian
  • The likelihood $p(\mathbf{z}_k \mid \mathbf{x}_k)$ is Gaussian
  • The initial distribution $p(\mathbf{x}_0)$ is Gaussian

Two facts about Gaussian distributions make the recursion close in closed form.

Linear transformations of Gaussians are Gaussian. If $\mathbf{x} \sim \mathcal{N}(\boldsymbol{\mu}, \boldsymbol{\Sigma})$ and $\mathbf{y} = \mathbf{A}\mathbf{x} + \mathbf{b}$, then $\mathbf{y} \sim \mathcal{N}(\mathbf{A}\boldsymbol{\mu} + \mathbf{b}, \mathbf{A}\boldsymbol{\Sigma}\mathbf{A}^T)$. This makes the prediction step in equation (3.3) analytically tractable. The integral over a Gaussian transition density and a Gaussian prior produces a Gaussian result.

Conditional distributions of jointly Gaussian variables are Gaussian. If $\mathbf{x}$ and $\mathbf{z}$ are jointly Gaussian, then $p(\mathbf{x} \mid \mathbf{z})$ is Gaussian with mean and covariance given by the Gaussian conditioning formula (Appendix B). This makes the update step in equation (3.4) analytically tractable. The ratio of Gaussian likelihood and Gaussian prior produces a Gaussian posterior.

The combination of these two facts means that if the posterior at time $k-1$ is Gaussian, then the posterior at time $k$ is also Gaussian. The recursion preserves the form. Because the initial distribution is Gaussian by assumption, the posterior at every time step is Gaussian.

A Gaussian distribution is fully characterized by its mean and covariance. The recursion therefore reduces from "propagate a distribution" to "propagate two quantities: a mean vector and a covariance matrix." That reduction is what makes the Kalman filter computationally tractable.

The sections that follow derive the explicit recursions for the mean and covariance. Section 4 derives the prediction step. Section 5 derives the update step. Section 6 establishes the optimality properties that hold under the linear-Gaussian assumptions.

4. The Prediction Step

4.1 Propagating the mean

The prediction step propagates the state estimate forward through the dynamics. Starting from the posterior at time $k-1$:

$$ p(\mathbf{x}_{k-1} \mid \mathbf{z}_{1:k-1}) = \mathcal{N}(\hat{\mathbf{x}}_{k-1}, \mathbf{P}_{k-1}) \tag{4.1} $$

The predicted state at time $k$ is the conditional expectation of $\mathbf{x}_k$ given the measurements through time $k-1$:

$$ \hat{\mathbf{x}}_k^- = \mathbb{E}[\mathbf{x}_k \mid \mathbf{z}_{1:k-1}] \tag{4.2} $$

Substituting the process model from equation (2.1):

$$ \hat{\mathbf{x}}_k^- = \mathbb{E}[\mathbf{F}_k \mathbf{x}_{k-1} + \mathbf{w}_k \mid \mathbf{z}_{1:k-1}] \tag{4.3} $$

Expectation is linear. The process noise $\mathbf{w}_k$ is independent of all prior measurements and has zero mean. Equation (4.3) reduces to:

$$ \hat{\mathbf{x}}_k^- = \mathbf{F}_k , \mathbb{E}[\mathbf{x}_{k-1} \mid \mathbf{z}_{1:k-1}] + \mathbb{E}[\mathbf{w}_k] \tag{4.4} $$

The first expectation is the posterior mean from time $k-1$. The second expectation is zero by the noise assumption in equation (2.3). The predicted mean is therefore:

$$ \boxed{\hat{\mathbf{x}}_k^- = \mathbf{F}_k \hat{\mathbf{x}}_{k-1}} \tag{4.5} $$

The predicted state is the previous state estimate propagated through the state transition matrix. No measurement information has been used. The dynamics alone determine the predicted state.

4.2 Propagating the covariance

The predicted covariance is the conditional covariance of $\mathbf{x}_k$ given the measurements through time $k-1$:

$$ \mathbf{P}_k^- = \mathbb{E}\left[(\mathbf{x}_k - \hat{\mathbf{x}}_k^-)(\mathbf{x}_k - \hat{\mathbf{x}}_k^-)^T \mid \mathbf{z}_{1:k-1}\right] \tag{4.6} $$

Substituting the process model and the predicted mean from equation (4.5):

$$ \mathbf{x}_k - \hat{\mathbf{x}}_k^- = \mathbf{F}_k \mathbf{x}_{k-1} + \mathbf{w}_k - \mathbf{F}_k \hat{\mathbf{x}}_{k-1} = \mathbf{F}_k(\mathbf{x}_{k-1} - \hat{\mathbf{x}}_{k-1}) + \mathbf{w}_k \tag{4.7} $$

The error in the predicted state has two contributions: the propagated error from the previous estimate, and the new process noise. Both are zero-mean. Substituting equation (4.7) into equation (4.6) and expanding the outer product:

$$ \mathbf{P}_k^- = \mathbb{E}\left[\mathbf{F}_k(\mathbf{x}_{k-1} - \hat{\mathbf{x}}_{k-1})(\mathbf{x}_{k-1} - \hat{\mathbf{x}}_{k-1})^T \mathbf{F}_k^T \mid \mathbf{z}_{1:k-1}\right] $$

  • \mathbb{E}\left[\mathbf{F}k(\mathbf{x}{k-1} - \hat{\mathbf{x}}_{k-1}) \mathbf{w}k^T \mid \mathbf{z}{1:k-1}\right] $$ $$
  • \mathbb{E}\left[\mathbf{w}k (\mathbf{x}{k-1} - \hat{\mathbf{x}}_{k-1})^T \mathbf{F}k^T \mid \mathbf{z}{1:k-1}\right] $$ $$
  • \mathbb{E}\left[\mathbf{w}_k \mathbf{w}k^T \mid \mathbf{z}{1:k-1}\right] \tag{4.8} $$

The two cross terms vanish. The process noise $\mathbf{w}k$ is independent of $\mathbf{x}{k-1}$ and of all prior measurements; the expectation of the product factors, and one factor is zero. The remaining terms are the propagated covariance and the process noise covariance:

$$ \boxed{\mathbf{P}_k^- = \mathbf{F}_k \mathbf{P}_{k-1} \mathbf{F}_k^T + \mathbf{Q}_k} \tag{4.9} $$

The predicted covariance has two terms. The first term, $\mathbf{F}k \mathbf{P}{k-1} \mathbf{F}_k^T$, is the prior covariance transformed by the dynamics. The second term, $\mathbf{Q}_k$, is the additional uncertainty introduced by the process noise. The covariance always grows in the prediction step. Time without measurement increases uncertainty.

4.3 The predicted state distribution

Equations (4.5) and (4.9) give the mean and covariance of the predicted state. Because the dynamics are linear and the noise is Gaussian, the predicted state distribution is itself Gaussian:

$$ p(\mathbf{x}_k \mid \mathbf{z}_{1:k-1}) = \mathcal{N}(\hat{\mathbf{x}}_k^-, \mathbf{P}_k^-) \tag{4.10} $$

This is the result of the marginalization integral in equation (3.3) under linear-Gaussian assumptions. The integral does not need to be evaluated explicitly. The Gaussian closure property from Section 3.3 guarantees the result.

The predicted distribution carries everything the dynamics know about the state at time $k$ before the measurement $\mathbf{z}_k$ is processed. The mean $\hat{\mathbf{x}}_k^-$ is the best estimate of the state under the dynamics alone. The covariance $\mathbf{P}_k^-$ quantifies the uncertainty in that estimate.

The next section incorporates the measurement.

5. The Update Step

5.1 The innovation

When the measurement $\mathbf{z}_k$ arrives, it can be compared to what the predicted state suggests it should be. The expected measurement under the predicted state is $\mathbf{H}_k \hat{\mathbf{x}}_k^-$. The difference between the actual measurement and the expected measurement is the innovation:

$$ \boldsymbol{\nu}_k = \mathbf{z}_k - \mathbf{H}_k \hat{\mathbf{x}}_k^- \tag{5.1} $$

The innovation captures what the measurement adds beyond what was already predicted. If the predicted state were exact, the innovation would be zero. The innovation is nonzero because of two sources: the error in the predicted state, and the measurement noise.

Substituting the measurement model from equation (2.2) and the predicted state from Section 4:

$$ \boldsymbol{\nu}_k = \mathbf{H}_k \mathbf{x}_k + \mathbf{v}_k - \mathbf{H}_k \hat{\mathbf{x}}_k^- = \mathbf{H}_k(\mathbf{x}_k - \hat{\mathbf{x}}_k^-) + \mathbf{v}_k \tag{5.2} $$

The first term is the prediction error projected through the measurement matrix. The second term is the measurement noise. Both are zero-mean, so the innovation has zero mean under the model assumptions:

$$ \mathbb{E}[\boldsymbol{\nu}_k] = \mathbf{0} \tag{5.3} $$

The innovation covariance is:

$$ \mathbf{S}_k = \mathbb{E}[\boldsymbol{\nu}_k \boldsymbol{\nu}_k^T] \tag{5.4} $$

Substituting equation (5.2) and expanding:

$$ \mathbf{S}_k = \mathbb{E}\left[\mathbf{H}_k(\mathbf{x}_k - \hat{\mathbf{x}}_k^-)(\mathbf{x}_k - \hat{\mathbf{x}}_k^-)^T \mathbf{H}_k^T\right] + \mathbb{E}[\mathbf{v}_k \mathbf{v}_k^T] \tag{5.5} $$

The cross terms vanish because the measurement noise is independent of the prediction error. The first remaining term is $\mathbf{H}_k \mathbf{P}_k^- \mathbf{H}_k^T$ by the definition of $\mathbf{P}_k^-$ in equation (4.6). The second remaining term is $\mathbf{R}_k$ by the noise assumption in equation (2.4). The innovation covariance is:

$$ \boxed{\mathbf{S}_k = \mathbf{H}_k \mathbf{P}_k^- \mathbf{H}_k^T + \mathbf{R}_k} \tag{5.6} $$

The innovation covariance has two terms. The first term is the prediction uncertainty projected into measurement space. The second term is the measurement noise. Both are uncertainty contributions to the difference between predicted and observed measurements.

The properties of the innovation sequence are central to filter diagnostics. Section 7 develops them in detail.

5.2 The Kalman gain

The updated state estimate is a linear combination of the predicted state and the innovation:

$$ \hat{\mathbf{x}}_k = \hat{\mathbf{x}}_k^- + \mathbf{K}_k \boldsymbol{\nu}_k \tag{5.7} $$

The matrix $\mathbf{K}_k \in \mathbb{R}^{n \times m}$ is the Kalman gain. It determines how much weight the new measurement receives relative to the prediction. The derivation that follows establishes the specific form of $\mathbf{K}_k$ that minimizes the updated covariance.

The updated covariance is:

$$ \mathbf{P}_k = \mathbb{E}\left[(\mathbf{x}_k - \hat{\mathbf{x}}_k)(\mathbf{x}_k - \hat{\mathbf{x}}_k)^T\right] \tag{5.8} $$

Substituting the update equation (5.7):

$$ \mathbf{x}_k - \hat{\mathbf{x}}_k = \mathbf{x}_k - \hat{\mathbf{x}}_k^- - \mathbf{K}_k \boldsymbol{\nu}_k \tag{5.9} $$

Substituting the innovation from equation (5.2):

$$ \mathbf{x}_k - \hat{\mathbf{x}}_k = (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k)(\mathbf{x}_k - \hat{\mathbf{x}}_k^-) - \mathbf{K}_k \mathbf{v}_k \tag{5.10} $$

The updated error is a linear combination of the prediction error and the measurement noise. Substituting into the covariance definition in equation (5.8) and expanding the outer product:

$$ \mathbf{P}_k = (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k) \mathbf{P}_k^- (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k)^T + \mathbf{K}_k \mathbf{R}_k \mathbf{K}_k^T \tag{5.11} $$

The cross terms vanish by the independence of the measurement noise and prediction error. Equation (5.11) gives the updated covariance as a function of the choice of $\mathbf{K}_k$. This is the Joseph form of the covariance update, developed further in Section 5.4.

The Kalman gain is the choice of $\mathbf{K}_k$ that minimizes the trace of $\mathbf{P}_k$. The trace is the sum of the variances of the state estimate components; minimizing it minimizes the total mean squared error in the updated estimate.

Differentiating the trace of equation (5.11) with respect to $\mathbf{K}_k$ and setting the derivative to zero:

$$ \frac{\partial , \text{tr}(\mathbf{P}_k)}{\partial \mathbf{K}_k} = -2(\mathbf{I} - \mathbf{K}_k \mathbf{H}_k) \mathbf{P}_k^- \mathbf{H}_k^T + 2 \mathbf{K}_k \mathbf{R}_k = \mathbf{0} \tag{5.12} $$

Rearranging:

$$ \mathbf{K}_k (\mathbf{H}_k \mathbf{P}_k^- \mathbf{H}_k^T + \mathbf{R}_k) = \mathbf{P}_k^- \mathbf{H}_k^T \tag{5.13} $$

The matrix on the left is the innovation covariance $\mathbf{S}_k$ from equation (5.6). Solving for the gain:

$$ \boxed{\mathbf{K}_k = \mathbf{P}_k^- \mathbf{H}_k^T \mathbf{S}_k^{-1}} \tag{5.14} $$

The Kalman gain weighs the prediction uncertainty against the innovation uncertainty. The numerator $\mathbf{P}_k^- \mathbf{H}_k^T$ is the cross-covariance between the state and the measurement. The denominator $\mathbf{S}_k^{-1}$ is the inverse of the innovation covariance. Large prediction uncertainty pushes the gain higher, giving more weight to the measurement. Large measurement noise pushes the gain lower, giving more weight to the prediction. The filter balances the two automatically.

5.3 The updated state distribution

Substituting the optimal gain from equation (5.14) into the covariance expression in equation (5.11) and simplifying produces the standard form of the updated covariance:

$$ \boxed{\mathbf{P}_k = (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k) \mathbf{P}_k^-} \tag{5.15} $$

The updated covariance is the predicted covariance reduced by the information contributed by the measurement. The factor $\mathbf{I} - \mathbf{K}_k \mathbf{H}_k$ projects out the dimensions of the state that the measurement informs about.

The updated state estimate from equation (5.7), the updated covariance from equation (5.15), and the Gaussian closure property from Section 3.3 give the posterior distribution at time $k$:

$$ p(\mathbf{x}_k \mid \mathbf{z}_{1:k}) = \mathcal{N}(\hat{\mathbf{x}}_k, \mathbf{P}_k) \tag{5.16} $$

This is the same result that would be obtained by applying the Gaussian conditioning formula (Appendix B) directly to the joint distribution of $\mathbf{x}_k$ and $\mathbf{z}_k$. The minimum-variance derivation and the conditional-distribution derivation produce the same answer under linear-Gaussian assumptions. The agreement is expected. Section 6 develops the optimality properties that explain why.

5.4 The Joseph form

Equation (5.11) gave the updated covariance as a function of any gain $\mathbf{K}_k$, including non-optimal choices. With the optimal gain from equation (5.14), equation (5.11) simplifies to the standard form in equation (5.15). Algebraically the two are equivalent.

Numerically they are not.

The standard form $\mathbf{P}_k = (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k) \mathbf{P}_k^-$ involves a subtraction. In finite-precision arithmetic, subtraction of near-equal quantities introduces error. Over many filter iterations, this error can compound, causing the updated covariance $\mathbf{P}_k$ to lose symmetry or even positive-definiteness. A covariance matrix that is not positive-definite is a structural failure of the filter.

The Joseph form preserves positive-definiteness:

$$ \mathbf{P}_k = (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k) \mathbf{P}_k^- (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k)^T + \mathbf{K}_k \mathbf{R}_k \mathbf{K}_k^T \tag{5.17} $$

Both terms in equation (5.17) are symmetric positive-semi-definite products. Their sum is symmetric positive-semi-definite by construction. The Joseph form is numerically robust even when the gain is slightly suboptimal due to rounding error or covariance specification mismatch.

The reference implementation uses the Joseph form. The standard form is mathematically cleaner; the Joseph form is what production estimation deploys.

6. Optimality

6.1 MMSE optimality under linear-Gaussian assumptions

The Kalman gain in equation (5.14) was derived by minimizing the trace of the updated covariance. The trace is the sum of the variances of the state estimate components. Minimizing it minimizes the total mean squared error of the estimate:

$$ \text{MSE}(\hat{\mathbf{x}}_k) = \mathbb{E}\left[(\mathbf{x}_k - \hat{\mathbf{x}}_k)^T (\mathbf{x}_k - \hat{\mathbf{x}}_k)\right] = \text{tr}(\mathbf{P}_k) \tag{6.1} $$

The Kalman filter is the minimum mean squared error (MMSE) estimator within the class of linear estimators. No linear combination of the prediction and the measurement produces a smaller mean squared error.

Under the linear-Gaussian assumptions of Section 2, a stronger statement holds. The Kalman filter is the MMSE estimator within the class of all estimators, linear or otherwise. The linear restriction is not a restriction at all when the underlying distributions are Gaussian.

The reason is the Gaussian conditioning property. When $\mathbf{x}_k$ and $\mathbf{z}_k$ are jointly Gaussian, the conditional mean $\mathbb{E}[\mathbf{x}k \mid \mathbf{z}{1:k}]$ is a linear function of the measurements. Nonlinear estimators cannot do better than linear ones because the optimal estimator is already linear. The Kalman filter computes this linear conditional mean exactly.

This is the strongest optimality claim available. Under linear-Gaussian assumptions, the Kalman filter is not just the best linear estimator. It is the best estimator.

6.2 The orthogonality principle

The minimum-variance derivation in Section 5.2 produces the optimal gain through calculus. The orthogonality principle produces the same gain through a geometric argument. Both paths arrive at the same result. Seeing both clarifies the structural reason for the optimality.

The orthogonality principle states that the optimal estimation error is uncorrelated with the data used to form the estimate. Formally, for the optimal $\hat{\mathbf{x}}_k$:

$$ \mathbb{E}\left[(\mathbf{x}_k - \hat{\mathbf{x}}_k) \mathbf{z}_j^T\right] = \mathbf{0} \quad \text{for all } j \leq k \tag{6.2} $$

The error of the optimal estimate is orthogonal to every measurement used in computing the estimate. If the error were correlated with any measurement, that correlation could be exploited to reduce the error further. The estimate would not be optimal. Optimality requires the correlation to vanish.

Applying the orthogonality condition to the innovation:

$$ \mathbb{E}\left[(\mathbf{x}_k - \hat{\mathbf{x}}_k) \boldsymbol{\nu}_k^T\right] = \mathbf{0} \tag{6.3} $$

Substituting the updated estimate from equation (5.7) and the innovation from equation (5.2) into equation (6.3), then solving for $\mathbf{K}_k$, produces the same gain as the minimum-variance derivation. The two paths are different views of the same structure.

The orthogonality principle is geometric. The optimal estimate $\hat{\mathbf{x}}_k$ is the projection of $\mathbf{x}_k$ onto the subspace spanned by the measurements. The error $\mathbf{x}_k - \hat{\mathbf{x}}_k$ is the component of $\mathbf{x}_k$ orthogonal to that subspace. Projection minimizes distance. The Kalman filter is performing orthogonal projection in the Hilbert space of square-integrable random variables.

This view does not change the equations. It changes how the equations are understood.

6.3 What "optimal" means and what it doesn't

The optimality of the Kalman filter is conditional. The conditions are the assumptions of Section 2:

  • Linear dynamics and measurements
  • Gaussian noise with known covariance
  • Process noise and measurement noise mutually independent and white
  • Initial state Gaussian with known mean and covariance

When all four conditions hold, the filter is the MMSE estimator and no improvement is possible. When any condition is weakened, the optimality claim weakens with it.

Nonlinear dynamics or measurements. The Extended Kalman Filter linearizes the dynamics or measurement model around the current estimate and applies the linear Kalman update to the linearized system. The result is an approximation. It is no longer optimal in any formal sense, though it is often the best tractable approximation available. The Unscented Kalman Filter uses a deterministic sampling approach to propagate the mean and covariance through the nonlinear functions directly. It typically achieves better accuracy than the EKF for strongly nonlinear systems, but it too is an approximation.

Non-Gaussian noise. When the noise distribution is non-Gaussian, the Kalman filter remains the best linear estimator but is no longer the best estimator overall. Particle filters and other Monte Carlo methods can outperform the Kalman filter in this regime by representing the non-Gaussian posterior explicitly.

Misspecified covariance. The filter assumes $\mathbf{Q}_k$ and $\mathbf{R}_k$ are known. In practice they are modeled, estimated, or guessed. When they are wrong, the filter produces a suboptimal estimate whose error is larger than the covariance suggests. The innovation sequence statistics developed in Section 7 detect this failure mode.

Correlated noise. When the process and measurement noises are correlated, or when either is colored rather than white, the standard Kalman filter is no longer optimal. Modified forms of the filter exist that handle these cases when the correlation structure is known.

Optimality is a property of the model, not a property of reality. The filter is optimal given the model. Whether the model fits reality is a separate question. The diagnostics in Section 7 provide tools to assess that question quantitatively.

This is the correct framing for production deployment. The filter is not magic. It is the best estimator under specific assumptions. The work of an estimation engineer is to determine whether those assumptions hold, to characterize how they fail when they do, and to deploy the filter with full knowledge of both.

7. Properties of the Innovation Sequence

The innovation sequence ${\boldsymbol{\nu}_1, \boldsymbol{\nu}_2, \ldots, \boldsymbol{\nu}_k}$ is the most informative signal the filter produces. Under the assumptions of Section 2, the innovation has three statistical properties: zero mean, known covariance, and whiteness. These properties hold by mathematical necessity when the model is correct. They fail in specific, characterizable ways when the model is wrong. The properties are therefore both a verification that the filter is operating correctly and a diagnostic that surfaces when it is not.

7.1 Whiteness under correct model

The innovation sequence is white. Innovations at different times are mutually uncorrelated:

$$ \mathbb{E}[\boldsymbol{\nu}_k \boldsymbol{\nu}_j^T] = \mathbf{0} \quad \text{for all } k \neq j \tag{7.1} $$

The whiteness property follows from the optimality of the Kalman filter. The orthogonality principle from Section 6.2 states that the optimal estimation error is uncorrelated with every measurement used in computing the estimate. The innovation at time $k$ is a linear function of the prediction error at time $k$ and the measurement noise at time $k$. The prediction error at time $k$ depends on the estimate at time $k-1$, which was constructed to be orthogonal to all measurements through time $k-1$. The measurement noise at time $k$ is independent of all prior measurement noise by assumption.

Combining these: the innovation at time $k$ is uncorrelated with every prior measurement and with every prior measurement noise. Since prior innovations are linear combinations of prior measurements and prior measurement noise, the innovation at time $k$ is uncorrelated with every prior innovation.

The whiteness of the innovation sequence is the signature of optimal filtering. A correctly specified filter produces a white innovation sequence. A misspecified filter does not.

7.2 Zero mean

The innovation has zero mean:

$$ \mathbb{E}[\boldsymbol{\nu}_k] = \mathbf{0} \tag{7.2} $$

This was established in equation (5.3). The derivation is direct. The innovation is the difference between the actual measurement and the predicted measurement:

$$ \boldsymbol{\nu}_k = \mathbf{z}_k - \mathbf{H}_k \hat{\mathbf{x}}_k^- = \mathbf{H}_k(\mathbf{x}_k - \hat{\mathbf{x}}_k^-) + \mathbf{v}_k \tag{7.3} $$

Both terms have zero mean. The prediction error has zero mean because the predicted state is the conditional mean of $\mathbf{x}_k$ given prior measurements. The measurement noise has zero mean by the assumption in equation (2.3). Their sum has zero mean.

A nonzero innovation mean indicates a bias in the filter's predictions. The cause is structural: either the dynamics model $\mathbf{F}_k$ is incorrect, the measurement model $\mathbf{H}_k$ is incorrect, or one of the noise processes has a nonzero mean that has not been absorbed into the deterministic part of the model. None of these conditions can be corrected by tuning the noise covariances. They require correcting the model itself.

7.3 Known covariance

The innovation covariance equals $\mathbf{S}_k$ exactly:

$$ \mathbb{E}[\boldsymbol{\nu}_k \boldsymbol{\nu}_k^T] = \mathbf{S}_k = \mathbf{H}_k \mathbf{P}_k^- \mathbf{H}_k^T + \mathbf{R}_k \tag{7.4} $$

This was established in equation (5.6). The innovation covariance is fully determined by the prediction covariance $\mathbf{P}_k^-$ and the measurement noise covariance $\mathbf{R}_k$. It does not require empirical estimation. The filter computes it as part of its normal operation.

The empirical innovation covariance, computed from the observed innovation sequence, should agree with the filter-computed $\mathbf{S}_k$ to within sampling variation. When the empirical covariance is systematically larger than $\mathbf{S}_k$, the filter is underestimating its own uncertainty. When the empirical covariance is systematically smaller, the filter is overestimating its uncertainty. Both cases indicate covariance misspecification.

The normalized innovation squared (NIS) statistic operationalizes this comparison:

$$ \epsilon_{\nu, k} = \boldsymbol{\nu}_k^T \mathbf{S}_k^{-1} \boldsymbol{\nu}_k \tag{7.5} $$

Under the model assumptions, $\epsilon_{\nu, k}$ follows a chi-squared distribution with $m$ degrees of freedom, where $m$ is the measurement dimension. A time-averaged NIS that falls outside the expected chi-squared range indicates that the empirical innovation covariance disagrees with the filter-computed $\mathbf{S}_k$. The disagreement is detectable from a finite sample.

7.4 Why these properties matter

The three innovation properties are whiteness, zero mean, and known covariance. They are jointly necessary and sufficient for the filter to be operating correctly under its specified model. When all three hold in the observed data, the model and the filter agree. When any fails, the model and the data disagree in a specific way:

Non-white innovations indicate that the filter is not extracting all available information from the measurements. The dynamics model or measurement model is missing structure that the data contains. This manifests as autocorrelation in the innovation sequence at one or more lags. The lag structure of the autocorrelation often points to the missing model component.

Non-zero mean innovations indicate bias in the filter's predictions. The structural causes are model error in $\mathbf{F}_k$ or $\mathbf{H}_k$, or unmodeled deterministic forcing. Tuning noise covariances does not correct bias.

Innovation covariance mismatch indicates that $\mathbf{Q}_k$ or $\mathbf{R}_k$ are misspecified. The filter's uncertainty estimate is inconsistent with the actual error magnitude. NIS statistics quantify this mismatch and identify whether the filter is over-confident or under-confident.

The innovation sequence is therefore a self-diagnosing structure. It contains within itself the information needed to detect when the filter is operating outside its assumptions. The diagnostics that extract this information are developed in the implementation in diagnostics.py.

This is the structural reason for caring about innovation properties in production estimation. A filter that runs without monitoring its innovation sequence is operating blind. A filter that monitors its innovation sequence knows when it is performing correctly and when it is not. The work of an estimation engineer is to deploy the latter.

8. The Complete Algorithm

The Kalman filter is a recursive estimator with three operational phases: initialization, prediction, and update. Sections 4 and 5 derived the equations for prediction and update. Section 7 established the statistical properties that govern correct operation. This section consolidates the result into a complete algorithm.

8.1 Initialization

The filter requires an initial state estimate and an initial covariance:

$$ \hat{\mathbf{x}}_0, \quad \mathbf{P}_0 \tag{8.1} $$

The initial estimate $\hat{\mathbf{x}}_0$ represents the best prior knowledge of the state before any measurements have been processed. The initial covariance $\mathbf{P}_0$ represents the uncertainty in that prior knowledge.

When prior information is strong, $\hat{\mathbf{x}}_0$ should be set to the known mean and $\mathbf{P}_0$ should reflect the known uncertainty. When prior information is weak, $\mathbf{P}_0$ should be set large enough that the first few measurements dominate the estimate. A large $\mathbf{P}_0$ models an uninformative prior.

The filter is sensitive to initialization for the first several measurements. After enough measurements have been processed, the effect of the initialization fades. The filter converges to estimates that depend primarily on the dynamics, the measurements, and the noise covariances.

The process and measurement noise covariances $\mathbf{Q}_k$ and $\mathbf{R}_k$ must also be specified before the filter operates. These are model parameters, not internal filter state. They reflect the modeler's assumptions about the system, not estimates the filter produces.

8.2 Predict

The prediction step propagates the state estimate and its covariance forward through the dynamics:

$$ \hat{\mathbf{x}}_k^- = \mathbf{F}_k \hat{\mathbf{x}}_{k-1} \tag{8.2} $$

$$ \mathbf{P}_k^- = \mathbf{F}_k \mathbf{P}_{k-1} \mathbf{F}_k^T + \mathbf{Q}_k \tag{8.3} $$

The predicted state is the previous estimate propagated through the state transition matrix. The predicted covariance is the previous covariance propagated through the same matrix, plus the additional uncertainty introduced by the process noise.

The prediction step requires no measurement. It can be applied at each time step regardless of whether a measurement is available. When measurements arrive irregularly, the filter can predict forward through multiple steps before incorporating the next measurement.

8.3 Update

The update step incorporates a new measurement into the state estimate. It begins with the innovation:

$$ \boldsymbol{\nu}_k = \mathbf{z}_k - \mathbf{H}_k \hat{\mathbf{x}}_k^- \tag{8.4} $$

The innovation covariance follows from the prediction covariance and the measurement noise covariance:

$$ \mathbf{S}_k = \mathbf{H}_k \mathbf{P}_k^- \mathbf{H}_k^T + \mathbf{R}_k \tag{8.5} $$

The Kalman gain weighs the prediction uncertainty against the innovation uncertainty:

$$ \mathbf{K}_k = \mathbf{P}_k^- \mathbf{H}_k^T \mathbf{S}_k^{-1} \tag{8.6} $$

The updated state estimate combines the prediction with the innovation, weighted by the gain:

$$ \hat{\mathbf{x}}_k = \hat{\mathbf{x}}_k^- + \mathbf{K}_k \boldsymbol{\nu}_k \tag{8.7} $$

The updated covariance uses the Joseph form for numerical robustness:

$$ \mathbf{P}_k = (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k) \mathbf{P}_k^- (\mathbf{I} - \mathbf{K}_k \mathbf{H}_k)^T + \mathbf{K}_k \mathbf{R}_k \mathbf{K}_k^T \tag{8.8} $$

The update step requires a measurement. When no measurement is available at time $k$, the update step is skipped and the filter output is the predicted state and covariance.

8.4 Pseudocode form

The complete algorithm in pseudocode:

Inputs:
    F[k], Q[k]                       # process model at each time step
    H[k], R[k]                       # measurement model at each time step
    z[k]                             # measurements at each time step
    x_hat_0, P_0                     # initial state and covariance

Initialize:
    x_hat = x_hat_0
    P = P_0

For k = 1, 2, ...:

    # Predict
    x_pred = F[k] @ x_hat
    P_pred = F[k] @ P @ F[k].T + Q[k]

    if measurement z[k] is available:

        # Compute innovation and its covariance
        nu = z[k] - H[k] @ x_pred
        S = H[k] @ P_pred @ H[k].T + R[k]

        # Compute Kalman gain
        K = P_pred @ H[k].T @ inverse(S)

        # Update state estimate
        x_hat = x_pred + K @ nu

        # Update covariance (Joseph form)
        I_KH = I - K @ H[k]
        P = I_KH @ P_pred @ I_KH.T + K @ R[k] @ K.T

        # Diagnostic: record innovation and its statistics
        record(nu, S)

    else:

        # No measurement; carry prediction forward
        x_hat = x_pred
        P = P_pred

    output(x_hat, P)

The reference implementation in kalman_filter.py follows this structure directly. The diagnostic recording step feeds the innovation sequence analysis tools in diagnostics.py developed from the properties in Section 7.

The complete recursion is contained in equations (8.2) through (8.8). Given the inputs and initial conditions, the filter produces estimates that are optimal under the assumptions of Section 2 and self-diagnosing through the innovation sequence developed in Section 7.