The Kalman Filter, Derived and Tuned

The linear-Gaussian world

Assume the dynamics are linear, the observation is linear, and both noises are Gaussian. Then a beautiful thing happens: the posterior is always Gaussian, so it is fully described by just its mean \hat{x}_k and covariance P_k. The intractable Bayes integrals from guide 1 reduce to closed-form matrix algebra on those two objects. That algebra is the Kalman filter — the recursive, optimal (minimum mean-squared-error) linear estimator.

The predict step

Predicting means pushing the previous estimate through the dynamics and letting uncertainty grow by the process noise W:

\hat{x}_k^{-} = A\,\hat{x}_{k-1}, \qquad P_k^{-} = A\,P_{k-1}\,A^{\top} + W

The minus superscript marks the a-priori (pre-measurement) estimate.

Push last step's best estimate forward through the dynamics to guess the new state, and grow the uncertainty — because prediction alone always adds a bit of doubt.

\hat{x}_k^{-}: The predicted (a-priori) state, before this step's measurement.
P_k^{-}: The predicted uncertainty (covariance) of that estimate.
A: The dynamics matrix carrying the state forward.
W: The process noise added at each step, which inflates P.

This is the predict half of the Kalman filter loop.

The update step: gain, innovation, correction

Now fold in the measurement. The innovation r_k is the surprise — how far the actual neural reading is from what the prediction expected. The Kalman gain K_k decides how much of that surprise to believe, and the state is nudged accordingly.

\begin{aligned} \text{innovation:}\quad & r_k = y_k - C\,\hat{x}_k^{-} \\ \text{gain:}\quad & K_k = P_k^{-} C^{\top}\big(C\,P_k^{-} C^{\top} + Q\big)^{-1} \\ \text{update:}\quad & \hat{x}_k = \hat{x}_k^{-} + K_k\,r_k, \qquad P_k = (I - K_k C)\,P_k^{-} \end{aligned}

The gain optimally balances prediction uncertainty P against measurement noise Q.

See how wrong the prediction was (the innovation), then nudge the estimate toward the measurement by an amount set by the gain — trust the sensor more when it is reliable, the prediction more when it is not.

r_k: The innovation: measurement minus predicted measurement.
K_k: The Kalman gain — how far to move toward the data.
Q: The measurement noise covariance; large Q shrinks the gain.
P_k: The updated, shrunken uncertainty after correction.

The Kalman gain balances prediction uncertainty P against measurement noise Q; see the innovation sequence.

# one Kalman step (linear-Gaussian)
def kalman_step(x, P, y, A, W, C, Q):
    # predict
    x = A @ x
    P = A @ P @ A.T + W
    # update
    r = y - C @ x                    # innovation (surprise)
    S = C @ P @ C.T + Q              # innovation covariance
    K = P @ C.T @ np.linalg.inv(S)   # Kalman gain
    x = x + K @ r
    P = (np.eye(len(x)) - K @ C) @ P
    return x, P

The entire real-time decoder is these ~8 lines running per time bin.

The velocity Kalman filter

The version that actually drives implanted cursors is the velocity Kalman filter. Its state holds position and velocity, but the crucial design choice — learned the hard way in the BrainGate era — is to have neurons encode velocity, and to obtain position by integrating velocity rather than decoding it directly.

Batch cousins: the Wiener filter and OLE

If you ignore dynamics entirely and just regress kinematics on a sliding window of neural history, you get the Wiener filter / optimal linear estimator — the batch, non-recursive decoder that many labs still use as a baseline.

L^{\star} = \arg\min_{L}\; \mathbb{E}\big\lVert x_k - L\,\mathbf{y}_k \big\rVert^{2} = \Sigma_{xy}\,\Sigma_{yy}^{-1}, \qquad \hat{x}_k = L^{\star}\mathbf{y}_k

y_k stacks a lag window of neural features; the Wiener solution is a regularized covariance ratio.

The batch cousin of the Kalman filter: instead of recursing, directly fit one matrix that maps a window of neural features to the state, by regressing the state on those features. The answer is a covariance ratio.

\mathbf{y}_k: A stacked lag window of neural features.
L^{\star}: The optimal linear map from features to state.
\Sigma_{xy}: Cross-covariance between the state and the features.
\Sigma_{yy}^{-1}: Inverse feature covariance, which whitens the inputs.

This is the Wiener filter / optimal linear estimator — regularize \Sigma_{yy}^{-1} when features outnumber samples.

The Kalman filter is essentially the recursive, dynamics-aware upgrade of this idea: instead of a fixed regression it maintains a running estimate whose confidence adapts moment to moment. In practice the two often decode similarly offline — the Kalman filter's advantage shows up in closed loop, where its smoothness and integrator behaviour matter (guide 4).

Tuning the filter in practice

Fitting is system identification: estimate A, W from the autocorrelation of the training kinematics, and C, Q by regressing neural features on kinematics (next guide). But W and Q are effectively hyperparameters — many teams deliberately inflate Q or shrink W beyond their fitted values to buy extra smoothness, because a filter that is statistically optimal on open-loop data is often not the one that feels best to a user closing the loop.