Multilevel approximation of Gaussian random fields: Covariance compression, estimation and spatial prediction

Abstract

Centered Gaussian random fields (GRFs) indexed by compacta such as smooth, bounded domains in Euclidean space or smooth, compact and orientable manifolds are determined by their covariance operators. We consider centered GRFs given sample-wise as variational solutions to coloring operator equations driven by spatial white noise, with pseudodifferential coloring operator being elliptic, self-adjoint and positive from the Hörmander class. This includes the Matérn class of GRFs as a special case. Using microlocal tools and biorthogonal multiresolution analyses on the manifold, we prove that the precision and covariance operators, respectively, may be identified with bi-infinite matrices and finite sections may be diagonally preconditioned rendering the condition number independent of the dimension p of this section. We prove that a tapering strategy by thresholding as e.g. in [Bickel, P.J. and Levina, E. Covariance regularization by thresholding, Ann. Statist., 36 (2008), 2577-2604] applied on finite sections of the bi-infinite precision and covariance matrices results in optimally numerically sparse approximations. Numerical sparsity signifies that only asymptotically linearly many nonzero matrix entries are sufficient to approximate the original section of the bi-infinite covariance or precision matrix using this tapering strategy to arbitrary precision. This tapering strategy is non-adaptive and the locations of these nonzero matrix entries are known a priori. The tapered covariance or precision matrices may also be optimally diagonal preconditioned. Analysis of the relative size of the entries of the tapered covariance matrices motivates novel, multilevel Monte Carlo (MLMC) oracles for covariance estimation, in sample complexity that scales log-linearly with respect to the number p of parameters. This extends [Bickel, P.J. and Levina, E. Regularized Estimation of Large Covariance Matrices, Ann. Stat., 36 (2008), pp. 199-227] to estimation of (finite sections of) pseudodifferential covariances for GRFs by this fast MLMC method. Assuming at hand sections of the bi-infinite covariance matrix in wavelet coordinates, we propose and analyze a novel compressive algorithm for simulating and kriging of GRFs. The complexity (work and memory vs. accuracy) of these three algorithms scales near-optimally in terms of the number of parameters p of the sample-wise approximation of the GRF in Sobolev scales.

The SPDE approach for Gaussian and non-Gaussian fields: 10 years and still running

“…Considering the generalised Whittle-Matérn fields (7) for a convex polytope D ⊂ R d , d ∈ {1, 2, 3}, we know from that there exists a unique solution to the SPDE given that α > d/2 (which corresponds to ν > 0 in the stationary case), κ is an essentially bounded function, κ ∈ L ∞ (D), and H is a sufficiently nice (Lipschitz continuous on D and uniformly positive definite) matrix valued function. Herrmann et al (2020) extended this existence result by showing that it holds also when considering the SPDE (7) on a closed, connected, orientable, smooth, compact 2-surface in R 3 , under the assumption that κ, H are smooth, and Harbrecht et al (2021) derived similar results for the model on more general manifolds without boundaries. Cox and Kirchner (2020) generalized the case D ⊂ R d further by only requiring that the domain D has a Lipschitz boundary, and by relaxing the requirement on H to only assume essential boundedness and uniformly positive definiteness.…”

“…The spectral representations are linked to the eigenfunctions and eigenvalues of the Laplacian and its manifold versions. The methods for fractional operators discussed in Section 4.1 have recently been extended using high order numerical methods for PDEs on Riemannian manifolds by Lang and Pereira (2021) and Harbrecht et al (2021), involving polynomial and wavelet basis expansions. Other approaches to constructing valid models on spheres can be found in Porcu et al (2016).…”

Galerkin--Chebyshev approximation of Gaussian random fields on compact Riemannian manifolds

“…A quadrature approximation allowed them to derive numerical approximations of such fields defined on bounded Euclidean domains (Bolin et al, 2018a,b) and even compact metric spaces (Herrmann et al, 2020). This approach requires to solve multiple (large but sparse) linear systems in order to generate samples of the random fields, and work has been done to find suitable and efficient preconditionners to tackle them (Harbrecht et al, 2021). In comparison, our approach does not assume the existence of a SPDE defining the random field but still includes Whittle-Matérn fields as a particular case.…”

“…As such, Z can be seen as an instance of a Whittle-Matérn random field on a manifold, as popularized by Lindgren et al (2011) for compact Riemannian manifolds. This class of random fields were studied by Jansson et al (2021) for the particular case where the manifold is a sphere, and for compact Riemannian manifolds by Herrmann et al (2020) and Harbrecht et al (2021).…”

Regularity theory for a new class of fractional parabolic stochastic evolution equations

“…This relation gave rise to the SPDE approach proposed by Lindgren, Rue, and Lindström [49], where the SPDE (1.2) is considered on a bounded domain D R d and augmented with Dirichlet or Neumann boundary conditions. Besides enabling the applicability of efficient numerical methods available for (S)PDEs, such as finite element methods [12,14,15,22,40,49] or wavelets [17,39], this approach has the advantage of allowing for (a) nonstationary or anisotropic generalizations by replacing the operator κ 2 −∆ in (1.2) with more general strongly elliptic second-order differential operators,…”