Log In Sign Up

Stochastic uncertainty analysis of gravity gradient tensor components and their combinations

by   Pejman Shamsipour, et al.

Full tensor gravity (FTG) devices provide up to five independent components of the gravity gradient tensor. However, we do not yet have a quantitative understanding of which tensor components or combinations of components are more important to recover a subsurface density model by gravity inversion. This is mainly because different components may be more appropriate in different scenarios or purposes. Knowledge of these components in different environments can aid with selection of optimal selection of component combinations. In this work, we propose to apply stochastic inversion to assess the uncertainty of gravity gradient tensor components and their combinations. The method is therefore a quantitative approach. The applied method here is based on the geostatistical inversion (Gaussian process regression) concept using cokriging. The cokriging variances (variance function of the GP) are found to be a useful indicator for distinguishing the gravity gradient tensor components. This approach is applied to the New Found dataset to demonstrate its effectiveness in real-world applications.


page 10

page 11

page 13

page 20

page 21

page 22


Understanding Deflation Process in Over-parametrized Tensor Decomposition

In this paper we study the training dynamics for gradient flow on over-p...

Tensor Basis Gaussian Process Models of Hyperelastic Materials

In this work, we develop Gaussian process regression (GPR) models of hyp...

Inductive Framework for Multi-Aspect Streaming Tensor Completion with Side Information

Low-rank tensor completion is a well-studied problem and has application...

Autoregressive Density Modeling with the Gaussian Process Mixture Transition Distribution

We develop a mixture model for transition density approximation, togethe...

Additive Tensor Decomposition Considering Structural Data Information

Tensor data with rich structural information becomes increasingly import...

Efficient AutoML Pipeline Search with Matrix and Tensor Factorization

Data scientists seeking a good supervised learning model on a new datase...

Lower Bounds for the Convergence of Tensor Power Iteration on Random Overcomplete Models

Tensor decomposition serves as a powerful primitive in statistics and ma...

1 Introduction

Gravity gradiometry permits single components and combinations of components to be applied in the interpretation and inversion of data. The effectiveness of choosing either the single components or the combinations highly depends on the information content of the components and their combinations. a structure’s trends, noise, orientation, among others, can significantly effect which components are to be selected. As there are a wide range of combinations and more are likely to be suggested in the future, a quantitative ranking of these values is vital.

The majority of quantitative analyses of tensor components were limited to the comparison of gravity measurements and gradiometry gradient tensor (Murphy and Brewster, 2007; Jordan, 1978; Vasco, 1989).

In other approaches, authors attempted to compare the 3D inversion results with known geology qualitatively (Martinez et al. (2012); Zhdanov et al. (2004); Li (2001)). Previous works (Pilkington, 2012, 2014; Mikhailov et al., 2007; Beiki and Pedersen, 2010) compared tensor components and their combinations in terms of their information content and error levels when inverting for simple models. In these cases the models used were either very simple or did not include the effects of full covariance errors. Here, we use cokriging to study different tensor components using more complex and realistic models.

Commonly, cokriging is used as a mathematical interpolation and extrapolation tool that uses the spatial correlation between primary and secondary variables to improve the estimation of the primary variable. However, in this study cokriging can be seen as an inversion method

(Shamsipour et al., 2010; Tchikaya et al., 2016). The benefit of performing inversion using cokriging is that we can use the cokriging error as an indicator to determine the uncertainty of estimated densities obtained from different gravity gradients or any combination of them. Even though the cokriging error does not directly depend on the data values used for estimation, it can be used as an indicator of uncertainty.

The cokriging variance is available for each prism, so we can determine which component gives a better density estimation at each individual location. After testing, we realized that the cokriging variance for different components or their combinations is distinguishable. In order to investigate the effectiveness of different components, we first assume a geometry with fixed covariance parameters. The only difference for the different components is the kernel, or Green’s function. However it is possible to have different geometries by changing the parameters of the covariance (lengthscales/ranges and rotation angles) and investigate which component gives the better estimate for each geometry. We can also consider the effect of the noise in each case by adding a nugget effect to the covariance model. To validate our experiments, we conduct stochastic inversion on gravity gradiometer data provided in the New Found dataset.


1.1 Stochastic inversion and error of estimated model

Shamsipour et al. (2010) presents the application of cokriging to invert potential field data. The method can be applied to any linear geophysical problem. Consider that there are data observations and parameters (physical properties) on rectangular prisms (), their relationship can be written in matrix form:


