Kalman Smoother Comparison: RTS vs Fixed-Lag vs Fixed-Point Smoothing in Python (Theory, Latency, and RMSE Benchmarks)

Compare RTS (fixed-interval), fixed-lag, and fixed-point Kalman smoothers in Python: recursions, complexity, and a from-scratch numpy implementation with numpy.linalg.inv, numpy.random.default_rng, and multivariate_normal. A 50-run Monte Carlo benchmark shows RMSE improving monotonically from filter to fixed-lag to RTS, with the lag-length N vs latency trade-off quantified.

Introduction: Filtering vs Smoothing

State estimation in state-space models comes in two flavors, depending on which observations you are allowed to use.

  • Filtering: estimation using observations up to time \(k\) → \(p(x_k \mid y_{1:k})\) . This is a causal operation: the estimate can be updated as each observation arrives, making it suitable for online, real-time processing.
  • Smoothing: estimation that also uses observations from the future → \(p(x_k \mid y_{1:T})\) with \(T > k\) . This is non-causal: exploiting future data always yields accuracy equal to or better than filtering, at the cost of latency.

Depending on “how far into the future you wait,” smoothers fall into three classes:

  1. Fixed-interval smoother: wait for the entire record \(y_{1:T}\) , then estimate \(\hat{x}_{k|T}\) for every \(k = 0, \ldots, T\) . The standard algorithm is the Rauch-Tung-Striebel (RTS) smoother. Fully offline.
  2. Fixed-lag smoother: always estimate the state \(N\) steps in the past, \(\hat{x}_{k-N|k}\) . Near-real-time processing with an acceptable delay of \(N\) steps.
  3. Fixed-point smoother: keep refining the estimate of one particular time \(j\) , \(\hat{x}_{j|k}\) , as new observations arrive. Used for initial-state or event-time estimation.

This article derives the recursions of all three smoothers, compares their complexity, latency, and use cases, and then benchmarks the Kalman filter, fixed-lag smoother, and RTS smoother with a from-scratch numpy implementation, measuring RMSE on the same simulated data.

For the underlying theory, see Kalman Filter: Theory and Python Implementation.

Comparison Table

Here is the summary up front. \(n\) is the state dimension, \(T\) the record length, \(N\) the lag.

AspectFixed-interval (RTS)Fixed-lagFixed-point
Estimate\(\hat{x}_{k\|T}\) for all \(k\)\(\hat{x}_{k-N\|k}\) (always \(N\) steps behind)\(\hat{x}_{j\|k}\) for one fixed \(j\)
Latencyfull record (offline)\(N\) steps (near real time)none (incremental updates)
Cost / step\(O(n^3)\) × forward + backward passes\(O(N n^3)\) (windowed RTS) / \(O((Nn)^3)\) (augmented)\(O(n^3)\)
Memory\(O(T n^2)\) (all intermediate results)\(O(N n^2)\)\(O(n^2)\)
Accuracybest (uses all observations)approaches RTS as \(N \to \infty\)approaches RTS at time \(j\) only
Typical usebatch analysis, trajectory reconstruction, EMcommunications, speech enhancement, trackinginitial-state / event-time refinement

Kalman Filter Review (Forward Pass)

All three smoothers consume the output of the Kalman filter forward pass. Define the linear Gaussian state-space model:

\[ x_k = F x_{k-1} + w_{k-1}, \quad w_{k-1} \sim \mathcal{N}(0, Q) \tag{1} \] \[ y_k = H x_k + v_k, \quad v_k \sim \mathcal{N}(0, R) \tag{2} \]

Prediction Step

\[ \hat{x}_{k|k-1} = F \hat{x}_{k-1|k-1} \tag{3} \] \[ P_{k|k-1} = F P_{k-1|k-1} F^T + Q \tag{4} \]

Update Step

\[ K_k = P_{k|k-1} H^T (H P_{k|k-1} H^T + R)^{-1} \tag{5} \] \[ \hat{x}_{k|k} = \hat{x}_{k|k-1} + K_k (y_k - H \hat{x}_{k|k-1}) \tag{6} \] \[ P_{k|k} = (I - K_k H) P_{k|k-1} \tag{7} \]

For smoothing, we must store \(\hat{x}_{k|k}, P_{k|k}, \hat{x}_{k|k-1}, P_{k|k-1}\) at every step. For nonlinear systems, replacing equations \((3)\) –\((7)\) with the corresponding EKF or UKF steps yields nonlinear smoothers within the same framework.

Fixed-Interval Smoother (RTS)

Fixed-interval smoothing computes the posterior \(p(x_k \mid y_{1:T})\) for every time step using the whole record. The RTS smoother solves this with a backward recursion from \(T\) down to \(0\) .

