Correcting the Smoothing Effect of Ordinary Kriging Estimates

The smoothing effect of ordinary kriging is a well-known dangerous effect associated with this estimation technique. Consequently kriging estimates do not reproduce both histogram and semivariogram model of sample data. A four-step procedure for correcting the smoothing effect of ordinary kriging estimates is shown to be efficient for the reproduction of histogram and semivariogram without loss of local accuracy. Furthermore, this procedure provides a unique map sharing both local and global accuracies. Ordinary kriging with a proper correction for smoothing effect can be revitalized as a reliable estimation method that allows a better use of the available information.     

Estimation of semivariogram parameters and evaluation of the effects of data sparsity

Semivariogram parameters are estimated by a weighted least-squares method and a jackknife kriging method. The weighted least-squares method is investigated by differing the lag increment and maximum lag used in the fit. The jackknife kriging method minimizes the variance of the jackknifing error as a function of semivariogram parameters. The effects of data sparsity and the presence of a trend are investigated by using 400-, 200-, and 100-point synthetic data sets. When the two methods yield significantly different results, more data may be needed to determine reliably the semivariogram parameters, or a trend may be present in the data.     

Correcting the Smoothing Effect of Estimators: A Spectral Postprocessor

The postprocessing algorithm introduced by Yao for imposing the spectral amplitudes of a target covariance model is shown to be efficient in correcting the smoothing effect of estimation maps, whether obtained by kriging or any other interpolation technique. As opposed to stochastic simulation, Yao's algorithm yields a unique map starting from an original, typically smooth, estimation map. Most importantly it is shown that reproduction of a covariance/semivariogram model (global accuracy) is necessarily obtained at the cost of local accuracy reduction and increase in conditional bias. When working on one location at a time, kriging remains the most accurate (in the least squared error sense) estimator. However, kriging estimates should only be listed, not mapped, since they do not reflect the correct (target) spatial autocorrelation. This mismatch in spatial autocorrelation can be corrected via stochastic simulation, or can be imposed a posteriori via Yao's algorithm.     

Fitting the Linear Model of Coregionalization by Generalized Least Squares

In geostatistical studies, the fitting of the linear model of coregionalization (LMC) to direct and cross experimental semivariograms is usually performed with a weighted least-squares (WLS) procedure based on the number of pairs of observations at each lag. So far, no study has investigated the efficiency of other least-squares procedures, such as ordinary least squares (OLS), generalized least squares (GLS), and WLS with other weighing functions, in the context of the LMC. In this article, we compare the statistical properties of the sill estimators obtained with eight least-squares procedures for fitting the LMC: OLS, four WLS, and three GLS. The WLS procedures are based on approximations of the variance of semivariogram estimates at each distance lag. The GLS procedures use a variance–covariance matrix of semivariogram estimates that is (i) estimated using the fourth-order moments with sill estimates (GLS₁), (ii) calculated using the fourth-order moments with the theoretical sills (GLS₂), and (iii) based on an approximation using the correlation between semivariogram estimates in the case of spatial independence of the observations (GLS₃). The current algorithm for fitting the LMC by WLS while ensuring the positive semidefiniteness of sill matrix estimates is modified to include any least-squares procedure. A Monte Carlo study is performed for 16 scenarios corresponding to different combinations of the number of variables, number of spatial structures, values of ranges, and scale dependence of the correlations among variables. Simulation results show that the mean square error is accounted for mostly by the variance of the sill estimators instead of their squared bias. Overall, the estimated GLS₁ and theoretical GLS₂ are the most efficient, followed by the WLS procedure that is based on the number of pairs of observations and the average distance at each lag. On that basis, GLS₁ can be recommended for future studies using the LMC.     

On the equivalence of kriging and maximum entropy estimators

