7

I'm using Ordinary Kriging to interpolate a set of spatial data. I'm working on different sets which represent different days of measurement. I get good results, apart from one specific set. Will post some examples:

Good interpolations:

Good interpolation

Good interpolation

Strange interpolation:

enter image description here

What might be causing this behaviour? Why is it going light-blue on almost all the map, without spatially expanding lows and highs (which are pretty clustered!)? What parameters should I fix, or what could be the dataset error behind?

I see a pretty evident pattern, in that values tend to get high in the top-right corner, but it is not reflected by the interpolating surface.

Thoughts:

There could be some bug in my code, but as said, it goes well for all other sets except this one. I've tried applying different variograms but the results is still silly. On the other hand, other algorithms (e.g. inverse distance) do not show any issue.

Other info:

I'm using R, but could switch to GIS if needed. For R folks, I'm just calling:

kriging <- autoKrige(value~1,finalDataset,new_data=DEM)

from the automap package.

Sample data

I uploaded a working case and the problematic one (hopefully).

Improve this question
4
  • It seems to be a R-related issue, because interpolating in ArcGIS gives good results. Anyway I'm forced to go with R for work reasons, so I'd like to know what's wrong there.
    – natario
    May 18, 2015 at 16:12
  • 1
    My only guess is that that kind of "dipole" near the center (a low and a high tightly coupled) is telling R that each point has to have a very local influence. But then.. same set is working in ArcGIS.
    – natario
    May 18, 2015 at 16:23
  • 1
    Why don't you take the model variogram parameters from ArcGIS and apply them to your R model? You can then compare the actual interpolation results. You model looks overfit and I imagine that there is something in your data that is causing the autoKrige function to misspecify the model parameters. Perhaps autoKrige is fitting a pure nugget effect. Did you transform your data? Fitting Kriging models is a bit of an art and one should not rely entirely on default parameters. Take a look at the eyefit function in geoR. It will let you visualize how different parameters change your variogram.
    – Jeffrey Evans
    May 18, 2015 at 16:27
  • @Jeffrey ArcGIS can't be part of our daily workflow sadly. I need to define a Kriging model that could fit to present and future datasets, but right now I cannot imagine what parameters to change to make the dotted-figure more pleasing. Thank you for the eyefit suggestion, will check it out.
    – natario
    May 19, 2015 at 12:15

2 Answers 2

Reset to default
9

In a nutshell, the problem lies in a mismatch between data behavior and some (strong) assumptions you are implicitly making.

Diagnosis

The strongest of these is that the data are one realization of a second-order stationary process. They clearly are not, as you can tell by comparing the region near (450000, 5075000) in the upper "neck" (which I will call "Region X") to other regions within the map. In all other regions, the neighbors of any data support tend to have comparable values. In Region X, however, there are many wild swings between high and low values over short distances. (This is a manifestation of a strong change in the second-order structure.) These swings will inflate the estimates of semivariance at short lags, causing the automatic variogram estimator to infer that the range is short. This range is unlikely to be any greater than a typical distance between an upper and lower extreme point. It is approximately the radius of the colored dots appearing in the map.

One property of Ordinary Kriging is that when all sample points are beyond the range of the variogram, its estimate will be the neighborhood mean. Evidently you have used a fairly large neighborhood--perhaps the entire dataset. This results in many of the kriged values equaling the dataset mean, shown in cyan. The kriged values will vary from this only at grid cells within the range of a data point.

In addition, somehow the automatic fit included almost no nugget effect. This probably was due to mixing in the contributions from the rest of the dataset, where spatial variation is locally small and smooth looking. In particular, there are clusters of tightly spaced points (such as those along the middle left boundary) that would steer the variogram down towards zero at short lags. The absence of any appreciable nugget effect forces the kriged surface to pass through the data points.

These two characteristics of the variogram--short range and small or nonexistent nugget--cause the acned, multiple-bull's-eye pattern to emerge, since the kriged values must change smoothly from the cyan background to each data value within very short distances of each data support.

Recommendations

If you were to change the variogram parameters, then you would wind up smoothing over the variation in Region X. This would hide the interesting and possibly important change of data behavior in Region X.

A more sophisticated approach might combine some of the following elements:

  • Careful, supervised variography, including cross-validation of the variogram and use of directional and mixed variograms.

  • Partition of the study area into a priori identifiable regions of stationary behavior.

  • Co-kriging with the other variables.

  • Universal kriging to incorporate explanatory covariates.

  • Initial exploratory data analysis to check the validity of all the data, especially those in Region X.

  • Adopting more complex models of spatial covariation.

With the limited information available, it would be irresponsible to recommend any particular course of action, since it would depend on why you are kriging, what the data represent, and how you intend to interpret the maps.

Improve this answer
4
  • Thank you for your answer. I confirm I get better results by narrowing the neighborhood, but of course the variogram remains the same. These measures represent rainfall, and to now we are just trying to get eye-pleasing results. We are trying to merge rain gauge data with SWE (snow water equivalent) measures, using different techniques; I have implemented some initial data analysis to filter out invalid points, but some "dipoles" remain. Would like to hear some suggestion from you.
    – natario
    May 19, 2015 at 12:10
  • If you restrict the analysis to areas outside Region X, the variogram should change dramatically: in particular, its range should lengthen. If all you want is "eye-pleasing results," then you shouldn't be using a complicated, computing-intensive method like kriging. Consider a smoothing spline instead.
    – whuber
    May 19, 2015 at 12:13
  • I wanted to say there's a methodological issue behind, that lies in the fact each dataset is created from two different sources/physical phenomena (-> we get dipoles). Before using the maps I surely need to get a non-dotted representation. The problem I see with your suggestion is that I might encounter a new dataset where region X stands in a different place, while I need to define a model valid for all the sets we have. However, will think about it. Thank you.
    – natario
    May 19, 2015 at 12:26
  • I would suggest taking a different approach, one that has some valid statistical theory to support it: when you have different kinds of measurement, it is rarely appropriate or effective to lump them together into the same dataset. Instead, consider regression-based approaches or even cokriging.
    – whuber
    May 19, 2015 at 12:31
4

You said you used "R" but R is not a geostatistical program, it is more like a programming language. However there are geostatistical "packages" in R, e.g. gstat. Which package did you use?

Side notes 1. The smoothing spline can be replicated by using cokriging, Ordinary kriging & universal kriging are equivalent to Radial Basis function interpolation and interpolation with a spline is a special case of RBF

  1. You don't need or want more than 15-20 data locations in any search neigborhood, it can degrade the results and/or lead to ill-conditioned kriging matrices.

  2. If the variogram changes with different regions then you are violating a fundamental assumption of kriging, i.e. intrinsic stationarity. You may need to split the region into subregions and do the analyses separately for each subregion.

  3. It is likely you need to do more exploratory statistical analysis both before variogram modeling and kriging as well as afterwards. Relying on an automatic variogram fitting routine is not a good idea if you don't have a better understanding of your data.

Improve this answer
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.