where is the matrix of the geometric terms (kernel). Since and are linearly related, then, their covariance matrices are also linearly related:


where and are respectively the and covariance matrices for observations and parameters, and is the

data error covariance matrix which is usually modeled as a nugget effect (white noise). We also have:


where is the cross-covariance matrix between the observation data and parameters. Finally, the estimates of parameters (primary variable) are obtained from the observation data (secondary variable):


where is the optimal matrix of weighting coefficients obtained by minimization of the estimation variance. The estimation variances are found on the diagonal of the following matrix (Myers, 1982):


where is the matrix of weighting coefficients. Here and

are multidimensional random variables. Minimization of the above estimation variance yields the simple cokriging solution


Finally, the vector of cokriging variances is obtained from:


where diag(.) finds the diagonal elements of a matrix and put them in a vector. Note that the maximum variance equals the overall variance () and when the weights increase, this variance drops.

1.2 Gaussian processes for inversion and error estimation

Gaussian processes are non-parametric priors and a common machine learning framework that have been applied to both regression and classification tasks (

Rasmussen (2003)). Through this analysis we will observe the equivalence between cokriging inversion within the geostatistics domain and Gaussian process inversion.

Bayesian inversion corresponds to computing the posterior distribution of the parameters of interest given our observations . Reid et al. (2013) provide an overview of the application of inversion within Gaussian processes. Given an sensitivity matrix, , and as the covariance matrix, we may represent the covariance structures of our observations and parameters as:


where is the noise variance. These representations are identical to those of equations 2 and 3, respectively. From here, a predictive distribution can represented by a Gaussian prior as




Note that the matrix of weighting coefficients introduced in equation 4 is similarly represented by . The cokriging variance vector of equation 7 could be calculated as . Through these comparisons, we can clearly visualize the similarities between Gaussian process inversion and the cokriging inversion used in this work.

1.3 Forward modeling

The vertical component of the gravitational effect of a right rectangular prism with density ; limits , ; limits , ; and limits and is


where is the Newton’s gravitational constant, and


The tensor components for in equation (13) are given by (Forsberg, 1984; Li and Chouteau, 1998)




Please note that


which means that out of the six tensor components, only five are independent.

1.4 Gravity gradients inversion

Consider that there are n gravity data () and density values on rectangular prisms (), their relationship can be written in the matrix form as


where is the kernel and its elements are calculated using the formula in equation (13). For each tensor component of , similar expressions can be written. For example


where the elements of are calculated using the formula in equation (17). This means that stochastic inversion can be applied on each of the tensor components using the procedure explained in the previous section and error of the estimated model. Similarly, error of the estimated model can be calculated for each tensor component.

Combination of different tensor components can also be used by stochastic inversion. For example for combination of and we have


Please note that since based on equation (23) out of , and , only two are independent, and as such only two of them should be used in a combination.

2 Results & Dibscussion

We generate densities for m prisms by FFTMA simulation, which is designed for 3D Gaussian model with variable covariance model properties simulation (Le Ravalec et al., 2000). We applied a spherical variogram model with C kg/m and and m where C is the variogram sill, and are the variogram ranges in different directions. The 3D domain is divided into cubic prisms. The stochastic density model is shown in Figure 1. Using these density values, we calculate all gravity gradient components (Li and Chouteau (1998)) shown in Figure 2.

Figure 1: Stochastic density model simulated using FFTMA simulation method, with normalized density values. Exponential isotropic variogram with range of 10 m is utilized to generate the model.
Figure 2: Gravity tensor component quantities calculated from the stochastic density model.

Then cokriging variances for each of the components for all prisms are calculated using Eq. 7. These cokriging variances (error estimations) are shown in Figure 3 for a depth m. We can easily see the differences between estimation errors for all components. It should be noted the difference of cokriging variances between component layers close to the surface are minimal. However, it is still possible to distinguish between these tensor components. We should also point out that the patterns close to the edges are less liable because of the edge effect. The cokriging (error) variances for each case have a distinguishable structure. This results suggest that a gradient component, or combination of components, can be chosen lead to a more advantageous result in parameter estimation, depending on the given application. In Figure 4, we choose a grid in the center of the first layer below the surface of the horizontal layers in the z-direction. From these, the cokriging variances are drawn for several components and component combinations, as a function of depth for each prism. We can see that each component shares an upward trend in variance as the depth increases, however the and components show significantly reduced variance as compared to the other component combinations. These findings are even more prominent at shallower depths.