This study compares kriging and maximum entropy estimators for spatial estimation and monitoring network design. For second-order stationary random fields (a subset of Gaussian fields) the estimators and their associated interpolation error variances are identical. Simple lognormal kriging differs from the lognormal maximum entropy estimator, however, in both mathematical formulation and estimation error variances. Two numerical examples are described that compare the two estimators. Simple lognormal kriging yields systematically higher estimates and smoother interpolation surfaces compared to those produced by the lognormal maximum entropy estimator. The second empirical comparison applies kriging and entropy-based models to the problem of optimizing groundwater monitoring network design, using six alternative objective functions. The maximum entropy-based sampling design approach is shown to be the more computationally efficient of the two.     

On unbiased backtransform of lognormal kriging estimates

Lognormal kriging is an estimation technique that was devised for handling highly skewed data distributions. This technique takes advantage of a logarithmic transformation that reduces the data variance. However, backtransformed lognormal kriging estimates are biased because the nonbias term is totally dependent on a semivariogram model. This paper proposes a new approach for backtransforming lognormal kriging estimates that not only presents none of the problems reported in the literature but also reproduces the sample histogram and, consequently, the sample mean.     

The Second-Order Stationary Universal Kriging Model Revisited

Universal kriging originally was developed for problems of spatial interpolation if a drift seemed to be justified to model the experimental data. But its use has been questioned in relation to the bias of the estimated underlying variogram (variogram of the residuals), and furthermore universal kriging came to be considered an old-fashioned method after the theory of intrinsic random functions was developed. In this paper the model is reexamined together with methods for handling problems in the inference of parameters. The efficiency of the inference of covariance parameters is shown in terms of bias, variance, and mean square error of the sampling distribution obtained by Monte Carlo simulation for three different estimators (maximum likelihood, bias corrected maximum likelihood, and restricted maximum likelihood). It is shown that unbiased estimates for the covariance parameters may be obtained but if the number of samples is small there can be no guarantee of good estimates (estimates close to the true value) because the sampling variance usually is large. This problem is not specific to the universal kriging model but rather arises in any model where parameters are inferred from experimental data. The validity of the estimates may be evaluated statistically as a risk function as is shown in this paper.     

Quadratic covariance estimation and equivalence of predictions

This article illustrates the use of linear and nonlinear regression models to obtain quadratic estimates of covariance parameters. These models lead to new insights into the motivation behind estimation methods, the relationships between different methods, and the relationship of covariance estimation to prediction. In particular, we derive the standard estimating equations for minimum norm quadratic unbiased translation invariant estimates (MINQUEs) from an appropriate linear model. Connections between the linear model, minimum variance quadratic unbiased translation invariant estimates (MIVQUEs), and MINQUEs are examined and we provide a minimum norm justification for the use of one-step normal theory maximum likelihood estimates. A nonlinear regression model is used to define MINQUEs for nonlinear covariance structures and obtain REML estimates. Finally, the equivalence of predictions under various models is examined when covariance parameters are estimated. In particular, we establish that when using MINQUE, iterative MINQUE, or restricted maximum likelihood (REML) estimates, the choice between a stationary covariance function and an intrinsically stationary semivariogram is irrelevant to predictions and estimated prediction variances.     

An Alternative Measure of the Reliability of Ordinary Kriging Estimates

This paper presents an interpolation variance as an alternative to the measure of the reliability of ordinary kriging estimates. Contrary to the traditional kriging variance, the interpolation variance is data-values dependent, variogram dependent, and a measure of local accuracy. Natural phenomena are not homogeneous; therefore, local variability as expressed through data values must be recognized for a correct assessment of uncertainty. The interpolation variance is simply the weighted average of the squared differences between data values and the retained estimate. Ordinary kriging or simple kriging variances are the expected values of interpolation variances; therefore, these traditional homoscedastic estimation variances cannot properly measure local data dispersion. More precisely, the interpolation variance is an estimate of the local conditional variance, when the ordinary kriging weights are interpreted as conditional probabilities associated to the n neighboring data. This interpretation is valid if, and only if, all ordinary kriging weights are positive or constrained to be such. Extensive tests illustrate that the interpolation variance is a useful alternative to the traditional kriging variance.     

