Varying dimensional Bayesian acoustic waveform inversion for 1D semi-infinite heterogeneous media

Highlights

• Innovative methodology on probabilistic model selection for varying dimensional wave models.

• Bayesian inversion method to quantify uncertainty in geologic media for inversion process.

• Use of seismic forward model accommodating realistic far-field boundary conditions via PMLs.

Abstract

This paper introduces a methodology to infer the spatial variation of the acoustic characteristics of a 1D vertical elastic heterogeneous earth model via a Bayesian calibration approach, given a prescribed sequence of loading and the corresponding time history response registered at the ground level. This involves solving an inverse problem that maps the ground seismic response onto a random profile of the ground stratigraphy (i.e. a 1D continuous spatial random field). From a Bayesian point of view, the solution to an inverse problem is fully characterized by a posterior density function of the forward model random parameters, which explicitly overcomes the solution's non-uniqueness. This subsurface earth model is parameterized using a Bayesian partition model, where the number of soil layers, the location of the layers' interfaces, and their corresponding mechanical characteristics are defined as random variables. The partition model approach to an inverse problem is closely related to a Bayesian model selection problem, where the likely dimensionality of the inverse problem (number of unknowns) is inferred conditioned on the experimental observations. The main benefit of the proposed approach is that the explicit regularization of the inverted profile by global damping procedures is not required. A Reversible Jump Markov Chain Monte Carlo (RJMCMC) algorithm is used to sample the target posterior of varying dimension, dependent on the number of layers. A synthetic case study is provided to indicate the applicability of the proposed methodology.

Introduction

A subsurface earth model is composed of complex geophysical formations, which embodies a wide range of physical and mechanical heterogeneities. The aim of probabilistic inverse modeling is to reconstruct the random field structure of these subsurface properties, while accounting for various sources of uncertainty stemming from ground seismic geologic measurement errors, aleatory formations, and limited theoretical understanding about underground wave propagation.

In practice, one of the main goals of geophysical investigations is to identify the main geomorphological features of an unknown medium, meaning the spatial location and concentration of geological features such as the transition between materials, discontinuities and material concentrations [34]. In the case of a vertical 1D profile, this requires the definition of the location of the sharp transitions between layer interfaces (i.e. material properties), and the characterization of its corresponding mechanical properties.

In a horizontally stratified earth model, prior to making an inference about the likely variation of the elastic parameters within the geological layers, an assumption must be made concerning the number of layers in a certain depth range of interest. This assumption defines the dimensionality (i.e., the number of unknowns) of the inverse problem. In reality, however, such information is rarely available for the dimension and definition of the parameter space to be fixed. Consequently, fixing the number of layers based on an incorrect assumption results in an erroneous subsurface characterization.

To relax the hypothesis about the subsurface structure or spatial layering of the media's mechanical parameters, before the forward model is calibrated, it is proposed to define the number of layers, their locations, and their corresponding mechanical parameters, as random variables. From a Bayesian perspective, this set up is closely associated with probabilistic model selection, where a collection of models with varying number of parameters are presented for inversion, and the task is to select the models that most likely describe the experimental observations.

To illustrate the applicability of the proposed probabilistic calibration method, a one dimensional horizontally stratified medium is presented in terms of a Bayesian partition model [8]. The partition model divides the unknown material field into a number of non-overlapping regions, where each region represents a soil layer. Formulating the inverse medium problem in terms of a partition model may help reduce the dimensionality of the parameter space. Hence, regularizing the solution through specific prior distributions, which bears smoothness constraints (in a Bayesian inversion framework [45], [12], [20]), or regularization terms (in deterministic optimization problems [44], [36], [15]), is precluded.

A generalization of the simulation-based Markov Chain Monte Carlo methods, the so-called reversible jump [18], is used to sample the posterior distribution of varying dimensionality. In this setting, the Markov Chain is capable of undergoing dimension changes while moving among a number of candidate models. The key aspect of the reversible jump algorithm is the introduction of some auxiliary random variables to equalize the dimensionality of the parameter space across models. A series of one-to-one deterministic functions are defined to perform dimension matching such that the detailed balance condition is satisfied. Balance condition is the necessary condition for a Markov Chain to converge to the target density [7].

Since the introduction of Bayesian inference methods to the geophysical community, this has received a great deal of attention in a variety of applications [14], [16], [17], [45], [41], [43]. However, a limited number of studies have addressed the subsurface parameter estimation as a model selection problem, many of which resort to approximate methods to fulfill the model determination [10], [13]. The varying dimensional formulation was first introduced to the geophysics literature by Malinverno [31] in a 1D-DC resistivity sounding inversion, and later implemented in a number of geophysical probing inverse problems [39], [11], [1], [35].

The major impact of utilizing probabilistic calibration via a Bayesian approach is the systematic exploration of all combinations of the model parameters through a transparent definition of the impact of the participating uncertainty sources. During such exhaustive parameter exploration, a probability metric is defined to assess the likelihood of selecting sets of parameters that serve to approximate the experimental observations (likelihood); but also a probability density is defined to reflect the degree of knowledge on the model parameters (prior) before the model inversion. The combination of these two states of knowledge about the model of interest yields the following benefits: a transition from deterministic to probabilistic model parameters, assessment of the type and degree of correlation between the model parameters (e.g. linear or non-linear), measurement of the impact of the varying experimental observations (e.g. the effect of the number of observations on the prediction of confidence levels), assessment of the model performance, and most importantly, that among a number of competing models to choose from, selecting the best model which can describe the process that generated the observations. Key applications of the probabilistic subsurface imaging include integrated site investigation, since the recovery of geophysical mechanical parameters allows enhanced geomechanical characterization. [33], [34]. The varying parameter dimensionality is formulated through a Bayesian inversion, to populate likely configurations of a heterogeneous elastic medium occupying a semi-infinite domain.

A key defining characteristic of full waveform inversion is the numerical solution of the equations of motion. The governing forward physics consist of a 1D transient scalar acoustic wave propagation, where in order to model the semi-infinite extent of the physical domain, a perfectly matched layer (PML) is introduced at the truncation boundary to emulate the infiniteness of the earth structure [22]. A displacement–stress mixed finite element scheme is used for the numerical solution of the PML-augmented wave PDE.