Posterior contraction rates for non-parametric state and drift estimation

We consider a combined state and drift estimation problem for the linear stochastic heat equation. The infinite-dimensional Bayesian inference problem is formulated in terms of the Kalman-Bucy filter over an extended state space, and its long-time asymptotic properties are studied. Asymptotic posterior contraction rates in the unknown drift function are the main contribution of this paper. Such rates have been studied before for stationary non-parametric Bayesian inverse problems, and here we demonstrate the consistency of our time-dependent formulation with these previous results building upon scale separation and a slow manifold approximation.


Introduction
In this paper, we consider the combined state and drift estimation problem for the stochastic heat equation (1) du t = ∂ 2 x u t dt + f * dt + γ 1/2 dw t over the domain Ω = (0, π) and for t ≥ 0, with zero Dirichlet boundary conditions u t (0) = u t (π) = 0 and zero initial condition u 0 (x) = 0 for all x ∈ Ω, without loss of generality. We denote by w t a cylindrical Wiener process where γ > 0 characterises the strength of the model error. The unknown drift function f * is assumed to belong to a Sobolev space of appropriate regularity. We specify the properties of w t and f * more precisely in Section 2 below. Linear stochastic partial differential equations (SPDEs) of the form (1) are well understood from analytical and numerical perspectives. See, for example, Da Prato & Zabczyk (1992), Hairer (2009), Lord et al. (2014). In this paper, we specifically focus on the estimation of the states u t , t > 0, as well as the unknown drift function f * from noisy state measurements given by (2) dy t = u t dt + ρ 1/2 dv t , where y 0 (x) = 0 for all x ∈ Ω, ρ > 0 denotes the strength of the measurement noise, and v t is another cylindrical Wiener process independent of the model error w t . Our focus is on asymptotic posterior contraction rates, which have been widely studied in the context of non-parametric stationary inverse problems. See, for example, Knapik et al. (2011), Giné & Nickl (2016), Ghosal & van der Vaart (2017). We note that pure parameter estimation problems for SPDEs given exact observations of the Date: March 23, 2020. Institute of Mathematics, University of Potsdam, Karl-Liebknecht-Str. 24/25, D-14476 Potsdam, Germany, e-mail: sebastian.reich@uni-potsdam.de. states are also well-studied. See, for example, Cialenco (2018) for a recent survey. A first step towards non-parametric time-dependent inverse problems has been taken in Yan (2019). However, there they are framed as inference problems over a fixed time interval resulting in smoothing problems, whereas filtering problems are the focus of this paper. Furthermore, state and drift estimation are treated separately. In addition, parameter estimation for a class of linear SPDEs from noiseless local state measurements has been considered in Altmeyer & Reiß (2019), whereas we assume noisy observations of the full states throughout this paper.
In this paper, the infinite-dimensional, time-dependent combined state and drift estimation problem is formulated and analysed in terms of the well-known Kalman-Bucy filter equations for the posterior mean and covariance operator in an augmented state space. See, for example, Jazwinski (1970), Simon (2006), Curtain (1975) for an introduction to the Kalman-Bucy filter in finite and infinite dimensions. More recently, the Kalman-Bucy filter has been reformulated as a set of mean-field equations termed the ensemble Kalman-Bucy filter. See, for example Bergemann & Reich (2012), Wiljes et al. (2018), Nüsken et al. (2019). These mean-field equations allow for a concise formulation of our time-dependent estimation problem. Furthermore, although not considered in this paper, the Kalman-Bucy mean-field equations can be generalised to non-Gaussian estimation problems (Yang et al. 2013, Taghvaei et al. 2017. A key observation arising from the analysis of the Kalman-Bucy mean-field equations is the time-scale separation in the dynamics of the state and parameter variables, allowing us to apply the concept of slow manifolds (Verhulst 2007).
The remainder of this paper is structured as follows. A mathematical formulation of the combined state and drift estimation problem in terms of Fourier modes is provided in Section 2. Furthermore, the time-dependent estimation problem is formulated in terms of Kalman-Bucy mean-field equations. This formulation is applied to a stationary drift estimation problem in Section 3. Asymptotic posterior contraction rates are then derived and compared to existing results for a closely related non-parametric Bayesian inference problem. These results are extended to the combined state and drift estimation problem in Section 4. The essential analysis of the single-mode Kalman-Bucy filter equations is carried out in Section 4.1. A numerical exploration of the combined state and parameter estimation problem is carried out in Section 5. The paper concludes with Section 6.

Mathematical problem formulation
In this paper, we analyse the state and drift estimation problem using a spectral representation in terms of Fourier modes. We provide the essential background in this section. It is well known that solutions of (1) can be expanded in Fourier modes {U t (k)} k≥1 , that is, Because of (1), the Fourier modes obey the stochastic differential equations (SDEs) Here F * (k), k ≥ 1, denote the Fourier coefficients of f * , and W t (k) independent standard Brownian motions. The initial conditions are U 0 (k) = 0 for all k ≥ 1. Therefore (1) can be viewed as a stochastic evolution equation on the separable Hilbert space of square integrable functions u : Ω → R of sufficient spatial regularity with zero boundary conditions. The measurement model (2) can also be transformed into Fourier space, yielding with V t (k) representing independent standard Brownian motions, independent of the model errors W t (k ′ ) for all k ′ . Throughout this paper, we will exclusively work with the formulation (4) and (5) of our state and drift estimation problem in Fourier space. We now proceed to formulate a continuous-time Bayesian inference framework for this problem, which is used to derive asymptotic posterior contraction rates in Section 4. The posterior Bayesian random variables for the state and drift at time t, conditioned on the observed data Y (0,t] (k), k ≥ 1, are denoted by U t (k) and F t (k), respectively. Here Y (0,t] (k) represents all the data up to time t > 0, where is an observation at time s ∈ (0, t]. The mean-field Kalman-Bucy equations (Wiljes et al. 2018, Nüsken et al. 2019) for the combined state and drift estimation problem in the variables (U t (k), F t (k)) are then given by with the innovation I t (k) given by and Kalman gains Denoting the distribution of the joint random variable where σ u t (k) denotes the variance of U t (k). The prior assumptions are that U 0 (k) = 0 almost surely, and that the F 0 (k), k ≥ 1, are independent Gaussian random variables with mean zero and variance (4) and (5) are linear in the unknowns, U t (k) and F t (k) remain Gaussian under the evolution equations (7), and we investigate their behaviour in terms of their mean and covariance for t ≫ 1. We further assume that the true drift function f * has Sobolev regularity β > 0, that is for any δ > 0 and any sequence of coefficients c k with bounded |c k | as k → ∞. These coefficients can, for example, be realisations of i.i.d. uniform random variables from the interval [−1, 1].
In addition to an analysis of the continuous-time Bayesian inference problem (7), we also study the (linear) dependence of the estimators (Bayesian posterior means) , on the random measurement process Y (0,t] in Section 4. We denote the expectation values of (17) with respect to the observation process . These two quantities characterise the systematic bias in the estimators (17), while the implied variances are a measure of the 'frequentist' uncertainty of the estimators (17). According to Lemma 8.2 from Ghosal & van der Vaart (2017), the posterior contraction rate ε t (k) in the kth Fourier mode of the drift estimator F t (k), that is, The main contributions from our subsequent analysis consist of providing bounds for the three quantities σ f t (k), p f t (k), and m f t (k) − F * in terms of the statistical models (4) and (5), and additionally studying the asymptotic behaviour of the l 2 asymptotic contraction rate ε t , defined by in terms of the Sobolev regularity β of the true drift function f * and the variance of the prior characterised by α > 0 in (14). Here, the l 2 asymptotic contraction rate ε t is to be understood in the sense that where f t denotes the inverse Fourier transform of F t . We first discuss a stationary formulation of the drift estimation problem in Section 3, which is closely linked to available results for non-parametric Bayesian inverse problems (Knapik et al. 2011, Giné & Nickl 2016, Ghosal & van der Vaart 2017. The combined state and drift estimation problem is analysed in Section 4. Remark 2.1. A more general class of time-scaled prior variances over the unknown set of static parameters F (k) has been considered in Knapik et al. (2011). This more general class of prior distributions could be incorporated into our analysis as well. However, we restrict the subsequent analysis to the simpler case (14).

Stationary drift estimation problem
Before addressing the full time-dependent problem, we first consider the stationary drift estimation problem in which γ = 0 and dU t (k) = 0 in (4), giving rise to the forward model and observations for t ≥ 0. The forward model (26) leads to a mildly ill-posed inverse problem. We follow a Bayesian approach by placing a Gaussian mean-zero prior with variance given by (14) on the Fourier coefficients of the unknown drift function. Let us summarise our assumptions for ease of reference in the following: Assumption 3.1. We assume that the Fourier coefficients of the true drift function satisfy (15). The prior coefficients, denoted by F 0 (k), are independent Gaussian random variables with mean zero, that is, F 0 (k) = 0, and with variance σ f 0 (k) satisfying (14). Let us denote the Bayesian posterior estimate of the drift at time t > 0 by F t (k), k ≥ 1. We recall that F t (k) is Gaussian for all t > 0 provided F 0 (k) is Gaussian. Additionally, recall that we denote the mean of F t (k) by F t (k) and the variance by σ f t (k). The Kalman-Bucy evolution equations (7) reduce to the following equations for the mean and the variance of the kth Fourier mode: where K f t (k) denotes the Kalman gain factor The initial conditions are provided by Assumption 3.1.
Remark 3.2. Alternatively, the Kalman-Bucy equations (27) for the mean and the variance can be reformulated as mean-field equations in F t (k) directly, that is, These equations are a special case of (7), and provide the starting point for numerical implementation in the form of the ensemble Kalman-Bucy filter (Bergemann & Reich 2012), as well as extensions to nonlinear and non-Gaussian Bayesian inference problems (Yang et al. 2013, Taghvaei et al. 2017).
The evolution equation (27b) for the posterior variance has the closed form solution The trace of the covariance of the full joint posterior process thus satisfies the asymptotic estimate with the initial variances satisfying (14). This result follows from asymptotic bounds for infinite sums. See Lemma K.7 in Ghosal & van der Vaart (2017) in particular. Next, we carry out a 'frequentist' analysis of the mean F t (k) in terms of its bias with respect to the true F * (k), and its variance p f t (k) with respect to the measurement noise V t (k). We rewrite (27a) as Denoting the expectation value of F t (k) with respect to the measurement noise V t (k) by m f t (k), we obtain with initial condition m f 0 (k) = 0. Equation (33) has the closed-form solution

POSTERIOR CONTRACTION RATES FOR NON-PARAMETRIC STATE AND DRIFT ESTIMATION 7
Hence, based on (16), the l 2 -norm of the frequentist bias satisfies the asymptotic estimate This result follows again from asymptotic bounds for infinite sums. We finally investigate the time evolution of the frequentist variance p f t (k). Starting from the stochastic differential equation The initial variance is p f 0 (k) = 0 for all k ≥ 1. In order to analyse the solution behaviour of (37), we introduce ∆p f t (k) := σ f t (k)− p f t (k) and find the associated evolution equation In fact, the explicit solution of (37) is The 'frequentist' uncertainty is characterised by , so the probability of the event has probability tending to zero for any nonnegative function M t with M t → ∞. Here N(m, σ) denotes the Gaussian distribution with mean m and variance σ.
Theorem 3.3. Under Assumption 3.1, the posterior contraction rate (22) is given by Proof. According to Lemma 8.2 from Ghosal & van der Vaart (2017), the posterior contraction rate ε t (k) in each Fourier mode is provided by Because of (31) and (35) together with p f t (k) < σ f t (k), the result (43) follows. Remark 3.4. The following white noise forward model has been investigated in Knapik et al. (2011): where Ξ n (k) are i.i.d. Gaussian random variables with mean zero and variance one. The associated inference problem corresponds to a sequence of measurements with measurement error decreasing as 1/n. In our continuous-time problem, we obtain the same asymptotic rates as in the case considered above with p = 2 under the formal equivalence t = n → ∞.

Time-dependent state and drift estimation
We now return to the full dynamic model (4) subject to observations (5). We primarily wish to estimate the drift function F * . However, because of the stochastic model errors in (4), we also need to estimate the states U * t . We start with a careful analysis of the single mode system for both small-and large-k Fourier modes.

4.1.
Analysis of the single-mode filtering problem. In this section, we conduct a careful analysis of the single-mode filtering and parameter estimation problem. We suppress the dependence on the mode number k, and introduce the parameter ǫ = k −2 . The signal process is therefore given by (46) dU t = −ǫ −1 U t dt + F * dt + γ 1/2 dW t with given initial condition U 0 = 0 almost surely. Observations of the process are given by (47) dY t = U t dt + ρ 1/2 dV t with Y 0 = 0. The complete observation path up to time t is denoted by Y (0,t] . Recall that the mean-field Kalman-Bucy equations are given by (7). We again drop the dependence on the Fourier mode number k in the subsequent analysis. The meanfield equations (7) give rise to the following evolution equations in the conditional means U t := π t [u] and F t := π t [f ]: The deviations ∆U t := U t − U t and ∆F t := F t − F t thus satisfy We can further decompose (48) by introducing the 'frequentist' expectation values m u t = E * [U t ] and m f t = E * [F t ] with respect to the observation process (47), and the deviations ∆U t := U t − m u t and ∆F t : and introducing the shorthand µ t = E * [U t ], one is left with the ordinary differential equations d dt for the mean values (m u t , m f t ), as well as for the deviations (∆U t , ∆F t ). The expectation value µ t satisfies the evolution equation Both equations follow from (46). We finally combine (52a) and (54) into a single equation for the new variable U t := ∆U t − ∆U t , and replace (52) by The time-dependent linear stochastic differential equations (55) can also be analysed in terms of the variances Lemma 4.1. It holds that and the fact that K f t ∼ c/t for t sufficiently large with an appropriate constant c > 0. In particular, Remark 4.2. We note that the combined state and drift estimation problem under the slow manifold approximation behaves exactly like the stationary drift estimation problem from Section 3. In particular, compare equations (30) and (34). The numerical experiments in Section 5 reveal that the slow manifold approximation already holds for rather small-k Fourier modes. (c) k = 4 (d) k = 8 Figure 1. We display the time evolution of the bias |m f t (k)−F * (k)| and the variances σ f t (k) and p f t (k) for increasing values of k. We compare the results from the combined state and parameter estimation problem (labelled 'dynamic') with those from the associated direct inference problem in the parameter alone (labelled 'stationary'). The delay in the onset of the asymptotic t −1 regime as k increases can be clearly seen as well as that p f t (k) < σ f t (k) for all t > 0. Furthermore, as theoretically predicted, the dynamic and stationary problem formulations behave almost identically in terms of their drift estimates. drift does not deteriorate the asymptotic contraction rates. In fact, the data does not affect the posterior uncertainty in the states as k → ∞, which is entirely determined by the equilibrium distribution of the associated Ornstein-Uhlenbeck process. We verify this behaviour through a simple numerical experiment in the following section.

Numerical exploration
We numerically implemented the evolution equations (51), (58), and (61) for the Bayesian means, the Bayesian variances, and the frequentist variances, respectively, for different Fourier modes k. We also implemented the corresponding equations from Section 3 for the pure drift estimation problem. Simulations were run with ρ = γ = 1. The true reference value was set to F * (k) = k −β−1/2 with β = 1/2, and we used α = 1/2 for the prior variance (14).
The results can be found in Figure 1. They reveal a very similar behaviour for the combined state and drift estimation and the pure drift estimation problems. The delay in the onset of the asymptotic t −1 regime as k increases can also be seen. Overall, this simple numerical experiment confirms our theoretical investigations with regard to the dynamical behaviour in each Fourier mode from Section 4.

Conclusions
We have provided an analysis of the infinite dimensional Kalman-Bucy filter meanfiled equations (7) for a combined state and drift estimation problem defined in spectral space by (4) and (5). The derived asymptotic posterior contraction rates in the unknown drift function f * from Theorem 4.5 agree with those derived for an associated stationary problem formulation in Theorem 3.3. These theoretical findings imply that the required additional estimation of the states does not lead to a deterioration of the contraction rates, which has been confirmed numerically in Section 5. An extension to non-Gaussian estimation problems for SPDEs and different types of observation processes can be envisioned through nonlinear extensions of the Kalman-Bucy mean-field equations. See, for example, Nüsken et al. (2019). Furthermore, our results can be combined with those from Götze et al. (2019) to study the coverage probabilities of Bayesian credible sets in a non-asymptotic regime. One could also discuss adaptive choices for the prior F 0 (k). See, for example, Knapik et al. (2016).