Geophysical Inverse Problems Research Articles

The inverse dynamic problem of seismic sounding low-frequency regularization

A wave process induced by a plane wave in a vertically inhomogeneous half-space is considered. A traditional geophysical inverse problem is that of recovering the acoustic impedance if the shape of incident wave and the corresponding surface seismogram are known. The classical algorithm for solving this problem based on the reconstruction of the Sturm–Liouville differential operator by using its spectral function is very sensitive to errors in the low-frequency range. We present a regularization method based on approximate reconstruction of the medium with the use of the phase velocities of Love waves.

Invariance, groups, and non-uniqueness: The discrete case

Summary Lie group methods provide a valuable tool for examining invariance and non-uniqueness associated with geophysical inverse problems. The techniques are particularly well suited for the study of non-linear inverse problems. Using the infinitesimal generators of the group it is possible to move within the null space in an iterative fashion. The key computational step in determining the symmetry groups associated with an inverse problem is the singular value decomposition of a sparse matrix. I apply the methodology to the eikonal equation and examine the possible solutions associated with a crosswell tomographic experiment. Results from a synthetic test indicate that it is possible to vary the velocity model significantly and still fit the reference arrival times. The approach is also applied to data from crosswell surveys conducted before and after a CO2 injection at the Lost Hills field in California. The results highlight the fact that a fault cross-cutting the region between the wells may act as a conduit for the flow of water and CO2.

Improved hydrogeophysical characterization using joint inversion of cross‐hole electrical resistance and ground‐penetrating radar traveltime data

Appropriate regularizations of geophysical inverse problems and joint inversion of different data types improve geophysical models and increase their usefulness in hydrogeological studies. We have developed an efficient method to calculate stochastic regularization operators for given geostatistical models. The method, which combines circulant embedding and the diagonalization theorem of circulant matrices, is applicable for stationary geostatistical models when the grid discretization, in each spatial direction, is uniform in the volume of interest. We also used a structural approach to jointly invert cross‐hole electrical resistance and ground‐penetrating radar traveltime data in three dimensions. The two models are coupled by assuming, at all points, that the cross product of the gradients of the two models is zero. No petrophysical relationship between electrical conductivity and relative permittivity is assumed but is instead obtained as a by‐product of the inversion. The approach has been applied to data collected in a U.K. sandstone aquifer in order to improve characterization of the vadose zone hydrostratigraphy. By analyzing scatterplots of electrical conductivity versus relative permittivity together with petrophysical models a zonation could be obtained with corresponding estimates of the electrical formation factor, the water content, and the effective grain radius of the sediments. The approach provides greater insight into the hydrogeological characteristics of the subsurface than by using conventional geophysical inversion methods.

Linear inverse Gaussian theory and geostatistics

Inverse problems in geophysics require the introduction of complex a priori information and are solved using computationally expensive Monte Carlo techniques (where large portions of the model space are explored). The geostatistical method allows for fast integration of complex a priori information in the form of covariance functions and training images. We combine geostatistical methods and inverse problem theory to generate realizations of the posterior probability density function of any Gaussian linear inverse problem, honoring a priori information in the form of a covariance function describing the spatial connectivity of the model space parameters. This is achieved using sequential Gaussian simulation, a well-known, noniterative geostatisticalmethod for generating samples of a Gaussian random field with a given covariance function.This work is a contribution to both linear inverse problem theory and geostatistics. Our main result is an efficient method to generate realizations, actual solutions rather than the conventional least-squares-based approach, to any Gaussian linear inverse problem using a noniterative method. The sequential approach to solving linear and weakly nonlinear problems is computationally efficient compared with traditional least-squares-based inversion. The sequential approach also allows one to solve the inverse problem in only a small part of the model space while conditioned to all available data. From a geostatistical point of view, the method can be used to condition realizations of Gaussian random fields to the possibly noisy linear average observations of the model space.

Trans-dimensional inverse problems, model comparison and the evidence

