Sebastian Reich

Sebastian Reich

Many applications, such as intermittent data assimilation, lead to a recursive application of Bayesian inference within a Monte Carlo context. Popular data assimilation algorithms include sequential Monte Carlo methods and ensemble Kalman filters (EnKFs). These methods differ in the way Bayesian inference is implemented. Sequential Monte Carlo methods rely on importance sampling combined with a resampling step while EnKFs utilize a linear transformation of Monte Carlo samples based on the classic Kalman filter. While EnKFs have proven to be quite robust even for small ensemble sizes, they are not consistent since their derivation relies on a linear regression ansatz. In this paper, we propose another transform method, which does not rely on any a prior assumptions on the underlying prior and posterior distributions. The new method is based on solving an optimal transportation problem for discrete random variables.

Javier Amezcua, Kayo Ide, Eugenia Kalnay et al.

Two recent works have adapted the Kalman-Bucy filter into an ensemble setting. In the first formulation, BR10, the full ensemble is updated in the analysis step as the solution of single set of ODEs in pseudo-BGR09, the ensemble of perturbations is updated by the solution of an ordinary differential equation (ODE) in pseudo-time, while the mean is updated as in the standard KF. In the second formulation, BR10, the full ensemble is updated in the analysis step as the solution of single set of ODEs in pseudo-time. Neither requires matrix inversions except for the frequently diagonal observation error covariance. We analyze the behavior of the ODEs involved in these formulations. We demonstrate that they stiffen for large magnitudes of the ratio of background to observational error covariance, and that using the integration scheme proposed in both BGR09 and BR10 can lead to failure. An integration scheme that is both stable and is not computationally expensive is proposed. We develop transform-based alternatives for these Bucy-type approaches so that the integrations are computed in ensemble space where the variables are weights (of dimension equal to the ensemble size) rather than model variables. Finally, the performance of our ensemble transform Kalman-Bucy implementations is evaluated using three models: the 3-variable Lorenz 1963 model, the 40-variable Lorenz 1996 model, and a medium complexity atmospheric general circulation model (AGCM) known as SPEEDY. The results from all three models are encouraging and warrant further exploration of these assimilation techniques.

Nawinda Chustagulprom, Sebastian Reich, Maria Reinhardt

Data assimilation is the task to combine evolution models and observational data in order to produce reliable predictions. In this paper, we focus on ensemble-based recursive data assimilation problems. Our main contribution is a hybrid filter that allows one to adaptively bridge between ensemble Kalman and particle filters. While ensemble Kalman filters are robust and applicable to strongly nonlinear systems even with small and moderate ensemble sizes, particle filters are asymptotically consistent in the large ensemble size limit. We demonstrate numerically that our hybrid approach can improve the performance of both Kalman and particle filters at moderate ensemble sizes. We also show how to implement the concept of localization into a hybrid filter, which is key to its applicability to spatially extended systems.

Alastair Gregory, Colin Cotter, Sebastian Reich

This paper extends the Multilevel Monte Carlo variance reduction technique to nonlinear filtering. In particular, Multilevel Monte Carlo is applied to a certain variant of the particle filter, the Ensemble Transform Particle Filter. A key aspect is the use of optimal transport methods to re-establish correlation between coarse and fine ensembles after resampling; this controls the variance of the estimator. Numerical examples present a proof of concept of the effectiveness of the proposed method, demonstrating significant computational cost reductions (relative to the single-level ETPF counterpart) in the propagation of ensembles.

C. J. Cotter, D. A. Ham, C. C. Pain et al.

We introduce a new mixed discontinuous/continuous Galerkin finite element for solving the 2- and 3-dimensional wave equations and equations of incompressible flow. The element, which we refer to as P1dg-P2, uses discontinuous piecewise linear functions for velocity and continuous piecewise quadratic functions for pressure. The aim of introducing the mixed formulation is to produce a new flexible element choice for triangular and tetrahedral meshes which satisfies the LBB stability condition and hence has no spurious zero-energy modes. We illustrate this property with numerical integrations of the wave equation in two dimensions, an analysis of the resultant discrete Laplace operator in two and three dimensions, and a normal mode analysis of the semi-discrete wave equation in one dimension.

Colin J. Cotter, Sebastian Reich

This survey paper is written with the intention of giving a mathematical introduction to filtering techniques for intermittent data assimilation, and to survey some recent advances in the field. The paper is divided into three parts. The first part introduces Bayesian statistics and its application to statistical inference and estimation. Basic aspects of Markov processes, as they typically arise from scientific models in the form of stochastic differential and/or difference equations, are covered in the second part. The third and final part describes the filtering approach to estimation of model states by assimilation of observational data into scientific models. While most of the material is of survey type, very recent advances in the field of nonlinear data assimilation covered in this paper include a discussion of Bayesian inference in the context of optimal transportation and coupling of random variables, as well as a discussion of recent advances in ensemble transform filters. References and sources for further reading material will be listed at the end of each section.