Figure 3: Cokriging variances shown in section Z = 5 m depth for six components and three component combinations. The variance values have been normalized.
Figure 4: Cokriging variances for several components and component combinations, as a function of depth for one cell in centre grid of each layer.

2.1 Evaluation of Hyperparameters

We discuss the effects of hyperparameter tuning on the cokriging variances for each kernel component. This includes analysis of the range (lengthscale), as well as rotation angles,

. Further investigations focus on the effects of manipulating the error covariance, or nugget effect.

2.1.1 Effects of Range

In order to investigate more, we select a grid at the center of the first horizontal layer below the surface. We change the variogram range, ax, from 1 meter to 37 meters, and calculate variances for several components and component combinations. The variances for different variogram ranges are shown in Figure 5. It is seen that different gradients have different behaviours. For example, if we look at and components, at lower ranges they are similar but at the higher ranges starts to show lower error variance. Note that for the combination of components, the effects of changing of changing range are less significant as compared to single components.

Figure 5: Cokriging variances of several components and component combinations as a function of increasing range. In this example, the increasing range is in the x direction.

2.1.2 Geometric Rotation

To further study the effects of geometry, the outcomes of changing the coordinate system rotation angle, , for a given prism have been measured. This experiment has been visualized in Figure 6, where changes about the z-axis at regular intervals between 0 and 90 degrees. The noise of single gravitational tensor component is sensitive to the angle of rotation. As expected, tensor components that do not contain elements of the rotated axis (, , ) remain unaffected. This sensitivity is significantly reduced when applied to a greater combination of components. These combinations of tensor components also yield significantly lower variance than any single component.

Figure 6: Cokriging variances of several components and component combinations as a function of increasing coordinate system rotation angle.

2.1.3 Study of Noise

Switzer (1993) discussed how the kriging error reflects two scales of variability:

  • Local modulation: simply indicates the spatial configuration of the data locally. For example, errors are smaller close to data points, clustered data carry less information, and so on.

  • Regional modulation: as we apply a variogram for the whole domain, where data allows, the variogram parameters are adjusted regionally, and the kriging variance, without being conditional, becomes an useful indicator of uncertainty.

In this manuscript we take advantage of both the regional modulation effect and the effect of the different tensor kernels to analyze the errors from gravity gradients.

For the same prism, the effect of adding noise using nugget effect is studied. The nugget effect is a phenomenon present in many regionalized variables and represents short scale randomness or noise in the regionalized variable. It can be seen graphically in the variogram plot as a discontinuity at the origin of the function.

Before we can continue with these experiments, we must account for any preprocessing that may amplify the noise of any affected gravitational tensor components. Pilkington (2014) studied the affect of noise and determined that regardless of the preprocessing method and input error, that tensors will maintain a static relative noise as compared to the component. Using a ratio of 1:0.59:0.37:0.7 for the noise levels for the , /, , and /, respectively, the postprocessing errors should be equivalent.

Using these weighted relative errors, the results are shown in Figure 7. As we expect, for all gradient components and combinations, the error variance increases by increasing noise. The component and some of the component combinations show more resistance against noise increment. This figure also suggests that nugget noise increment, in general, increases the variance in the recovered models.

Figure 7: Cokriging variances for several components and component combinations, as a function of increasing noise. Noise is added as a nugget effect.

3 Case Study: New Found Dataset

Here we showcase the same inversion concepts and apply them to real-world solutions. In 2020, HeliFALCON Airborne Gravity Gradiometer (AGG) survey data was collected and provided by New Found Gold Corp. (NFGC)’s Queensway project in Newfoundland and Labrador, Canada, henceforth referred to in this paper as the New Found dataset. Of this data, we selected a km area from the centre of this selection for use in this study. Using the stochastic inversion methodology as described above, we conduct the same analysis on this real-world data. Figure 8, shows the observed gravity gradient components for the selected area.

Figure 8: Component quantities observed in a selected area in New Found dataset.

Figure 9 shows the cokriging variances of different gravity kernels at the slice taken from a depth of 5 m. From these results we can clearly see the variance of a kernel combination is lower than that of a single component, most notably in the larger combinations. These results remain consistent with the observed results from the synthetic data in Figure 3.

Figure 9: Cokriging variances shown in section Z = 5 m depth for six components and three component combinations for the selected area in the New Found dataset.

Figure 10 shows the density of prisms in the fifth layer below the surface in Z direction.

Figure 10: The fifth layer of the recovered density model in Z direction on the selected area in the New Found dataset.