Summary In most geophysical inverse problems the properties of interest are parametrized using a fixed number of unknowns. In some cases arguments can be used to bound the maximum number of parameters that need to be considered. In others the number of unknowns is set at some arbitrary value and regularization is used to encourage simple, non-extravagant models. In recent times variable or self-adaptive parametrizations have gained in popularity. Rarely, however, is the number of unknowns itself directly treated as an unknown. This situation leads to a transdimensional inverse problem, that is, one where the dimension of the parameter space is a variable to be solved for. This paper discusses trans-dimensional inverse problems from the Bayesian viewpoint. A particular type of Markov chain Monte Carlo (MCMC) sampling algorithm is highlighted which allows probabilistic sampling in variable dimension spaces. A quantity termed the evidence or marginal likelihood plays a key role in this type of problem. It is shown that once evidence calculations are performed, the results of complex variable dimension sampling algorithms can be replicated with simple and more familiar fixed dimensional MCMC sampling techniques. Numerical examples are used to illustrate the main points. The evidence can be difficult to calculate, especially in high-dimensional non-linear inverse problems. Nevertheless some general strategies are discussed and analytical expressions given for certain linear problems.

Two ways to quantify uncertainty in geophysical inverse problems

We present two approaches to invert geophysical measurements and estimate subsurface properties and their uncertainties when little is known a priori about the size of the errors associated with the data. We illustrate these approaches by inverting first-arrival traveltimes of seismic waves measured in a vertical well to infer the variation of compressional slowness in depth. First, we describe a Bayesian formulation based on probability distributions that define prior knowledge about the slowness and the data errors. We use an empirical Bayes approach, where hyperparameters are not well known ahead of time (e.g., the variance of the data errors) and are estimated from their most probable value, given the data. The second approach is a non-Bayesian formulation that we call spectral, in the sense that it uses the power spectral density of the traveltime data to constrain the inversion (e.g., to estimate the variance of the data errors). In the spectral approach, we vary assumptions made about the characteristics of the slowness signal and evaluate the resulting slowness estimates and their uncertainties. This approach is computationally simple and starts from a few assumptions. These assumptions can be checked during the analysis. On the other hand, it requires evenly spaced traveltime measurements, and it cannot be extended easily (e.g., to data that have gaps). In contrast, the Bayesian framework is based on a general theory that can be generalized immediately, but it is more involved computationally. Despite the conceptual and practical differences, we find that the two approaches give the same results when they start from the same assumptions: The allegiance to a Bayesian or non-Bayesian formulation matters less than what one is willing to assume when solving the inverse problem.

Solutions, algorithms and inter-relations for local minimization search geophysical inversion