Smoother gain:

\[ G_k = P_{k|k} F^T P_{k+1|k}^{-1} \tag{8} \]

Smoothed state and covariance:

\[ \hat{x}_{k|T} = \hat{x}_{k|k} + G_k (\hat{x}_{k+1|T} - \hat{x}_{k+1|k}) \tag{9} \] \[ P_{k|T} = P_{k|k} + G_k (P_{k+1|T} - P_{k+1|k}) G_k^T \tag{10} \]

The recursion starts from \(\hat{x}_{T|T}, P_{T|T}\) (the final forward-pass result) and runs for \(k = T-1, \ldots, 0\) . Since \(P_{k+1|T} \preceq P_{k+1|k}\) , we always have \(P_{k|T} \preceq P_{k|k}\) : smoothed uncertainty never exceeds filtered uncertainty. See the RTS smoother article for a full derivation.

Fixed-Lag Smoother

A fixed-lag smoother outputs the smoothed estimate \(\hat{x}_{k-N|k}\) — the state \(N\) steps in the past — every time a new observation \(y_k\) arrives. You keep streaming, yet each output has effectively “seen” \(N\) steps of the future.

Implementation 1: State Augmentation

Stack the last \(N+1\) states into an augmented state vector

\[ z_k = \begin{bmatrix} x_k^T & x_{k-1}^T & \cdots & x_{k-N}^T \end{bmatrix}^T \tag{11} \]

and run a standard Kalman filter on the augmented system whose transition matrix shifts the blocks down by one (the Biswas–Mahalanabis fixed-lag smoother). Simple to implement, but the state dimension becomes \((N+1)n\) , so the cost grows as \(O((Nn)^3)\) .

Implementation 2: Windowed RTS (used below)

At each time \(k\) , run the RTS backward recursion (equations \((8)\) –\((9)\) ) only over the sliding window \([k-N, k]\) and output the estimate at the head of the window, \(\hat{x}_{k-N|k}\) . The per-step cost is \(O(N n^3)\) and only the last \(N\) forward-pass results need to be kept, so memory is \(O(N n^2)\) .

The Lag-Length vs Accuracy Trade-off

Corrections from future observations propagate backward through the smoother gain \(G_k\) . With the steady-state gain \(\bar{G}\) , the contribution of an observation \(N\) steps ahead decays roughly like \(\bar{G}^N\) . Consequently:

  • Increasing \(N\) moves the accuracy toward RTS (which corresponds to \(N = \infty\) ), but the improvement saturates exponentially
  • Latency and memory grow linearly in \(N\)

In practice, pick the smallest \(N\) for which \(\bar{G}^N\) is negligible — typically a few multiples of the system time constant. The experiment below shows this saturation numerically.

Fixed-Point Smoother

The fixed-point smoother refines the estimate of one particular time \(j\) (e.g., the initial state or an event time), updating \(\hat{x}_{j|k}\) whenever a new observation \(y_k\) with \(k > j\) arrives.

Define the product of smoother gains from \(j\) to \(k-1\) :

\[ B_k = \prod_{i=j}^{k-1} G_i = B_{k-1} G_{k-1}, \quad B_j = I \tag{12} \]

Then, with the filter innovation \(\nu_k = y_k - H \hat{x}_{k|k-1}\) and Kalman gain \(K_k\) :

\[ \hat{x}_{j|k} = \hat{x}_{j|k-1} + B_k K_k \nu_k \tag{13} \] \[ P_{j|k} = P_{j|k-1} + B_k (P_{k|k} - P_{k|k-1}) B_k^T \tag{14} \]

This is an \(O(n^3)\) -per-step incremental update with fixed \(O(n^2)\) memory, since only time \(j\) is tracked. The gain product \(B_k\) shrinks as \(k\) grows, so the updates fade out and \(\hat{x}_{j|k}\) converges to the fixed-interval estimate \(\hat{x}_{j|T}\) .

Python Implementation: RMSE Benchmark on a Random Walk

We compare the Kalman filter, fixed-lag smoother, and RTS smoother on the simplest possible system: a 1-D random walk observed in noise.

Model Definition

\[ x_k = x_{k-1} + w_{k-1}, \quad y_k = x_k + v_k, \quad w \sim \mathcal{N}(0, 0.05), \; v \sim \mathcal{N}(0, 1) \tag{15} \]
import numpy as np

# ---- State-space model: random walk + observation noise ----
T = 300  # number of time steps

F = np.array([[1.0]])   # state transition (random walk)
H = np.array([[1.0]])   # observation matrix
Q = np.array([[0.05]])  # process noise variance
R = np.array([[1.0]])   # observation noise variance

n, m = 1, 1
x0 = np.array([0.0])
P0 = np.diag([1.0])

