Prediction

If prediction locations are defined in the command file, gstat chooses a prediction method depending on the model defined by the complete set of commands in a command file.

When no variograms are specified, inverse distance weighted interpolation is the default action (Fig. 2.1, example 3).

When variograms are specified the default prediction method is ordinary kriging [13,5] (example 4 and example 8).

Simple kriging is the default action when in addition for each variable the simple kriging mean (sk_mean or b) is set (section 4.2; example 5, universal kriging or uncorrelated linear model prediction is used when a model for the trend is defined ,section 2.7). Multiple prediction, multivariable prediction, and stratified prediction are described in section 3.1-3.4. Prediction of block averages is described in section 3.5.

If the prediction locations are specified as a mask map with the command

mask: 'file';

then predictions and prediction variances are written to output maps only when these maps are specified explicitly (section 4.1; example 5).

As an alternative to prediction on grid map locations, prediction on non-gridded locations is the default action when these locations are specified with the

data(): ... ;

command (note the absence of an identifier between the parentheses). In this case, output is written in ascii table or simplified GeoEAS format to the file defined by the command set output=' file'; (example 4, or defined with the command line option -o, section 5.2).

Figure 2.3: Local neighbourhood selections. Lines indicate selected points (+). Lower right: variable with anisotropic variogram having a strong north-south correlation

Local Neighbourhoods

By default, gstat uses global prediction, meaning that for each prediction all data values are used. However, it is often desirable to use not all data values, but only a subset in a (spatial) neighbourhood around the prediction (simulation) location, for either computational reasons or the wish to assume first-order stationarity only locally. Gstat allows local neighbourhood selections to be based on distance (radius), number of data points (max, min), variogram distance (vdist), and number of data points per octant (3D) or quadrant (2D) (omax). The options are explained below (see also Fig. 2.3, section 4.2 examples in chapter 6).

The quadtree-based algorithm used to obtain data points in a local search neighbourhood is described in [11], and is found at http://www.cs.umd.edu/ brabec/quadtree/index.html (Bucket PR Quadtree demo).

radius = 10

select all data points within 10 (euclidian) distance units from the prediction location

max = 8

select the 8 data points that are closest (in euclidian distance) to the prediction location (or take all data points if less than 8 are available)

Some options should be combined, and permitted combinations are explained below. (Combinations not mentioned might result in unexpected or undesired results.)

radius = 10, max = 8

after selecting all data points at (euclidian) distances from the prediction location less or equal to 10, choose the 8 closest when more than 8 are found

radius = 10, max = 8, min = 4

in addition to the previous selection, generate a missing value if less than 4 points are found within the search radius 10

radius = 10, max = 8, min = 4, force

in addition to the previous selection, if less than 4 data points are found in the search radius, instead of generating a missing value, select (force) the 4 nearest (in euclidian distance) data points, regardless their distance

radius = 10, omax = 2

after selecting all data points at distances less or equal to 10, choose the 2 closest data points in each octant (3D), quadrant (2D) or secant (1D)

radius = 10, vdist, ...

after the radius selection, decide what the nearest data points are on the base of point-to-point semivariance of the data variable instead of euclidian distances (``semivariance distance'': in case of anisotropy this allows the prevalence of more correlated points over the, in the euclidian sense, nearest points)

Indicator kriging

Basically, indicator kriging is equivalent to simple or ordinary kriging of indicator-transformed data. However, resulting estimates of indicator values are not guaranteed to satisfy order relations. During indicator kriging, gstat will do order relation violation correction for independent, cumulative or categorical (disjunct) indicators only if the order is to one of the values in Table 2.1, order and [8,10].