Uncertainty quantification for inverse problems with weak partial-differential-equation constraints

Abstract

In statistical inverse problems, the objective is a complete statistical description of unknown parameters from noisy observations to quantify uncertainties in unknown parameters. We consider inverse problems with partial-differential-equation (PDE) constraints, which are applicable to many seismic problems. Bayesian inference is one of the most widely used approaches to precisely quantify statistics through a posterior distribution, incorporating uncertainties in observed data, modeling kernel, and prior knowledge of parameters. Typically when formulating the posterior distribution, the PDE constraints are required to be exactly satisfied, resulting in a highly nonlinear forward map and a posterior distribution with many local maxima. These drawbacks make it difficult to find an appropriate approximation for the posterior distribution. Another complicating factor is that traditional Markov chain Monte Carlo (MCMC) methods are known to converge slowly for realistically sized problems. To overcome these drawbacks, we relax the PDE constraints by introducing an auxiliary variable, which allows for Gaussian errors in the PDE and yields a bilinear posterior distribution with weak PDE constraints that is more amenable to uncertainty quantification because of its special structure. We determine that for a particular range of variance choices for the PDE misfit term, the new posterior distribution has fewer modes and can be well-approximated by a Gaussian distribution, which can then be sampled in a straightforward manner. Because it is prohibitively expensive to explicitly construct the dense covariance matrix of the Gaussian approximation for problems with more than [Formula: see text] unknowns, we have developed a method to implicitly construct it, which enables efficient sampling. We apply this framework to 2D seismic inverse problems with 1800 and 92,455 unknown parameters. The results illustrate that our framework can produce comparable statistical quantities with those produced by conventional MCMC-type methods while requiring far fewer PDE solves, which are the main computational bottlenecks in these problems.

Ensemble-based uncertainty estimation in Full Waveform Inversion

“…Though this cost is indeed reasonably low, it does not include the computational cost of the reverse time migration they are using to precondition their sampling. Finally, the number of PDE solved to sample the posterior covariance in Fang et al (2018) proposition, is the number of sources plus the number of receivers per frequencies (not including the number of PDE to solve the inverse problem). Besides, this method does seem to display challenging memories limitation as it requires to store the optimal wavefields in memory for each frequency bands, which may become challenging for large scale 3D application.…”

“…Zhu et al (2016), Eliasson & Romdhane (2017) and Liu & Peter (2019) are relying on the randomized Singular Value Decomposition (SVD) to estimate a truncated Hessian in a tractable way. Finally, using the Wavefield Reconstruction Inversion (WRI) to relax the inverse problem formulation, Fang et al (2018) demonstrate that the WRI cost function is particularly suited for the quadratic approximation that is assumed in all the methods mentioned above and therefore justifies the assumption of Gaussianity of the posterior distribution. To estimate uncertainty, they approximate the Gauss-Newton Hessian, from which the square root makes it possible to sample the posterior covariance once the approximate Gauss-Newton Hessian is computed and stored.…”

Seismic wavefield imaging of Earth’s interior across scales

“…Ideally, the currently used quasi-Newton methods would be abandoned in favour of Bayesian inference 105,233 , which provides the complete, possibly multi-modal, a posteriori model distribution, but for this to happen, the computational cost of FWI would have to be reduced dramatically. Meanwhile, there have been successful attempts at sampling the a posteriori model distribution in the vicinity of the global minimum based on random probing 234,235 , Kalman filtering 236,237 and the square-root variable-metric method 238,239 . Synthetic examples of these techniques have focused on the Marmousi model 240 , as well as a real data set from the Valhall oil field 236 .…”

Bayes Meets Tikhonov: Understanding Uncertainty Within Gaussian Framework for Seismic Inversion

“…It remains challenging to solve inverse problems practically due to limitations in data acquisition, measurement uncertainties and the non-uniqueness of the solution. Several advancements in geophysics have been proposed to tackle these challenges [8,9,10,11,12]. Despite these rapid advancements, the research on uncertainty analysis for the seismic inversion solutions is progressing slowly due to limitations of computational power and algorithms advancement.…”