Kalman Filter (Forward Pass)

def kalman_filter(y_obs, F, H, Q, R, x0, P0):
    """Kalman filter forward pass (stores all intermediates for smoothing)"""
    T = len(y_obs)
    n = len(x0)
    x_filt = np.zeros((T + 1, n))       # x_{k|k}
    P_filt = np.zeros((T + 1, n, n))    # P_{k|k}
    x_pred = np.zeros((T, n))           # x_{k|k-1}
    P_pred = np.zeros((T, n, n))        # P_{k|k-1}
    x_filt[0], P_filt[0] = x0, P0
    for k in range(T):
        # Prediction (eqs. 3, 4)
        x_pred[k] = F @ x_filt[k]
        P_pred[k] = F @ P_filt[k] @ F.T + Q
        # Update (eqs. 5, 6, 7)
        S = H @ P_pred[k] @ H.T + R
        K = P_pred[k] @ H.T @ np.linalg.inv(S)
        x_filt[k + 1] = x_pred[k] + K @ (y_obs[k] - H @ x_pred[k])
        P_filt[k + 1] = (np.eye(n) - K @ H) @ P_pred[k]
    return x_filt, P_filt, x_pred, P_pred

RTS Smoother (Backward Pass)

def rts_smoother(x_filt, P_filt, x_pred, P_pred, F):
    """RTS (fixed-interval) smoother backward pass"""
    T = len(x_pred)
    n = x_filt.shape[1]
    x_smooth = np.zeros((T + 1, n))
    P_smooth = np.zeros((T + 1, n, n))
    x_smooth[T], P_smooth[T] = x_filt[T], P_filt[T]
    for k in range(T - 1, -1, -1):
        G = P_filt[k] @ F.T @ np.linalg.inv(P_pred[k])  # eq. 8
        x_smooth[k] = x_filt[k] + G @ (x_smooth[k + 1] - x_pred[k])  # eq. 9
        P_smooth[k] = P_filt[k] + G @ (P_smooth[k + 1] - P_pred[k]) @ G.T  # eq. 10
    return x_smooth, P_smooth

Fixed-Lag Smoother (Windowed RTS)

def fixed_lag_smoother(y_obs, F, H, Q, R, x0, P0, N):
    """Fixed-lag smoother: at each time k, run the backward recursion over
    the last N steps only and output x_{k-N|k} (windowed-RTS implementation)"""
    T = len(y_obs)
    n = len(x0)
    x_filt, P_filt, x_pred, P_pred = kalman_filter(y_obs, F, H, Q, R, x0, P0)
    x_lag = np.zeros((T + 1, n))
    x_lag[0] = x0
    for k in range(1, T + 1):
        j0 = max(0, k - N)  # time index k-N to output
        xs = x_filt[k].copy()
        for j in range(k - 1, j0 - 1, -1):  # backward recursion over [j0, k]
            G = P_filt[j] @ F.T @ np.linalg.inv(P_pred[j])
            xs = x_filt[j] + G @ (xs - x_pred[j])
        x_lag[j0] = xs
    # The last N points have no future observations yet; fill them with
    # the backward recursion of the final window
    xs = x_filt[T].copy()
    x_lag[T] = xs
    for j in range(T - 1, max(0, T - N) - 1, -1):
        G = P_filt[j] @ F.T @ np.linalg.inv(P_pred[j])
        xs = x_filt[j] + G @ (xs - x_pred[j])
        x_lag[j] = xs
    return x_lag

Monte Carlo Evaluation (50 Runs)

A single realization can shuffle the ranking due to sampling noise, so we generate 50 sequences with numpy.random.default_rng and compare the average RMSE.

def rmse(a, b):
    return np.sqrt(np.mean((a - b) ** 2))

n_trials = 50
lags = [1, 2, 3, 5, 10, 20]
rmse_kf_list, rmse_rts_list = [], []
rmse_lag_lists = {N: [] for N in lags}

rng = np.random.default_rng(42)
for trial in range(n_trials):
    # Generate the true states and observations
    x_true = np.zeros((T + 1, n))
    y_obs = np.zeros((T, m))
    for k in range(T):
        x_true[k + 1] = F @ x_true[k] + rng.multivariate_normal(np.zeros(n), Q)
        y_obs[k] = H @ x_true[k + 1] + rng.multivariate_normal(np.zeros(m), R)

    x_filt, P_filt, x_pred, P_pred = kalman_filter(y_obs, F, H, Q, R, x0, P0)
    x_rts, P_rts = rts_smoother(x_filt, P_filt, x_pred, P_pred, F)

    rmse_kf_list.append(rmse(x_true[1:, 0], x_filt[1:, 0]))
    rmse_rts_list.append(rmse(x_true[1:, 0], x_rts[1:, 0]))
    for N in lags:
        x_lag = fixed_lag_smoother(y_obs, F, H, Q, R, x0, P0, N)
        rmse_lag_lists[N].append(rmse(x_true[1:, 0], x_lag[1:, 0]))

