Efficient sampling approaches based on generalized Golub-Kahan methods for large-scale hierarchical Bayesian inverse problems
Elle Buser, Julianne Chung
TL;DR
This work tackles uncertainty quantification for large-scale hierarchical Bayesian inverse problems by developing Metropolis-Hastings independence samplers within a Gibbs framework. It introduces two genGK-based proposal strategies: a low-rank genGK approximation to the conditional covariance and a preconditioned Lanczos method for exact square-root sampling of the conditional covariance, both designed to circumvent expensive high-dimensional Gaussian draws when the prior covariance %5Cmathbf{Q}%5C is large or ill-conditioned. The authors provide explicit formulations for acceptance ratios and show how genGK subspaces can be reused across hyperparameters, enabling scalable sampling across dynamic, seismic, and atmospheric inverse problems. Numerical results on seismic tomography, atmospheric transport, and dynamic photoacoustic tomography demonstrate competitive acceptance rates, effective sample sizes, and good equilibrium behavior, validating the methods’ applicability to real-world, large-scale UQ tasks.
Abstract
Uncertainty quantification for large-scale inverse problems remains a challenging task. For linear inverse problems with additive Gaussian noise and Gaussian priors, the posterior is Gaussian but sampling can be challenging, especially for problems with a very large number of unknown parameters (e.g., dynamic inverse problems) and for problems where computation of the square root and inverse of the prior covariance matrix are not feasible. Moreover, for hierarchical problems where several hyperparameters that define the prior and the noise model must be estimated from the data, the posterior distribution may no longer be Gaussian, even if the forward operator is linear. Performing large-scale uncertainty quantification for these hierarchical settings requires new computational techniques. In this work, we consider a hierarchical Bayesian framework where both the noise and prior variance are modeled as hyperparameters. Our approach uses Metropolis-Hastings independence sampling within Gibbs where the proposal distribution is based on generalized Golub-Kahan based methods. We consider two proposal samplers, one that uses a low rank approximation to the conditional covariance matrix and another that uses a preconditioned Lanczos method. Numerical examples from seismic imaging, dynamic photoacoustic tomography, and atmospheric inverse modeling demonstrate the effectiveness of the described approaches.