Angwenyi David, Jana de Wiljes, Sebastian Reich

Simultaneous state and parameter estimation arises from various applicational areas but presents a major computational challenge. Most available Markov chain or sequential Monte Carlo techniques are applicable to relatively low dimensional problems only. Alternative methods, such as the ensemble Kalman filter or other ensemble transform filters have, on the other hand, been successfully applied to high dimensional state estimation problems. In this paper, we propose an extension of these techniques to high dimensional state space models which depend on a few unknown parameters. More specifically, we combine the ensemble Kalman-Bucy filter for the continuous-time filtering problem with a generalized ensemble transform particle filter for intermittent parameter updates. We demonstrate the performance of this two stage update filter for a wave equation with unknown wave velocity parameter.

Manfred Opper, Sebastian Reich

This paper outlines a novel extension of the classical Pontryagin minimum (maximum) principle to stochastic optimal control problems. Contrary to the well-known stochastic Pontryagin minimum principle involving forward-backward stochastic differential equations, the proposed formulation is deterministic and of mean-field type. We denote it by the McKean-Pontryagin minimum principle. The Hamiltonian structure of the proposed McKean-Pontryagin minimum principle is achieved via the introduction of a pair of auxiliary functions. A gauge freedom in the choice of one of these two functions can be used to decouple the forward and reverse time equations; hence simplifying the solution of the underlying boundary value problem. We also consider infinite horizon discounted cost optimal control problems. In this case, the mean-field formulation allows one to convert the computation of the desired optimal control law into solving a pair of forward mean-field ordinary differential equations. The McKean-Pontryagin minimum principle is tested numerically for a controlled inverted pendulum, a controlled Lorenz-63 system, and a controlled Lorenz-96 system. Although the focus is on linear-quadratic control problems, the proposed methodology is extendable to more general problems including mean-field type control formulations.

Sebastian Reich

This note outlines a mean-field approach to dynamic optimal transport problems based on the recently proposed McKean-Pontryagin maximum principle. Key aspects of the proposed methodology include i) avoidance of sampling over stochastic paths, ii) a fully variational approach leading to constrained Hamiltonian equations of motion, and iii) a unified treatment of deterministic and stochastic optimal transport problems. We also discuss connections to well-known dynamic formulations in terms of forward-backward stochastic differential equations and extensions beyond classical entropic-regularized transport problems.

Zhang Jiang, Marios Andreou, Sebastian Reich et al.

Data assimilation (DA) integrates observational information with model predictions to improve state estimation in complex systems. While filtering provides the basis for online forecasts by using only past and present observations, it can exhibit delays and biases when the underlying dynamics evolve rapidly or undergo regime transitions. Smoothing, which additionally incorporates future observations, provides a natural pipeline for hindcasting and reanalysis that yields an uncertainty reduction beyond the filter. This paper introduces an ensemble Kalman-Bucy smoother (EnKBS) for continuous-time DA of nonlinear dynamical systems, where the smoother's conditional distributions are reconstructed using ensemble moments. The result is a derivative-free framework that does not require explicit computation of tangent-linear or adjoint models, which converges to the exact smoother solution at the infinite-ensemble limit for a wide class of complex systems. Incorporating standard regularization techniques for high-dimensional systems, such as covariance localization and inflation, the skill of the EnKBS is demonstrated in various important scientific problems. By integrating future observations, which reveal the underlying causal mechanisms for retrospective state updates, the EnKBS is used for Bayesian-based inference of causal relationships and their temporal influence range in a dyadic trigger-feedback model and the development of a causality-driven iterative learning algorithm that identifies the structure and recovers the hidden parameters of a nonlinear reduced-order model mimicking midlatitude atmospheric circulation. Notably, both tasks remain effective with an ensemble size of $O(10)$ under partial observations, suggesting that EnKBS can support the instantaneous discovery of high-dimensional complex systems over time.

James Matthew Young, Paula Cordero-Encinar, Sebastian Reich et al.

We develop a diffusion-based sampler for target distributions known up to a normalising constant. To this end, we rely on the well-known diffusion path that smoothly interpolates between a (simple) base distribution and the target distribution, widely used in diffusion models. Our approach is based on a practical implementation of diffusion-annealed Langevin Monte Carlo, which approximates the diffusion path with convergence guarantees. We tackle the score estimation problem by developing an efficient sequential Monte Carlo sampler that evolves auxiliary variables from conditional distributions along the path, which provides principled score estimates for time-varying distributions. We further develop novel control variate schedules that minimise the variance of these score estimates. Finally, we provide theoretical guarantees and empirically demonstrate the effectiveness of our method on several synthetic and real-world datasets.

