Comparison between Tikhonov regularization and truncated SVD in gravity data inversion



In this paper the 3D inversion of gravity data using two different regularization methods, namely Tikhonov regularization and truncated singular value decomposition (TSVD), is considered. The earth under the survey area is modeled using a large number of rectangular prisms, in which the size of the prisms are kept fixed during the inversion and the values of densities of the prisms are the model parameters to be determined. A depth weighting matrix is used to counteract the natural decay of the kernel, so the inversion obtains reliable information about the source distribution with respect to depth. To generate a sharp and focused model, the minimum support (MS) constraint is used, which minimizes the total area with non zero departure of the model parameters from a given a priori model. Then, the application of iteratively reweighted least square algorithm is required to deal with non-linearity introduced by MS constraint. At each iteration of the inversion, a priori variable weighting matrix is updated using model parameters obtained at the previous iteration. We use the singular value decomposition (SVD) for computing Tikhonov solution, which also helps us to compare the results with the solution obtained by TSVD. Thus, the algorithms presented here are suitable for small to moderate size problems, where it is feasible to compute the SVD. In Tikhonov regularization method, the optimal regularization parameter at each iteration is obtained by application of the parameter-choice method. The method is based on the statistical distribution of the minimum of the Tikhonov function. For weighting of the data fidelity by a known Gaussian noise distribution on the measured data and, when the regularization term is considered to be weighted by unknown inverse covariance information on the model parameters, the minimum of the Tikhonov functional becomes a random variable that follows a  distribution. Then, a Newton root-finding algorithm can be used to find the regularization parameter. For truncated SVD regularization, the Picard plot is used to find a suitable value of truncation index. In math literature, a plot of singular values together with SVD and solution coefficients is often referred to as Picard plot. To test the algorithms, a density model which consists of a dipping dike embedded in a uniform half-space is used. The surface gravity anomaly produced by this model is contaminated with three different noise levels, and are used as input for introduced inversion algorithms. The results indicate that the algorithms are able to recover the geometry and density distribution of the original model. In general, the reconstructed model is more sparse using TSVD method as compare with Tikhonov solution. This especially happens for high noise level, where there is an important difference between two solutions. In this case, while TSVD produces a sparse model, the solution of Tikhonov regularization is not sparse. Furthermore, the number of iterations, which is required to terminate the algorithms, is more for TSVD as compare with Tikhonov method. This feature, along with automatic determination of regularization parameter, makes the implementation of the Tikhonov regularization method faster than TSVD. The inversion methods are used on real gravity data acquired over the Gotvand dam site in the south-west of Iran. Tertiary deposits of the Gachsaran formation are the dominant geological structure in this area, and it is mainly comprised of marl, gypsum, anhydrite and halite. There are several solution cavities in the area so that relative negative anomalies are distinguishable in the residual map. A window of residual map consists of 640 gridded data, which includes three negative anomalies, that is selected for modeling. The reconstructed models are shown and compare with results obtained by bore holes.


Main Subjects

Ardestani, V. E., 2013, Detecting, delineating and modeling the connected solution cavities in a dam site via microgravity data, Acta Geodaetica and Geophysica, 48, 123-138.
Blakely, R. J., 1996, Potential theory in gravity and magnetic applications, Cambridge University Press, Cambridge.
Boulanger, O. and Chouteau M., 2001, Constraint in 3D gravity inversion, Geophysical Prospecting, 49, 265-280.
Farquharson, C. G. and Oldenburg, D. W., 2004, A comparison of automatic techniques for estimating the regularization parameter in non-linear inverse problems, Geophys. J. Int., 156, 411-425.
Hansen, P. C., 1998, Rank-deficient and discrete Ill-posed problems: numerical aspects of linear inversion, Monographs on Mathematical Modeling and Computation, 4, SIAM.
Last, B. J. and Kubik, K., 1983, Compact gravity inversion, Geophysics, 48, 713-721.
Li, Y. and Oldenburg, D. W., 1996, 3D inversion of magnetic data, Geophysics, 61, 394-408.
Pilkington, M., 1997, 3D magnetic imagingusing conjugate gradients, Geophysics, 62, 1132-1142.
Portniaguine, O. and Zhdanov, M. S. 1999, Focusing geophysical inversion images, Geophysics, 64, 874-887.
Vatankhah, S., Renaut, R. A. and Ardestani, V. E., 2014, Regularization parameter estimation for underdetermined problems by the principle with application to 2D focusing gravity inversion, Inverse Problems, 30, 085002
Vatankhah, S., Ardestani, V. E., and Renaut R. A., 2015, Application of the  principle and unbiased predictive risk estimator for determining the regularization parameter in 3-D focusing gravity inversion, Geophys. J. Int., 200, 265-277