Combining Robustness with Efficiency in?the?Estimation of the Variogram

In the present paper, we propose a new method for the estimation of the variogram, which combines robustness with efficiency under intrinsic stationary geostatistical processes. The method starts by using a robust estimator to obtain discrete estimates of the variogram and control atypical observations that may exist. When the number of points used in the fit of a model is the same as the number of parameters, ordinary least squares and generalized least squares are asymptotically equivalent. Therefore, the next step is to fit the variogram by ordinary least squares, using just a few discrete estimates. The procedure is then repeated several times with different subsets of points and this produces a sequence of variogram estimates. The final estimate is the median of the multiple estimates of the variogram parameters. The suggested estimator will be called multiple variograms estimator. This procedure assures a global robust estimator, which is more efficient than other robust proposals. Under the assumed dependence structure, we prove that the multiple variograms estimator is consistent and asymptotically normally distributed. A simulation study confirms that the new method has several advantages when compared with other current methods.     

Empirical Maximum Likelihood Kriging: The General Case

Although linear kriging is a distribution-free spatial interpolator, its efficiency is maximal only when the experimental data follow a Gaussian distribution. Transformation of the data to normality has thus always been appealing. The idea is to transform the experimental data to normal scores, krige values in the “Gaussian domain” and then back-transform the estimates and uncertainty measures to the “original domain.” An additional advantage of the Gaussian transform is that spatial variability is easier to model from the normal scores because the transformation reduces effects of extreme values. There are, however, difficulties with this methodology, particularly, choosing the transformation to be used and back-transforming the estimates in such a way as to ensure that the estimation is conditionally unbiased. The problem has been solved for cases in which the experimental data follow some particular type of distribution. In general, however, it is not possible to verify distributional assumptions on the basis of experimental histograms calculated from relatively few data and where the uncertainty is such that several distributional models could fit equally well. For the general case, we propose an empirical maximum likelihood method in which transformation to normality is via the empirical probability distribution function. Although the Gaussian domain simple kriging estimate is identical to the maximum likelihood estimate, we propose use of the latter, in the form of a likelihood profile, to solve the problem of conditional unbiasedness in the back-transformed estimates. Conditional unbiasedness is achieved by adopting a Bayesian procedure in which the likelihood profile is the posterior distribution of the unknown value to be estimated and the mean of the posterior distribution is the conditionally unbiased estimate. The likelihood profile also provides several ways of assessing the uncertainty of the estimation. Point estimates, interval estimates, and uncertainty measures can be calculated from the posterior distribution.     

Interpolation and mapping of probabilities for geochemical variables exhibiting spatial intermittency

In monitoring a minor geochemical element in groundwater or soils, a background population of values below the instrumental detection limit is frequently present. When those values are found in the monitoring process, they are assigned to the detection limit which, in some cases, generates a probability mass in the probability density function of the variable at that value (the minimum value that can be detected). Such background values could distort both the estimation of the variable at nonsampled locations and the inference of the spatial structure of variability of the variable. Two important problems are the delineation of areas where the variable is above the detection limit and the estimation of the magnitude of the variables inside those areas. The importance of these issues in geochemical prospecting or in environmental sciences, in general related with contamination and environmental monitoring, is obvious. In this paper the authors describe the two-step procedure of indicator kriging and ordinary kriging and compare it with empirical maximum likelihood kriging. The first approach consists of using a binary indicator variable for estimating the probability of a location being above the detection limit, plus ordinary kriging conditional to the location being above the detection limit. An estimation variance, however, is not available for that estimator. Empirical maximum likelihood kriging, which was designed to deal with skew distributions, can also deal with an atom at the origin of the distribution. The method uses a Bayesian approach to kriging and gives intermittency in the form of a probability map, its estimates providing a realistic assessment of their estimation variance. The pros and cons of each method are discussed and illustrated using a large dataset of As concentration in groundwater. The results of the two methods are compared by cross-validation.     