This paper presents in a general form the most popular local minimization search solutions for geophysical inverse problems—the Tikhonov regularization solutions, the smoothest model solutions and the subspace solutions, from which the inter-relationships between these solutions are revealed. For the Tikhonov regularization solution, a variety of forms exist—the general iterative formula, the iterative linearized scheme, the Levenberg–Marquardt version, the conjugate gradient solver (CGS) and local-search quadratic approximation CGS. It is shown here that the first three solutions are equivalent and are just a specified form of the gradient solution. The local-search quadratic approximation CGS is shown to be a more general form from which these three solutions emerge. The smoothest model solution (Occam's inversion) is a subset of the Tikhonov regularization solutions, for which a practical algorithm generally comprises two parts: (1) determine the subset of the Tikhonov regularization solutions which satisfy the desired tolerance of data fit, (2) choose the solution which best fits the data in the subset. The solution procedure is actually an adaptive Tikhonov regularization solution. The subspace solution is mathematically no more than a dimension-decreased transform of model parameterization in the Tikhonov regularization solution or the smoothest model solution. The crucial step is to choose the basis vectors which span the whole model space in the parameterization transform. A simple form of the transform is developed which yields greater flexibility for geophysical inverse applications.

Characterising the parameter space of a highly nonlinear inverse problem

In inverse problems, often there is no available analytical expression relating the physical quantities of interest and the available data. In these cases, one resorts to using a numerical model with a finite number of parameters, resulting in a discrete problem. Also, many discrete inverse problems involve a highly nonlinear mapping between the model parameters and the simulation of the data by the model. Algorithms exist for estimating the model parameters in nonlinear discrete inverse problems. However, one needs to investigate how these estimated models relate to the true structure of the studied system (i.e., the truth model). This is known as model appraisal and it is greatly affected by three sources of uncertainty: misleading search, non-uniqueness and errors. In this work, we aim at analysing the impact of the first two sources of uncertainty in model appraisal (misleading search and non-uniqueness) by characterising the parameter space of a highly nonlinear geophysical inverse problem. Our approach is to use an error-free synthetic data problem so that the truth model is known in advance and therefore the results can be compared against a known truth. The characterisation requires an exceptionally high quality sampling of the parameter space. A new Real-parameter Genetic Algorithm called Genetic Sampler (GS) is used for this purpose. The GS successfully characterises the parameter space by locating multiple unconnected optimal regions, including that of the truth model, in several instances of the inverse problem. It is argued that, given the characteristics of the parameter space, many search methods could misidentify a local optimum as the global optimum. Such a misleading search will, in this case, have terrible consequences, as all the local optima had no predictive value. Another important result is related to the issue of non-uniqueness. It seems that, for this model and in the absence of any error, there is a single global optimum, which constitutes a unique solution to this nonlinear inverse problem. Lastly, we investigated the effect of including more data in the inversion process. This resulted in the parameter space having fewer and less prominent locally optimal regions, which makes the search for the global optimum easier.

Quantum Annealing Inversion and Its Implementation

AbstractQuantum annealing (QA) is a new method to solve nonlinear inverse problems based on quantum tunneling effect of state transitions. Compared with simulated annealing, quantum annealing has a potential advantage in some aspects such as annealing convergent rate and avoiding being trapped in a local minimum. So the study of quantum annealing inversion is of important academic significance and practical value. In this paper, the inversion tests of some models and actual data show that QA has higher efficiency and accuracy. It is the first time to apply QA method to the geophysical inverse problem.

Joint analysis of refractions with surface waves: An inverse solution to the refraction-traveltime problem

We describe a possible solution to the inverse refraction-traveltime problem (IRTP) that reduces the range of possible solutions (nonuniqueness). This approach uses a reference model, derived from surface-wave shear-wave velocity estimates, as a constraint. The application of the joint analysis of refractions with surface waves (JARS) method provided a more realistic solution than the conventional refraction/tomography methods, which did not benefit from a reference model derived from real data. This confirmed our conclusion that the proposed method is an advancement in the IRTP analysis. The unique basic principles of the JARS method might be applicable to other inverse geophysical problems.

Nonlinear inversion of potential-field data using a hybrid-encoding genetic algorithm

Using a genetic algorithm to solve an inverse problem of complex nonlinear geophysical equations is advantageous because it does not require computer gradients of models or “good” initial models. The multi-point search of a genetic algorithm makes it easier to find the globally optimal solution while avoiding falling into a local extremum. As is the case in other optimization approaches, the search efficiency for a genetic algorithm is vital in finding desired solutions successfully in a multi-dimensional model space. A binary-encoding genetic algorithm is hardly ever used to resolve an optimization problem such as a simple geophysical inversion with only three unknowns. The encoding mechanism, genetic operators, and population size of the genetic algorithm greatly affect search processes in the evolution. It is clear that improved operators and proper population size promote the convergence. Nevertheless, not all genetic operations perform perfectly while searching under either a uniform binary or a decimal encoding system. With the binary encoding mechanism, the crossover scheme may produce more new individuals than with the decimal encoding. On the other hand, the mutation scheme in a decimal encoding system will create new genes larger in scope than those in the binary encoding. This paper discusses approaches of exploiting the search potential of genetic operations in the two encoding systems and presents an approach with a hybrid-encoding mechanism, multi-point crossover, and dynamic population size for geophysical inversion. We present a method that is based on the routine in which the mutation operation is conducted in the decimal code and multi-point crossover operation in the binary code. The mix-encoding algorithm is called the hybrid-encoding genetic algorithm (HEGA). HEGA provides better genes with a higher probability by a mutation operator and improves genetic algorithms in resolving complicated geophysical inverse problems. Another significant result is that final solution is determined by the average model derived from multiple trials instead of one computation due to the randomness in a genetic algorithm procedure. These advantages were demonstrated by synthetic and real-world examples of inversion of potential-field data.

Combined global/linear inversion of well-logging data in layer-wise homogeneous and inhomogeneous media

In the paper a combined inversion algorithm solving the nonlinear geophysical well-logging inverse problem is presented. We apply a successive combination of a float-encoded Genetic Algorithm as a global optimization method and the well-known linearized Marquardt algorithm forming a fast inversion procedure. The technique is able to decrease the CPU run time at least one order of magnitude compared to the Genetic Algorithm and gives the parameter estimation errors having a few linearized optimization steps at the end of the iteration process. We use depth-dependent tool response equations to invert all the data of a greater depth-interval jointly in order to determine petrophysical parameters of homogeneous or inhomogeneous layers in one inversion procedure. The so-called interval inversion method gives more accurate and reliable estimation for the petrophysical model parameters than the conventional point by point inversion methods. It also enables us to determine the layer-thicknesses that can not be extracted from the data set by means of conventional inversion techniques. We test the combined interval inversion method on synthetic data, and employ it to the interpretation of well logs measured in a Hungarian hydrocarbon exploratory borehole.

Bayesian inference in satellite gravity inversion

To solve a geophysical inverse problem means applying measurements to determine the parameters of the selected model. The inverse problem is formulated as the Bayesian inference. The Gaussian probability density functions are applied in the Bayes's equation. The CHAMP satellite gravity data are determined at the altitude of 400 km altitude over the South part of the Pannonian Basin. The model of interpretation is the right vertical cylinder. The parameters of the model are obtained from the minimum problem solved by the Simplex method.

Using a genetic algorithm for 3-D inversion of gravity data in Fuerteventura (Canary Islands)

The use of genetic algorithms in geophysical inverse problems is a relatively recent development and offers many advantages in dealing with the non-linearity inherent in such applications. We have implemented a genetic algorithm to efficiently invert a set of gravity data. Employing several fixed density contrasts, this algorithm determines the geometry of the sources of the anomaly gravity field in a 3-D context. The genetic algorithms, based on Darwin’s theory of evolution, seek the optimum solution from an initial population of models, working with a set of parameters by means of modifications in successive iterations or generations. This searching method traditionally consists of three operators (selection, crossover and mutation) acting on each generation, but we have added a further one, which smoothes the obtained models. In this way, we have designed an efficient inversion gravity method, confirmed by both a synthetic example and a real data set from the island of Fuerteventura. In the latter case, we identify crustal structures related to the origin and evolution of the island. The results show a clear correlation between the sources of gravity field in the model and the three volcanic complexes recognized in Fuerteventura by other geological studies.

The Inverse Problem of Refraction Travel Times, Part II: Quantifying Refraction Nonuniqueness Using a Three-layer Model

This paper is the second of a set of two papers in which we study the inverse refraction problem. The first paper, “Types of Geophysical Nonuniqueness through Minimization,” studies and classifies the types of nonuniqueness that exist when solving inverse problems depending on the participation of a priori information required to obtain reliable solutions of inverse geophysical problems. In view of the classification developed, in this paper we study the type of nonuniqueness associated with the inverse refraction problem. An approach for obtaining a realistic solution to the inverse refraction problem is offered in a third paper that is in preparation. The nonuniqueness of the inverse refraction problem is examined by using a simple three-layer model. Like many other inverse geophysical problems, the inverse refraction problem does not have a unique solution. Conventionally, nonuniqueness is considered to be a result of insufficient data and/or error in the data, for any fixed number of model parameters. This study illustrates that even for overdetermined and error free data, nonlinear inverse refraction problems exhibit exact-data nonuniqueness, which further complicates the problem of nonuniqueness. By evaluating the nonuniqueness of the inverse refraction problem, this paper targets the improvement of refraction inversion algorithms, and as a result, the achievement of more realistic solutions. The nonuniqueness of the inverse refraction problem is examined initially by using a simple three-layer model. The observations and conclusions of the three-layer model nonuniqueness study are used to evaluate the nonuniqueness of more complicated n-layer models and multi-parameter cell models such as in refraction tomography. For any fixed number of model parameters, the inverse refraction problem exhibits continuous ranges of exact-data nonuniqueness. Such an unfavorable type of nonuniqueness can be uniquely solved only by providing abundant a priori information. Insufficient a priori information during the inversion is the reason why refraction methods often may not produce desired results or even fail. This work also demonstrates that the application of the smoothing constraints, typical when solving ill-posed inverse problems, has a dual and contradictory role when applied to the ill-posed inverse problem of refraction travel times. This observation indicates that smoothing constraints may play such a two-fold role when applied to other inverse problems. Other factors that contribute to inverse-refraction-problem nonuniqueness are also considered, including indeterminacy, statistical data-error distribution, numerical error and instability, finite data, and model parameters.

The Inverse Problem of Refraction Travel Times, Part I:Types of Geophysical Nonuniqueness Through Minimization

In a set of two papers we study the inverse problem of refraction travel times. The purpose of this work is to use the study as a basis for development of more sophisticated methods for finding more reliable solutions to the inverse problem of refraction travel times, which is known to be nonunique. The first paper, “Types of Geophysical Nonuniqueness through Minimization,” emphasizes the existence of different forms of nonuniqueness in the realm of inverse geophysical problems. Each type of nonuniqueness requires a different type and amount of a priori information to acquire a reliable solution. Based on such coupling, a nonuniqueness classification is designed. Therefore, since most inverse geophysical problems are nonunique, each inverse problem must be studied to define what type of nonuniqueness it belongs to and thus determine what type of a priori information is necessary to find a realistic solution. The second paper, “Quantifying Refraction Nonuniqueness Using a Three-layer Model,” serves as an example of such an approach. However, its main purpose is to provide a better understanding of the inverse refraction problem by studying the type of nonuniqueness it possesses. An approach for obtaining a realistic solution to the inverse refraction problem is planned to be offered in a third paper that is in preparation. The main goal of this paper is to redefine the existing generalized notion of nonuniqueness and a priori information by offering a classified, discriminate structure. Nonuniqueness is often encountered when trying to solve inverse problems. However, possible nonuniqueness diversity is typically neglected and nonuniqueness is regarded as a whole, as an unpleasant “black box” and is approached in the same manner by applying smoothing constraints, damping constraints with respect to the solution increment and, rarely, damping constraints with respect to some sparse reference information about the true parameters. In practice, when solving geophysical problems different types of nonuniqueness exist, and thus there are different ways to solve the problems. Nonuniqueness is usually regarded as due to data error, assuming the true geology is acceptably approximated by simple mathematical models. Compounding the nonlinear problems, geophysical applications routinely exhibit exact-data nonuniqueness even for models with very few parameters adding to the nonuniqueness due to data error. While nonuniqueness variations have been defined earlier, they have not been linked to specific use of a priori information necessary to resolve each case. Four types of nonuniqueness, typical for minimization problems are defined with the corresponding methods for inclusion of a priori information to find a realistic solution without resorting to a non-discriminative approach. The above-developed stand-alone classification is expected to be helpful when solving any geophysical inverse problems.

A Homotopy Method for Solving the Impedance Inversion of 1‐D Wave Equation

AbstractBased on the inversion of one‐dimensional seismic wave equation, we study a homotopy method for determining the parameters of the earth. With this method, the problem of nonlinear system of equations is reduced to an initial value problem of ordinary differential equations. This algorithm is not only fast and stable but also globally convergent. Both synthetic data and field data confirm the effectiveness of the homotopy method as an inversion algorithm. It is especially appropriate for non‐linear geophysical inversion problems with multi‐extreme, and has widespread application prospect in non‐linear geophysical inversion.

Adjustment of regularization in ill‐posed linear inverse problems by the empirical Bayes approach

ABSTRACTRegularization is the most popular technique to overcome the null space of model parameters in geophysical inverse problems, and is implemented by including a constraint term as well as the data‐misfit term in the objective function being minimized. The weighting of the constraint term relative to the data‐fitting term is controlled by a regularization parameter, and its adjustment to obtain the best model has received much attention. The empirical Bayes approach discussed in this paper determines the optimum value of the regularization parameter from a given data set. The regularization term can be regarded as representing a priori information about the model parameters. The empirical Bayes approach and its more practical variant, Akaike's Bayesian Information Criterion, adjust the regularization parameter automatically in response to the level of data noise and to the suitability of the assumed a priori model information for the given data. When the noise level is high, the regularization parameter is made large, which means that the a priori information is emphasized. If the assumed a priori information is not suitable for the given data, the regularization parameter is made small. Both these behaviours are desirable characteristics for the regularized solutions of practical inverse problems. Four simple examples are presented to illustrate these characteristics for an underdetermined problem, a problem adopting an improper prior constraint and a problem having an unknown data variance, all frequently encountered geophysical inverse problems. Numerical experiments using Akaike's Bayesian Information Criterion for synthetic data provide results consistent with these characteristics. In addition, concerning the selection of an appropriate type of a priori model information, a comparison between four types of difference‐operator model – the zeroth‐, first‐, second‐ and third‐order difference‐operator models – suggests that the automatic determination of the optimum regularization parameter becomes more difficult with increasing order of the difference operators. Accordingly, taking the effect of data noise into account, it is better to employ the lower‐order difference‐operator models for inversions of noisy data.

Enhancing the resolving power of least‐squares inversion with active constraint balancing

Most geophysical inverse problems are solved using least‐squares inversion schemes with damping or smoothness constraints to improve stability and convergence rate. Since the Lagrangian multiplier controls resolution and stability of the inverse problem, we always want to use the optimum multiplier, which is not easy to get and is usually obtained by experience or a time‐consuming optimization process. We present a new regularization approach, in which the Lagrangian multiplier is set as a spatial variable at each parameterized block and automatically determined via the parameter resolution matrix and spread function analysis. For highly resolvable parameters, a small value of the Lagrangian multiplier is assigned, and vice versa. This approach, named “active constraint balancing” (ACB), tries to balance the constraints of the least‐squares inversion according to sensitivity for a given problem so that it enhances the resolution as well as the stability of the inversion process. We demonstrate the performance of the ACB by applying it to a two‐dimensional resistivity tomography problem, which results in a remarkable enhancement of the spatial resolution. Enhancement of the resolution is also verified in the application of resistivity tomography to a field data set acquired at a tunnel construction site.

Multiresolution parameterization for geophysical inverse problems

A model parameterization based on a multiresolution wavelet representation is proposed for geophysical inverse problems in this study. The matrix equation for the wavelet representation of the model is constructed by dual wavelet transforms of a readily built matrix based on the pixel parameterization. We demonstrate the merits of the proposed parameterization using the simple examples of geophysical data gridding and 3D seismic tomography, and comparing the results with those obtained by the conventional pixel parameterization. It is shown that for the conventional damped least‐squares solutions of pixel‐parameterized models, regularization by norm damping deteriorates the underlying model correlation due to sparse and heterogeneous sampling, whereas roughness damping schemes impose a priori arbitrary correlation length irrespective of the heterogeneous sampling. The scale hierarchy embedded in the proposed multiresolution parameterization facilitates the compromise between the dual spatial and scale resolution. Depending on the local richness of the data constraints across different scales at a site, model variation of various scale spectra can be robustly resolved through the scale hierarchy. There is thus no need to invoke additional smoothness regularization.