Using these stochastic inversion techniques, the generated density values can progress further through the processing pipeline for applications such as clustering and later structural interpretation.

4 Conclusion

We used cokriging variances from stochastic inversion (Shamsipour et al., 2010) as an indicator to quantitatively study different single tensor components. The indicator defines an error for each prism and the interpreter can easily select a single or combination of tensor components depending their purpose and available data. The method can simply be applied to study the effect of model geometry (in a stationary framework) or variation of the noise characteristics. For geometry studies, we analyze the effects of changing the variogram parameters: either ranges or rotations, as well as a study of the noise. We find that the combination tensor component errors are lower as compared to other components. We also find that error of combination tensor components remains independent of rotation, and continue to contain lower cokriging variances across all nugget effect values. These results remain consistent when the same inversion techniques are applied to the real New Found HeliFALCON AGG survey data.

5 Acknowledgment

We wish to thank Dr. Denis Marcotte for his valuable discussions throughout the preparation of this manuscript. The authors would also like to thank Dr. Mark Pilkington for his comments towards this work.


  • M. Beiki and L. B. Pedersen (2010) Eigenvector analysis of gravity gradient tensor to locate geologic bodies. Geophysics 75 (6), pp. I37–I49. Cited by: §1.
  • S. K. Jordan (1978) Moving-base gravity gradiometer surveys and interpretation. Geophysics 43 (1), pp. 94–101. Cited by: §1.
  • M. Le Ravalec, B. Noetinger, and L. Y. Hu (2000) The fft moving average (fft-ma) generator: an efficient numerical method for generating and conditioning gaussian simulations. Mathematical Geology 32 (6), pp. 701–723. Cited by: §2.
  • Y. Li (2001) 3-d inversion of gravity gradiometer data. In SEG Technical Program Expanded Abstracts 2001, pp. 1470–1473. Cited by: §1.
  • C. Martinez, Y. Li, R. Krahenbuhl, and M. A. Braga (2012) 3D inversion of airborne gravity gradiometry data in mineral exploration: a case study in the quadrilatero ferrifero, brazil. Geophysics 78 (1), pp. B1–B11. Cited by: §1.
  • V. Mikhailov, G. Pajot, M. Diament, and A. Price (2007) Tensor deconvolution: a method to locate equivalent sources from full tensor gravity data. Geophysics 72 (5), pp. I61–I69. Cited by: §1.
  • C. A. Murphy and J. Brewster (2007) Target delineation using full tensor gravity gradiometry data. ASEG extended abstracts 2007 (1), pp. 1–3. Cited by: §1.
  • D.E. Myers (1982) Matrix formulation of co-kriging. Mathematical Geology 14 (3), pp. 249–257. Cited by: §1.1.
  • M. Pilkington (2012) Analysis of gravity gradiometer inverse problems using optimal design measures. Geophysics 77 (2), pp. G25–G31. Cited by: §1.
  • M. Pilkington (2014) Evaluating the utility of gravity gradient tensor components. Geophysics 79 (1), pp. G1–G14. Cited by: §1, §2.1.3.
  • C. E. Rasmussen (2003) Gaussian processes in machine learning. In Summer school on machine learning, pp. 63–71. Cited by: §1.2.
  • A. S. Reid, S. O’Callaghan, E. Bonilla, L. McCalman, T. Rawling, and F. Ramos (2013) Bayesian joint inversions for the exploration of earth resources. In

    Twenty-Third International Joint Conference on Artificial Intelligence

    Cited by: §1.2.
  • P. Shamsipour, D. Marcotte, M. Chouteau, and P. Keating (2010) 3D stochastic inversion of gravity data using cokriging and cosimulation. Geophysics 75, pp. I1–I10. Cited by: §1, §4, §1.1.
  • P. Switzer (1993) The spatial variability of prediction errors. In Geostatistics Tróiaí92, pp. 261–272. Cited by: §2.1.3.
  • E. B. Tchikaya, M. Chouteau, P. Keating, and P. Shamsipour (2016) 3D unconstrained and geologically constrained stochastic inversion of airborne vertical gravity gradient data. Exploration Geophysics 47 (1), pp. 67–84. Cited by: §1.
  • D. Vasco (1989) Resolution and variance operators of gravity and gravity gradiometry. Geophysics 54 (7), pp. 889–899. Cited by: §1.
  • M. S. Zhdanov, R. Ellis, and S. Mukherjee (2004) Three-dimensional regularized focusing inversion of gravity gradient tensor component data. Geophysics 69 (4), pp. 925–937. Cited by: §1.