Georg A. Gottwald, Sebastian Reich

We consider the problem of sampling from an unknown distribution for which only a sufficiently large number of training samples are available. In this paper, we build on previous work combining Schrödinger bridges and plug & play Langevin samplers. A key bottleneck of these approaches is the exponential dependence of the required training samples on the dimension, $d$, of the ambient state space. We propose a localization strategy which exploits conditional independence of conditional expectation values. Localization thus replaces a single high-dimensional Schrödinger bridge problem by $d$ low-dimensional Schrödinger bridge problems over the available training samples. In this context, a connection to multi-head self attention transformer architectures is established. As for the original Schrödinger bridge sampling approach, the localized sampler is stable and geometric ergodic. The sampler also naturally extends to conditional sampling and to Bayesian inference. We demonstrate the performance of our proposed scheme through experiments on a high-dimensional Gaussian problem, on a temporal stochastic process, and on a stochastic subgrid-scale parametrization conditional sampling problem. We also extend the idea of localization to plug & play Langevin samplers using kernel-based denoising in combination with Tweedie's formula.

Georg A. Gottwald, Sebastian Reich

We present a supervised learning method to learn the propagator map of a dynamical system from partial and noisy observations. In our computationally cheap and easy-to-implement framework a neural network consisting of random feature maps is trained sequentially by incoming observations within a data assimilation procedure. By employing Takens' embedding theorem, the network is trained on delay coordinates. We show that the combination of random feature maps and data assimilation, called RAFDA, outperforms standard random feature maps for which the dynamics is learned using batch data.

Georg A. Gottwald, Sebastian Reich

Data-driven prediction and physics-agnostic machine-learning methods have attracted increased interest in recent years achieving forecast horizons going well beyond those to be expected for chaotic dynamical systems. In a separate strand of research data-assimilation has been successfully used to optimally combine forecast models and their inherent uncertainty with incoming noisy observations. The key idea in our work here is to achieve increased forecast capabilities by judiciously combining machine-learning algorithms and data assimilation. We combine the physics-agnostic data-driven approach of random feature maps as a forecast model within an ensemble Kalman filter data assimilation procedure. The machine-learning model is learned sequentially by incorporating incoming noisy observations. We show that the obtained forecast model has remarkably good forecast skill while being computationally cheap once trained. Going beyond the task of forecasting, we show that our method can be used to generate reliable ensembles for probabilistic forecasting as well as to learn effective model closure in multi-scale systems.

Caroline L. Wormell, Sebastian Reich

Diffusion maps is a manifold learning algorithm widely used for dimensionality reduction. Using a sample from a distribution, it approximates the eigenvalues and eigenfunctions of associated Laplace-Beltrami operators. Theoretical bounds on the approximation error are however generally much weaker than the rates that are seen in practice. This paper uses new approaches to improve the error bounds in the model case where the distribution is supported on a hypertorus. For the data sampling (variance) component of the error we make spatially localised compact embedding estimates on certain Hardy spaces; we study the deterministic (bias) component as a perturbation of the Laplace-Beltrami operator's associated PDE, and apply relevant spectral stability results. Using these approaches, we match long-standing pointwise error bounds for both the spectral data and the norm convergence of the operator discretisation. We also introduce an alternative normalisation for diffusion maps based on Sinkhorn weights. This normalisation approximates a Langevin diffusion on the sample and yields a symmetric operator approximation. We prove that it has better convergence compared with the standard normalisation on flat domains, and present a highly efficient algorithm to compute the Sinkhorn weights.

Kay Bergemann, Sebastian Reich

Ensemble Kalman filter techniques are widely used to assimilate observations into dynamical models. The phase space dimension is typically much larger than the number of ensemble members which leads to inaccurate results in the computed covariance matrices. These inaccuracies can lead, among other things, to spurious long range correlations which can be eliminated by Schur-product-based localization techniques. In this paper, we propose a new technique for implementing such localization techniques within the class of ensemble transform/square root Kalman filters. Our approach relies on a continuous embedding of the Kalman filter update for the ensemble members, i.e., we state an ordinary differential equation (ODE) whose solutions, over a unit time interval, are equivalent to the Kalman filter update. The ODE formulation forms a gradient system with the observations as a cost functional. Besides localization, the new ODE ensemble formulation should also find useful applications in the context of nonlinear observation operators and observations arriving continuously in time.