Fast Kalman Filtering and Forward–Backward Smoothing via a Low-Rank Perturbative Approach

Abstract

Kalman filtering-smoothing is a fundamental tool in statistical time-series analysis. However, standard implementations of the Kalman filter-smoother require O(d³) time and O(d²) space per time step, where d is the dimension of the state variable, and are therefore impractical in high-dimensional problems. In this article we note that if a relatively small number of observations are available per time step, the Kalman equations may be approximated in terms of a low-rank perturbation of the prior state covariance matrix in the absence of any observations. In many cases this approximation may be computed and updated very efficiently (often in just O(k²d) or O(k²d + kdlog d) time and space per time step, where k is the rank of the perturbation and in general k ≪ d), using fast methods from numerical linear algebra. We justify our approach and give bounds on the rank of the perturbation as a function of the desired accuracy. For the case of smoothing, we also quantify the error of our algorithm because of the low-rank approximation and show that it can be made arbitrarily low at the expense of a moderate computational cost. We describe applications involving smoothing of spatiotemporal neuroscience data. This article has online supplementary material.

Key Words:

• Covariance approximation
• Fast algorithm
• Low-rank methods
• Numerical analysis
• Tracking

View correction statement:

Correction

SUPPLEMENTARY MATERIAL

Appendix.pdf: Appendix including the omitted proofs and further technical results. Section A deals with extensions to non-Gaussian measurements. Section B provides further analysis and tighter bounds for the effective rank. Section C provides the proof for Theorem 4.1 and discusses the convergence properties of the LRBT algorithm when applied in an iterative fashion. Section D discusses how the LRBT can be extended to provide estimates of the smoothed covariance in O(K(d)T) time and O(dT) space.

Matlab code: Matlab code including routines for all the algorithms described in the paper (fast Kalman filter, LRBT, steepest descent LRBT, conjugate gradients with LRBT preconditioner). The wrapper code calls all the algorithms as well the exact Kalman filter and BT algorithms for the receptive field example presented in Section 5.1.

ACKNOWLEDGMENTS

LP is supported by a McKnight Scholar award, an NSF CAREER award, an NSF grant IIS-0904353, an NEI grant EY018003, and a DARPA contract N66001-11-1-4205. JH is supported by the Columbia College Rabi Scholars Program. The authors thank P. Jercog for kindly sharing his hippocampal data with us.

Notes

It is well known that the Woodbury formula can be numerically unstable when the observation covariance W is small (i.e., the high-SNR case). It should be possible to derive a low-rank square-root filter (Treebushny and Madsen 2005; Chandrasekar et al. 2008) to improve the numerical stability here, though we have not yet pursued this direction. Meanwhile, a crude but effective method to guarantee that C_t remains positive definite is to simply shrink Σ_t slightly if any negative eigenvalues are detected. This can be done easily in O(d) time by restricting attention to the subspace spanned by L_t.

Note that in a slight abuse of notation, we will recycle the names of some matrices (e.g., O_t and Q_t) that play a similar role in the LRBT approach as in the fast Kalman method described in the previous sections.