Fast sampling of parameterised Gaussian random fields

TL;DR

This paper introduces a fast sampling method for parameterized Gaussian random fields using reduced basis surrogate models to efficiently handle multiple hyperparameter configurations, especially in complex spatial models.

Contribution

It develops a reduced basis approach for Karhunen-Loève expansions, enabling efficient sampling of Gaussian fields with unknown parameters in high-dimensional PDE problems.

Findings

Reduced eigenpair approximation errors are acceptable for practical use.

The method significantly reduces computational costs in Bayesian inverse problems.

Numerical experiments demonstrate effectiveness in 2D spatial models.

Abstract

Gaussian random fields are popular models for spatially varying uncertainties, arising for instance in geotechnical engineering, hydrology or image processing. A Gaussian random field is fully characterised by its mean function and covariance operator. In more complex models these can also be partially unknown. In this case we need to handle a family of Gaussian random fields indexed with hyperparameters. Sampling for a fixed configuration of hyperparameters is already very expensive due to the nonlocal nature of many classical covariance operators. Sampling from multiple configurations increases the total computational cost severely. In this report we employ parameterised Karhunen-Lo\`eve expansions for sampling. To reduce the cost we construct a reduced basis surrogate built from snapshots of Karhunen-Lo\`eve eigenvectors. In particular, we consider Mat\'ern-type covariance operators with unknown correlation length and standard deviation. We suggest a linearisation of the covariance function and describe the associated online-offline decomposition. In numerical experiments we investigate the approximation error of the reduced eigenpairs. As an application we consider forward uncertainty propagation and Bayesian inversion with an elliptic partial differential equation where the logarithm of the diffusion coefficient is a parameterised Gaussian random field. In the Bayesian inverse problem we employ Markov chain Monte Carlo on the reduced space to generate samples from the posterior measure. All numerical experiments are conducted in 2D physical space, with non-separable covariance operators, and finite element grids with $\sim 10^4$ degrees of freedom.