rmse_kf = np.mean(rmse_kf_list)
rmse_rts = np.mean(rmse_rts_list)
print(f"KF  filter RMSE : {rmse_kf:.4f}")
for N in lags:
    r = np.mean(rmse_lag_lists[N])
    print(f"Fixed-lag N={N:2d}   : {r:.4f}  ({(1 - r / rmse_kf) * 100:.1f}% improvement)")
print(f"RTS smoother    : {rmse_rts:.4f}  ({(1 - rmse_rts / rmse_kf) * 100:.1f}% improvement)")

Results

KF  filter RMSE : 0.4519
Fixed-lag N= 1   : 0.4145  (8.3% improvement)
Fixed-lag N= 2   : 0.3890  (13.9% improvement)
Fixed-lag N= 3   : 0.3717  (17.7% improvement)
Fixed-lag N= 5   : 0.3521  (22.1% improvement)
Fixed-lag N=10   : 0.3397  (24.8% improvement)
Fixed-lag N=20   : 0.3387  (25.0% improvement)
RTS smoother    : 0.3388  (25.0% improvement)
MethodMean RMSEImprovement vs filterLatency
Kalman filter0.45190
Fixed-lag \(N=1\)0.41458.3%1 step
Fixed-lag \(N=2\)0.389013.9%2 steps
Fixed-lag \(N=3\)0.371717.7%3 steps
Fixed-lag \(N=5\)0.352122.1%5 steps
Fixed-lag \(N=10\)0.339724.8%10 steps
Fixed-lag \(N=20\)0.338725.0%20 steps
RTS (fixed-interval)0.338825.0%full record

As expected, accuracy improves monotonically filter < fixed-lag < RTS, and the fixed-lag smoother converges to RTS as \(N\) grows (at \(N=20\) the difference is within Monte Carlo error).

For this system the steady-state prediction variance is \(\bar{P}_{k|k-1} = 0.25\) and the steady-state smoother gain is \(\bar{G} = \bar{P}_{k|k} / \bar{P}_{k+1|k} = 0.8\) . Since the contribution of future observations decays like \(\bar{G}^N = 0.8^N\) , we get \(0.8^{10} \approx 0.11\) and \(0.8^{20} \approx 0.012\) — consistent with the observed behavior of “nearly saturated at \(N=10\) , indistinguishable from RTS at \(N=20\) .” A practical design rule: choose the lag from your latency budget, using the \(N\) at which \(\bar{G}^N\) drops below a few percent as an upper bound.

Beyond Linear Gaussian Models

This article was restricted to linear Gaussian systems, but the same three-way classification carries over directly to nonlinear and non-Gaussian settings. Besides EKF/UKF-based smoothers, particle smoothers (forward-filtering backward-smoothing, \(O(T M^2)\) with \(M\) particles) re-weight the samples of a particle filter backward in time, enabling smoothing of multimodal posteriors. More recently, deep state space models (Deep Kalman Filters, the S4/Mamba family) learn the transition and observation models with neural networks, and the same filtering-vs-smoothing distinction reappears in their variational posterior approximations. The latency and complexity intuitions of classical smoothers transfer directly to these modern methods.

Summary

  • Smoothing is non-causal estimation that uses future observations: it never does worse than filtering, but it introduces latency
  • RTS (fixed-interval) is the offline gold standard for all time steps; fixed-lag trades a delay of \(N\) steps for near-RTS accuracy in streaming settings; fixed-point incrementally refines one specific time step at minimal cost
  • A 50-run Monte Carlo on a random walk confirmed the monotone improvement: filter RMSE 0.4519 → fixed-lag \(N=5\) 0.3521 → RTS 0.3388
  • The benefit of the lag saturates exponentially as \(\bar{G}^N\) , so a lag of a few system time constants is enough in practice

References

  • Rauch, H. E., Tung, F., & Striebel, C. T. (1965). “Maximum likelihood estimates of linear dynamic systems.” AIAA Journal, 3(8), 1445-1450.
  • Biswas, K. K., & Mahalanabis, A. K. (1972). “An approach to fixed-lag smoothing problems.” IEEE Transactions on Aerospace and Electronic Systems, AES-8(5), 676-682.
  • Meditch, J. S. (1969). Stochastic Optimal Linear Estimation and Control. McGraw-Hill.
  • Särkkä, S. (2013). Bayesian Filtering and Smoothing. Cambridge University Press.