Effects of search parameters on kriged reserve estimates

Summary Reliable ore reserve estimates for deposits with highly skewed grade distributions are difficult tasks to perform. Although some recent geostatistical techniques are available to handle problems with these estimations, ordinary kriging or conventional interpolation methods are still widely used to estimate the ore reserves for such deposits. The estimation results can be very sensitive to the search parameters used during the interpolation of grades with these methods.This paper compares the ore reserve estimates from ordinary kriging using several cases in which certain search parameters are varied. The comparisons are extended to different mineralizations to show the changing effects of these parameters.     

Estimation or simulation? That is the question

The issue of smoothing in kriging has been addressed either by estimation or simulation. The solution via estimation calls for postprocessing kriging estimates in order to correct the smoothing effect. Stochastic simulation provides equiprobable images presenting no smoothing and reproducing the covariance model. Consequently, these images reproduce both the sample histogram and the sample semivariogram. However, there is still a problem, which is the lack of local accuracy of simulated images. In this paper, a postprocessing algorithm for correcting the smoothing effect of ordinary kriging estimates is compared with sequential Gaussian simulation realizations. Based on samples drawn from exhaustive data sets, the postprocessing algorithm is shown to be superior to any individual simulation realization yet, at the expense of providing one deterministic estimate of the random function.     

Kriging and Semivariogram Deconvolution in the Presence of Irregular Geographical Units

This paper presents a methodology to conduct geostatistical variography and interpolation on areal data measured over geographical units (or blocks) with different sizes and shapes, while accounting for heterogeneous weight or kernel functions within those units. The deconvolution method is iterative and seeks the point-support model that minimizes the difference between the theoretically regularized semivariogram model and the model fitted to areal data. This model is then used in area-to-point (ATP) kriging to map the spatial distribution of the attribute of interest within each geographical unit. The coherence constraint ensures that the weighted average of kriged estimates equals the areal datum.This approach is illustrated using health data (cancer rates aggregated at the county level) and population density surface as a kernel function. Simulations are conducted over two regions with contrasting county geographies: the state of Indiana and four states in the Western United States. In both regions, the deconvolution approach yields a point support semivariogram model that is reasonably close to the semivariogram of simulated point values. The use of this model in ATP kriging yields a more accurate prediction than a naïve point kriging of areal data that simply collapses each county into its geographic centroid. ATP kriging reduces the smoothing effect and is robust with respect to small differences in the point support semivariogram model. Important features of the point-support semivariogram, such as the nugget effect, can never be fully validated from areal data. The user may want to narrow down the set of solutions based on his knowledge of the phenomenon (e.g., set the nugget effect to zero). The approach presented avoids the visual bias associated with the interpretation of choropleth maps and should facilitate the analysis of relationships between variables measured over different spatial supports.     

Evaluation and comparison of spatial interpolators II

The performance of several variations on ordinary kriging and inverse distance estimators is evaluated. Mean squared errors (MSE) were calculated for estimates made on multiple resamplings from five exhaustive data bases representing two distinctly different types of estimation problem. Ordinary kriging, when performed with variograms estimated from the sample data, was more robust than inverse-distance methods to the type of estimation problem, and to the choice of estimation parameters such as number of neighbors.Notice: Although the research described in this article has been funded in part by the United States Environmental Protection Agency through Cooperative Agreement CR818526 to the Harry Reid Center for Environmental Studies, University of Nevada-Las Vegas, it has not been subjected to Agency review. Therefore it does not necessarily reflect the views of the Agency. Mention of trade names or commercial products does not constitute endorsement or recommendation for use.     

Kriging and Semivariogram Deconvolution in the Presence of Irregular Geographical Units

