Seismic Tomography With a Transdimensional Markov Chain
Thomas Bodin, Malcolm Sambridge
Research School of Earth Sciences, Australian National
University, Canberra, ACT 0200, Australia
Figure 1. Voronoi cells about 30 pseudo random points on the plane. The cell nuclei have been drawn from a 2-D uniform distribution over the spatial domain delimited by the red rectangle. The cell boundaries are defined as the perpendicular bisectors of pairs of nuclei. Any point inside a cell is closer to the nucleus of that cell than any other nucleus.
In seismic tomography, the Earth's interior must be parametrized in some fashion. This is typically done using uniform local cells, and the inversion process consists of finding seismic wave speeds within each cell. The number of cells (and size) has to be chosen according to a compromise between model resolution and model uncertainty. Usually, seismologists choose to have a large number of cells and then face the problem of non-uniqueness by imposing constraints on the solution that are independent of the observations, i.e. by employing regularization procedures like spatial smoothing and norm damping. The type and weighting applied to regularization terms often forms a subjective choice of the user. Another aspect is that the strength of damping and smoothing is determined globally which raises the possibility that, while the ill-constrained regions are being suitably damped, the well constrained regions may be over-smoothed with resulting loss of information from the data.
Our work this year has been devoted to using some new ideas in nonlinear inversion to determine the model dimension (i.e. the number of cells) during the inversion. Treating the number of unkowns as an unknown itself has received little attention in geophysics. However, for more than 10 years, Markov chain Monte Carlo (MCMC) methods that admit transitions between states of differing dimension have been actively developed in the area of Bayesian statistics.
We have developed an approach which uses Voronoi cells instead of a regular mesh for an Earth parametrization (see Figure 1). The Voronoi cells are defined by their centres which are able to move. That is, the number and the position of the cells defining the geometry of the velocity field, as well as the velocity field itself are unknowns in the inversion. The inversion is carried out with a fully non-linear parameter search method based on a trans-dimensional Markov chain.
At each step of the chain, a change from the current model is proposed : we either change the velocity or the position of one random cell. The algorithm also allows jumps between dimensions by adding or removing random cells. The forward problem is computed and provides new estimated travel times. The new misfit to observed travel times is compared to that of the current model. The proposed model is either accepted or rejected using a predefined probabilistic threshold.
Figure 2 : Upper left map shows the true model. The upper right map shows the ray geometry. The lower left map shows the model sampled with the best fit to the data and the lower right map shows the average estimated solution.
The Markov chain produces an ensemble of models with different dimensions which carries relevant statistic information about the unknown velocity field. The method takes as a solution the average of this family of models. Each model in the ensemble has a different parametrization but the average is continuous without obvious 'parametrization' artefacts. The standard deviation of the ensemble forms a continuous map and can be used as a proxy for the error for the solution model.
The method has been tested on synthetic situations where the ray coverage is not uniform and where the parametrisation is an issue (see Figure 2). A major advantage is that explicit regularisation of the model parameters is not required, thus avoiding global damping procedures and the need to find an optimal regularisation value. The technique has also been tested on real data and gives promising results.