Evaluating Data Assimilation Algorithms

Data assimilation leads naturally to a Bayesian formulation in which the posterior probability distribution of the system state, given the observations, plays a central conceptual role. The aim of this paper is to use this Bayesian posterior probability distribution as a gold standard against which to evaluate various commonly used data assimilation algorithms. A key aspect of geophysical data assimilation is the high dimensionality and low predictability of the computational model. With this in mind, yet with the goal of allowing an explicit and accurate computation of the posterior distribution, we study the 2D Navier-Stokes equations in a periodic geometry. We compute the posterior probability distribution by state-of-the-art statistical sampling techniques. The commonly used algorithms that we evaluate against this accurate gold standard, as quantified by comparing the relative error in reproducing its moments, are 4DVAR and a variety of sequential filtering approximations based on 3DVAR and on extended and ensemble Kalman filters. The primary conclusions are that: (i) with appropriate parameter choices, approximate filters can perform well in reproducing the mean of the desired probability distribution; (ii) however they typically perform poorly when attempting to reproduce the covariance; (iii) this poor performance is compounded by the need to modify the covariance, in order to induce stability. Thus, whilst filters can be a useful tool in predicting mean behavior, they should be viewed with caution as predictors of uncertainty. These conclusions are intrinsic to the algorithms and will not change if the model complexity is increased, for example by employing a smaller viscosity, or by using a detailed NWP model.

💡 Research Summary

The paper presents a rigorous benchmark for evaluating popular data assimilation algorithms by using the full Bayesian posterior distribution as a gold‑standard reference. Because computing the exact posterior is infeasible for realistic atmospheric models, the authors select the two‑dimensional incompressible Navier‑Stokes equations on a periodic domain as a testbed. This system is high‑dimensional (≈10³ degrees of freedom) and exhibits chaotic dynamics, yet it remains tractable for state‑of‑the‑art Markov chain Monte Carlo (MCMC) sampling.

First, the dynamical model is described together with the observation operator (linear, noisy, regularly spaced in time). A Gaussian prior is assumed, and synthetic observations are generated. Using a pre‑conditioned Hamiltonian Monte Carlo sampler with multiple parallel chains, the authors obtain a converged ensemble of posterior samples. From this ensemble they compute the true posterior mean, covariance, and higher‑order moments, which serve as the “gold‑standard” against which all other methods are judged.

Second, four representative assimilation schemes are implemented: (i) four‑dimensional variational assimilation (4D‑VAR), (ii) three‑dimensional variational assimilation (3D‑VAR) with a fixed background covariance, (iii) the extended Kalman filter (EKF), and (iv) several variants of the ensemble Kalman filter (EnKF), including stochastic EnKF, deterministic EnKF, and localized EnKF. All methods are supplied with the same observation sequence and the same prior information. For each method the authors tune algorithmic parameters (observation error variance, covariance inflation factor, localization radius, ensemble size) to achieve the best possible performance.

Third, performance is quantified by the relative error of the estimated mean and covariance with respect to the Bayesian gold‑standard. The results show a clear dichotomy. All methods are capable of reproducing the posterior mean with reasonable accuracy when parameters are well‑chosen; 4D‑VAR and EnKF typically achieve mean errors below 5 %. However, the reconstruction of the posterior covariance is dramatically poorer. 3D‑VAR, which uses a static background covariance, fails to capture any time‑varying uncertainty. EKF suffers from linearisation errors and numerical instability in high dimensions, often producing non‑positive‑definite covariances. EnKF, while conceptually able to estimate flow‑dependent covariances, requires a finite ensemble; with realistic ensemble sizes (tens to a few hundred) the sampled covariance is severely under‑dispersed. Inflation and localization improve stability and reduce mean‑square error, but they also distort the true uncertainty structure, leading to relative covariance errors ranging from 30 % to over 200 %.

Fourth, the authors explore the effect of increasing model complexity (reducing viscosity, refining the spatial grid) and find that the qualitative conclusions are unchanged. The inability of approximate filters to capture covariance information is intrinsic to the Gaussian/linear approximations they employ, not a consequence of the specific dynamical system.

The paper also discusses theoretical connections: 4D‑VAR can be interpreted as a maximum‑a‑posteriori (MAP) estimator, forward‑backward filtering corresponds to smoothing, and under linear‑Gaussian assumptions the Kalman filter provides the exact posterior. In practice, however, the approximations required for tractability break these equivalences.

In summary, the study provides strong evidence that while variational and ensemble‑based filters are useful tools for estimating the most likely state of a chaotic system, they must be used with caution when quantifying uncertainty. The authors recommend that operational centres complement mean‑state estimates with additional diagnostics (e.g., ensemble spread monitoring, rank histograms) and pursue research into non‑Gaussian, particle‑based, or hybrid variational‑filtering methods that can better approximate the full Bayesian posterior.