Hi Julio,

It's something that I also questioned once about. Finally I finished to implement the regularization. However it would be interesting to design an optimization scheme that dropout and reinsert hard data based on the condition number along the simulation (regulating in some way maybe the neighboorhood) ...and that at the same time solve also the bias issue introduced by the simulation path.

For the last problem, I suggest you this recent research:

https://goo.gl/42BCsm

https://goo.gl/tJrL86

Bests ! 

Hi Gionata, thank you for sharing the papers, but the question doesn't relate to the simulation path. Regarding the papers you linked, I did read them when they came out, it is a good work by Raphael and co-authors.

There are several possible reasons why the coefficient matrix (for the kriging equations) may not be invertible.

1. The sufficient condition for an invertible matrix is that the covariance function is strictly positive definite (not just semi-definite) or in the case of ordinary kriging that the variogram is strictly conditionally negative definite. If you are using one of the well known software programs this problem should not occur since they limit the covariance function/variogram choices to functions satisfying the required condition.

2. Although the Gaussian covariance function is strictly positive definite if you look at a graph of the function you will see that the function is very nearly constant for a considerable distance (until it starts to curve downward). This means if you have multiple pairs of data locations that are close together you have multiple entries in the coefficient matrix with the same or very nearly the same value (in variogram form you will have multiple zeros or multiple values very close to zero. On simple solution is to incorporate a small nugget effect.

Another solution is to delete one data location (i.e. suppress) in each pair of close locations.

Unlike IWD kriging not only takes into account the spatial pattern of data locations relative to the estimation location but also the spatial pattern of the data locations with respect to each other.

From the perspective of the kriging weights the kriging system will down weight one of the locations (depending on which is closer to the estimation location). The value at the second location (in the pair) is either duplicate information or it is dis-information. Based on the (Gaussian) correlation function the data values must be very highly correlated, i.e. nearly the same. If the data values are not essentially the same then  one of them is "dis-information" and maybe you need to look at the data to see why this might have occurred, e.g. one value is simply wrong.

Dear Prof. Myers,

Thank you for your detailed answer, I appreciate it. My main concern with adding a nugget effect, or equivalently, adding a term s^2 to the diagonal of the system, is that we loose the great interpolation property of the estimator. This addition is also equivalent to assuming noisy observations ~ Normal(0, s^2).

Since most traditional software for geostatistics such as GSLIB and SGeMS implement this preconditioning step, I will end up implement it as well in the next releases of the project. However, I hope to be able to revisit this issue in the future (there must be a clever and efficient way to avoid ill-conditioned Kriging systems based on the spatial configuration of the data inside of the neighborhood).

The comment about the Gaussian covariance specifically is also very pertinent. This is the model that produces the most ill-conditioned matrices from experimentation with the code.

Julio

You are correct that ONE interpretation of a nugget component in the covariance function/variogram is that the data is noisy. The noise term need not be normal and to write s^2 is a little misleading. This would normally correspond to the sample variance rather than to a model variance. More generally a nugget component appears because there are not enough (perhaps nearly none) pairs of data locations with very short lag distances to be able to adequately model the behavior of the covariance function/variogram close to zero. The variogram must theoretically be zero at zero distance (the covariance function must be equal to the sill at distance zero). . Theoretically a nugget corresponds to a discontinuity of the covariance function/variogram at zero. In particular it is not differentiable at zero.

If the data really is noisy and you have a separate way to estimate the error variance (e.g. the noise results from measurement/instrument error) then you can interpolate the non-noisy data using cokriging

You will find an extensive discussion of "nugget effects" in "Geostatistics: Modeling spatial uncertainty",  J-P. Chiles and P. Delfiner (2nd edition) J. Wiley & sons

If you look at the ordinary kriging equations with a variogram model incorporating a nugget component you will see that the nugget component changes the Lagrange multipler but does not change the kriging weights.

I am not sure which "preconditioning" you are referring to but I am dubious about any claims of "preconditioning" the data. e.g. transforming to a normal distribution. A linear transformation essentially has no effect, i.e. it is invertible.

Note particularly Chiles and Delfiner's discussion of the nugget effect and regularization and also the relationship to possible data location coordinate errors

Dear Prof. Myers,

I will definitely read the book again before I commit to a particular implementation. The preconditioning of the Kriging system, however; is due to *numerical* issues. They don't relate to the theory. I wish I could upload some visualizations and data here, but this website is very limited.

Other software such as GSLIB and SGeMS have this preconditioning step implemented. By adding the extra term to the diagonal of the covariance, the condition number of the matrix is decreased with the downside of approximating the solution weights. In practice, this approximation really doesn't matter, but it is still an approximation that I was trying to avoid.

Hi Jin,

How does your answer addresses my question? I don't see any connection, sorry.

Jin Li

How would you define :"predictive accuracy:"?  The most likely answer is you must compare the predicted value with the true value but if the true value is known you don't need a predictor or interpolator. Hence it is nonsensical to talk about "improving the predictive accuracy"

Jin Li

It is incorrect to use cross validation to quantify "predictive accuracy", you are not understanding the theory behind cross validation. The predictive accuracy algorithm listed in your SPM package is not cross-validation as it appears in the geostatistics literature. IDW is not a geostatistical method, it is not based on any theory or statistical assumptions. It does not reproduce the data values at the data locations. Random Forest can be used for spatial regression but it does not reproduce the data values.

1991, Myers,D.E., On Variogram Estimation in Proceedings of the First Inter. Conf. Stat. Comp., Cesme, Turkey, 30 Mar.-2 April 1987, Vol II, American Sciences Press, 261-281

Cross validation might be used to identify "unusual" data values/data locations, to compare different variogram model types, to examine the sensitivity of the kriging results to the choice of the search neighborhood parameters, to examine the sensitivity of the kriging results to changes in the variogram model parameters, examine whether the data values are from two different populations. But it does not quantify "predictive accuracy"

I suggest you look at the book "Geostatistics" by J.-P. Chiles and P. Delfiner (2nd edition), John Wiley & Sons.

Also "Mining Geostatstics" by A. Journel and Ch Huibregts , Academic Press (1978)

For something more elementary but which also discusses cross validation see

1991, Myers,D.E.Interpolation and Estimation with Spatially Located Data Chemometrics and Intelligent Laboratory Systems 11, 209-228

The original question pertained to "ill-conditioned coefficient kriging matrices" not to "predictive accuracy". The first condition that must be checked is whether the variogram model is strictly conditionally negative definite, if it is then the kriging matrix will be invertible. However one commonly used model, the Gaussian, is not only zero at distance zero it is very close ro zero for moderately non-zero distances (this function is extremely smooth near the origin). This is the situation mentioned by Julio. Actually Julio was asking about sequential simulation which uses simple kriging (and hence a covariance function instead of a variogram). Seq Simulation also is based on an assumption of a Gaussian distribution, not to be confused with a Gaussian model variogram/covariance. There is at least one geostatistics package that allows the user to make finer decisions about the kriging search neighborhoods and hence when there are two locations too close to each other to exclude one of them from the search neighborhood. See GeoEAS (US EPA), this software was used in the Chemometrics and Intelligent Lab Systems paper.