This paper presents a methodology to conduct geostatistical variography and interpolation on areal data measured over geographical units (or blocks) with different sizes and shapes, while accounting for heterogeneous weight or kernel functions within those units. The deconvolution method is iterative and seeks the pointsupport model that minimizes the difference between the theoretically regularized semivariogram model and the model fitted to areal data. This model is then used in area-to-point (ATP) kriging to map the spatial distribution of the attribute of interest within each geographical unit. The coherence constraint ensures that the weighted average of kriged estimates equals the areal datum.This approach is illustrated using health data (cancer rates aggregated at the county level) and population density surface as a kernel function. Simulations are conducted over two regions with contrasting county geographies: the state of Indiana and four states in the Western United States. In both regions, the deconvolution approach yields a point support semivariogram model that is reasonably close to the semivariogram of simulated point values. The use of this model in ATP kriging yields a more accurate prediction than a na?ve point kriging of areal data that simply collapses each county into its geographic centroid. ATP kriging reduces the smoothing effect and is robust with respect to small differences in the point support semivariogram model. Important features of the point-support semivariogram, such as the nugget effect, can never be fully validated from areal data. The user may want to narrow down the set of solutions based on his knowledge of the phenomenon (e.g., set the nugget effect to zero). The approach presented avoids the visual bias associated with the interpretation of choropleth maps and should facilitate the analysis of relationships between variables measured over different spatial supports.     

A parametric study of robustness of kriging variance as a function of range and relative nugget effect for a spherical semivariogram

In geostatistics, an estimation of blocks of a deposit is reported along with the variance of error made in their estimation. This calculation is based on the model chosen for the semivariogram of the deposit so that mistakes in its estimation can manifest themselves in the perception of accuracy with which blocks are known. Changes in kriging variance resulting from various amounts of error in modeling the relative nugget effect and range of the semivariogram are investigated for an extensive set of spherical semivariograms.     

Compensating for estimation smoothing in kriging

Smoothing is a characteristic inherent to all minimum mean-square-error spatial estimators such as kriging. Cross-validation can be used to detect and model such smoothing. Inversion of the model produces a new estimator—compensated kriging. A numerical comparison based on an exhaustive permeability sampling of a 4-ft² slab of Berea Sandstone shows that the estimation surface generated by compensated kriging has properties intermediate between those generated by ordinary kriging and stochastic realizations resulting from simulated annealing and sequential Gaussian simulation. The frequency distribution is well reproduced by the compensated kriging surface, which also approximates the experimental semivariogram well—better than ordinary kriging, but not as well as stochastic realizations. Compensated kriging produces surfaces that are more accurate than stochastic realizations, but not as accurate as ordinary kriging.     

A robust cylindrical fitting to point cloud data

Environmental, engineering and industrial modelling of natural features (e.g. trees) and man-made features (e.g. pipelines) requires some form of fitting of geometrical objects such as cylinders, which is commonly undertaken using a least-squares method that—in order to get optimal estimation—assumes normal Gaussian distribution. In the presence of outliers, however, this assumption is violated leading to a Gaussian mixture distribution. This study proposes a robust parameter estimation method, which is an improved and extended form of vector algebraic modelling. The proposed method employs expectation maximisation and maximum likelihood estimation (MLE) to find cylindrical parameters in case of Gaussian mixture distribution. MLE computes the model parameters assuming that the distribution of model errors is a Gaussian mixture corresponding to inlier and outlier points. The parameters of the Gaussian mixture distribution and the membership functions of the inliers and outliers are computed using an expectation maximisation algorithm from the histogram of the model error distribution, and the initial guess values for the model parameters are obtained using total least squares. The method, illustrated by a practical example from a terrestrial laser scanning point cloud, is novel in that it is algebraic (i.e. provides a non-iterative solution to the global maximisation problem of the likelihood function), is practically useful for any type of error distribution model and is capable of separating points of interest and outliers.