Efficient Marginalization-Based MCMC Methods for Hierarchical Bayesian Inverse Problems

Abstract

Hierarchical models in Bayesian inverse problems are characterized by an assumed prior probability distribution for the unknown state and measurement error precision, and hyperpriors for the prior parameters. Combining these probability models using Bayes' law often yields a posterior distribution that cannot be sampled from directly, even for a linear model with Gaussian measurement error and Gaussian prior, both of which we assume in this paper. In such cases, Gibbs sampling can be used to sample from the posterior [3], but problems arise when the dimension of the state is large. This is because the Gaussian sample required for each iteration can be prohibitively expensive to compute, and because the statistical efficiency of the Markov chain degrades as the dimension of the state increases. The latter problem can be mitigated using marginalization-based techniques, such as those found in [22,33,45], but these can be computationally prohibitive as well. In this paper, we combine the low-rank techniques of [8] with the marginalization approach of [45]. We consider two variants of this approach: delayed acceptance and pseudo-marginalization. We provide a detailed analysis of the acceptance rates and computational costs associated with our proposed algorithms, and compare their performances on two numerical test cases-image deblurring and inverse heat equation.

Computational methods for large-scale inverse problems: a survey on hybrid projection methods

“…In [23], the authors consider a fully Bayesian approach with assigned hyperpriors and use MCMC methods with Metropolis-Hastings independence sampling with a proposal distribution based on a low-rank approximation of the prior-preconditioned Hessian. In [188], these ideas are combined with marginalization. Theoretical results on approximating the posterior distribution and posterior covariance matrix can be found in [189,194].…”

Optimal experimental design for infinite-dimensional Bayesian inverse problems governed by PDEs: a review

“…Significant progress has been made in this area in recent years; the references [1,38,52,71,72,86,100,101] provide a small sample of literature from the last couple of decades. There have also been notable strides in computational methods for infinite-dimensional Bayesian inverse problems; see, e.g., [12,15,31,32,50,62,85,93,98].…”

Optimization-Based MCMC Methods for Nonlinear Hierarchical Statistical Inverse Problems

“…For linear inverse problems, [5] investigated the use of Gibbs sampling schemes that alternatively update the model parameters and hyperparameters, [1] analyzed the dimension scalability (w.r.t. the model parameters) of several Gibbs sampling schemes, [22] analyzed the consistency of the hyperparameter estimation, [25,49] investigated the use of the one-block-update of [48] and marginalization over model parameters to accelerate the sampling. The success of these developments commonly relies on two facts: there exists an analytic expression for the marginal posterior over the hyperparameters and one can the directly sample the conditional posterior over the model parameters for given hyperparameters.…”

“…This way, the RTO-PM is equivalent to a MH method that uses a proposal p RTO (u ♯ |θ ♯ ) q(θ ♯ , θ) to sample the joint posterior p(u, θ|y). Thus, for K = 1, the RTO-PM can also be viewed as the nonlinear extension of the one-block-update (see [25,48,49] for example) used in linear inverse problems. In the linear case, the conditional posterior p(u|y, θ) is Gaussian and can be directly sampled, whereas here we sample from p(u|y, θ) using an MH step with RTO proposal.…”

“…Dimension scalability. Algorithm 5.2 is the nonlinear analogue of sampling from the marginal density for linear inverse problems, as found in [6, Section 5.3] and [25,49], where it is shown that sampling from the marginal density removes the dimension scalability issues for the Gibbs sampler described at the end of Section 4. We expect the same result in the nonlinear case when using RTO-PM, since we do not simulate a Markov chain in the high-dimensional model parameter (u) space.…”