A Bayesian linear model for the high-dimensional inverse problem of seismic tomography

We apply a linear Bayesian model to seismic tomography, a high-dimensional inverse problem in geophysics. The objective is to estimate the three-dimensional structure of the earth's interior from data measured at its surface. Since this typically involves estimating thousands of unknowns or more, it has always been treated as a linear(ized) optimization problem. Here we present a Bayesian hierarchical model to estimate the joint distribution of earth structural and earthquake source parameters. An ellipsoidal spatial prior allows to accommodate the layered nature of the earth's mantle. With our efficient algorithm we can sample the posterior distributions for large-scale linear inverse problems and provide precise uncertainty quantification in terms of parameter distributions and credible intervals given the data. We apply the method to a full-fledged tomography problem, an inversion for upper-mantle structure under western North America that involves more than 11,000 parameters. In studies on simulated and real data, we show that our approach retrieves the major structures of the earth's interior as well as classical least-squares minimization, while additionally providing uncertainty assessments.