MVARIOGRAM procedure

Fits models to an experimental variogram (S.A. Harding, D.A. Murray & R. Webster).

Options

PRINT = string tokens Controls printed output from the fit (model, summary, estimates, correlations, fittedvalues, monitoring); default mode, summ, esti
MODELTYPE = string token Defines which model to fit (power, boundedlinear, circular, spherical, doublespherical, pentaspherical, exponential, besselk1, gaussian, affinepower, linear, cubic, stable, cardinalsine, matern); default powe
WEIGHTING = string token Method to be used for weighting (counts, cbyvar, equal); default coun
CONSTANT = string token How to treat the constant (estimate, omit); default esti
SMOOTHNESS = scalar Value of power parameter for the stable model, or ν parameter for the Matern model; default * i.e. estimate
ISOTROPY = string token Defines whether to fit an isotropic or geometrical anisotropic model (isotropic, geometrical); default isot
WINDOW = scalar Window in which to plot a graph; default 0 i.e. no graph
TITLE = text Title for the graph
XUPPER = scalar Upper limit for the x-axis in the graph
PENDATA = scalar Pen to be used to plot the data; default 1
PENMODEL = scalar Pen to be used to plot the model; default 2

Parameters

VARIOGRAM = variates or matrices Experimental variogram to which the model is to be fitted, as a variate if in only one direction or as a matrix if there are several
COUNTS = variates or matrices Counts for the points in each variogram (not required if WEIGHTING=equal)
DISTANCE = variates or matrices Mean lag distances for the points in each variogram
DIRECTION = variates Directions in which each variogram was computed
INITIAL = scalars or variates Scalar defining initial distance parameter for an isotropic model, or variate with two values for a double-spherical isotropic model, or a variate with three values for a geometrical anisotropic model
ESTIMATES = variates Estimated parameter values
FITTEDVALUES = variates Fitted values
EXIT = scalars Exit status from the nonlinear fitting
SAVE = pointers Saves the model name and estimates in a pointer that can be used in KRIGE

Description

Procedure MVARIOGRAM uses the directives FIT, FITCURVE and FITNONLINEAR to fit various models to the experimental variogram. Models must be authorized in the sense that they cannot give rise to negative variances when data are combined. Technically they are conditionally negative semi-definite (CNSD); see Webster & Oliver (1990, 2001), Journel & Huijbregts (1978) or Chiles & Delfiner (1999) for an explanation.

The MODELTYPE option specifies the model that is to be fitted. There are bounded isotropic models with finite ranges. These all take the value c + c0 for h ≥ a, and the following values for h < a

    boundedlinear c0 + ch/a
    circular c0 + c {1 – (2/π)arccos(h/a)
  + (2h/(πa))√(1-h2/a2)}
    spherical c0 + c {1.5h/a – 0.5(h/a)3 }
    doublespherical c0 + c1 {1.5h/a1 – 0.5(h/a1)3 }
  + c2 {1.5h/a2 – 0.5(h/a2)3 }                                                            for ha1
  c0 + c1 + c2 {1.5h/a2 – 0.5(h/a2)3}                                            for a1 < h < a2
  where c = c1 + c2
    pentaspherical c0 + c {1.875h/a – 1.25(h/a)3 + 0.375(h/a)5}
    cubic c0 + c {7(h/a)2 – 8.75(h/a)3 + 3.5(h/a)5
  – 0.75(h/a)7}

There are also bounded asymptotic models

    besselk1 c0 + c {1 – h/a K1(h/a) }
  (Whittle’s elementary correlation, Whittle 1954)
    exponential c0 + c {1 – exp(-h/a)}
    gaussian c0 + c {1 – exp(-h2/a2)}
    stable c0 + c {1 – exp(-(h/a)b))}
    matern c0 + c {1 – 1 / (2(ν-1) Γ(ν)) (h/a)ν Kν(h/a)}

unbounded models

    power c0 + g hα
  (power function with exponent α strictly between 0 and 2)
    linear c0 + c h
  which is a special case of the power function with exponent 1

and hole effect models

    cardinalsine c0 + c × (1 – a/h × sin(h/a)) .

Geometrically anisotropic models, i.e. ones that might be made isotropic by a simple linear transformation of the spatial coordinates, can be fitted by setting option ISOTROPY=geometrical. The following transformation is used:

Ω(θ) = √{ a2cos2(θ-φ) + b2sin2(θ-φ) }

where θ represents the direction (specified by the DIRECTION parameter) converted from degrees to radians. So, for example, a geometrical anisotropic power model would be

c0 + ( Ω(θ) h )power

(Note: this particular model can also be defined by setting MODELTYPE=affinepower; the ISOTROPY option is then ignored.)

In all these models, the intercept term (or nugget variance) c0 can be omitted by setting the CONSTANT option to omit; the default is estimate.

For the stable model the SMOOTHNESS option controls the power parameter for the model, and for the matern model it specifies the ν parameter. By default, the parameter is estimated. However, you can supply a value, to fix the parameter for the model fitting.

The data for the procedure can be taken directly from the FVARIOGRAM directive, with parameters DISTANCES, VARIOGRAMS and COUNTS corresponding to those with the same names in FVARIOGRAM. The data will be in variates if the variogram was calculated in only one direction. If it is in several, they can either be in matrices (as generated by FVARIOGRAM) or in variates. For ISOTROPY=geometrical, directions must be supplied, using the DIRECTIONS parameter. These should be in a variate with one value for each column if the other data are in matrices; alternatively, they should be in a variate of the same length as the other variates.

The WEIGHTING option controls the weights that are used when fitting the model. The default setting counts uses the values supplied by the COUNTS parameter, cbyvar uses the COUNTS divided by the values in VARIOGRAM, and equal uses equal weights (of one).

By default, MVARIOGRAM generates rough starting values for the parameters. If the solution does not converge there are two likely reasons. The model may be unsuited for the particular experimental variogram. For example, a bounded model is specified when the variogram is clearly unbounded, or vice versa. You should choose only models that have approximately the right shape. Alternatively, the starting values may be too far from a sensible solution. You should then supply initial values using the INITIAL parameter. For a double-spherical isotropic model, INITIAL must be set to a variate with two values representing the two distance parameters. For the other isotropic models it should be set to a scalar defining the initial distance parameter. Finally, for a geometrical anisotropic model, it should be set to a variate with three values, defining the initial values for φ, the maximum distance parameter and the minimum distance parameter.

Printed output is controlled by the PRINT option, and includes all the usual settings as in FIT, FITCURVE or FITNONLINEAR. You can also produce a high-resolution graph of the data and the fitted model, by setting the WINDOW option to the number of a suitable window. By default WINDOW is zero, and no graph is produced. The TITLE option can supply a title for the plot. Option XUPPER can define an upper value for the x-axis (i.e. distance), and PENDATA and PENMODEL can supply the numbers of the pens to be used to plot the experimental variogram and the fitted model respectively (by default 1 and 2). Alternatively, you can use the ESTIMATES parameter to save the parameter estimates, and plot the variogram and model later with the DVARIOGRAM procedure.

The SAVE parameter saves the parameter estimates and associated information required by the KRIGE directive. The FITTEDVALUES parameter saves the fitted values, and the EXIT parameter saves the exit “status code” from FIT, FITCURVE or FITNONLINEAR. A zero value indicates success (see the Guide to the Genstat Command Language, Part 2, Section 3.7.4).

Options: PRINT, MODELTYPE, WEIGHTING, CONSTANT, SMOOTHNESS, ISOTROPY, WINDOW, TITLE, XUPPER, PENDATA, PENMODEL.

Parameters: VARIOGRAM, COUNTS, DISTANCE, DIRECTION, INITIAL, ESTIMATES, FITTEDVALUES, EXIT, SAVE.

Method

The model is fitted using directives FIT, FITCURVE or FITNONLINEAR as appropriate.

Action with RESTRICT

If the data variates are restricted, only the units not excluded by the restriction will be used.

References

Chiles, J-P. & Delfiner, P. (1999). Geostatistics: Modeling Spatial Uncertainty. Wiley, Chichester.

Journel, A.G. & Huijbregts, C.J. (1978). Mining Geostatistics. Academic Press, London.

Webster, R. & Oliver, M.A. (1990). Statistical Methods in Soil and Land Resource Survey. Oxford University Press.

Webster, R. & Oliver, M.A. (2001). Geostatistics for Environmental Scientists. Wiley, Chichester.

Whittle, P. (1954). On stationary processes in the plane. Biometrika, 41, 434-449.

See also

Directives: FVARIOGRAM, MCOVARIOGRAM, KRIGE.

Procedure: DVARIOGRAM.

Commands for: Spatial statistics.

Example

CAPTION 'MVARIOGRAM example',\
        !t('Model a variogram levels of potassium levels at',\
        'Brooms Barn Experimental Station (see Webster & Oliver,',\
        '1990, Statistical Methods in Soil and Land Resource Survey,',\
        'Oxford University Press, pages 267-269).'); STYLE=meta,plain
VARIATE [VALUES=0,45,90,135] Angles
MATRIX  [COLUMNS=13; ROWS=Angles] LogKvar,Kcounts,Midpoints
READ    [SERIAL=yes] LogKvar,Kcounts,Midpoints
        *        *        *        *
 0.005988 0.008046 0.006741 0.008364
 0.008062 0.009615 0.008969 0.010654
 0.011550 0.011688 0.010114 0.012356
 0.013897 0.013190 0.011354 0.012891
 0.015261 0.016131 0.014271 0.013980
 0.016509 0.017591 0.016121 0.016478
 0.018149 0.017546 0.016925 0.018412
 0.019388 0.018019 0.017641 0.021293
 0.020312 0.017763 0.016811 0.022991
 0.020957 0.017235 0.017137 0.025295
 0.020331 0.014874 0.017571 0.026301
 0.019182 0.014330 0.017281 0.028884  :
    0    0    0    0
  396  374  399  376
  362 1024  375 1032
  971  612 1014  620
  890  859  968  875
 1336 1294 1490 1336
 1227  939 1407  989
 1106 1489 1323 1606
 1340 1307 1684 1461
 1157 1144 1599 1326
 1239 1150 1929 1381
 1060  896 1875 1087
  879 1238 1808 1522  :
         *         *         *         *
  1.000000  1.414214  1.000000  1.414214
  2.000000  2.425227  2.000000  2.426056
  3.106801  3.605560  3.106274  3.605559
  4.081167  4.395480  4.080601  4.395571
  5.189774  5.451577  5.187020  5.453624
  6.159935  6.554852  6.157439  6.555201
  7.138072  7.469445  7.136054  7.472196
  8.236724  8.515127  8.231771  8.518978
  9.211962  9.551049  9.207607  9.550405
 10.309932 10.421715 10.303755 10.424751
 11.282172 11.394851 11.277505 11.395464
 12.259038 12.353009 12.254770 12.356227  :
CALCULATE Kcounts=Kcounts*(Midpoints<11.75)
" fit various possible models "
FOR Mod='LINEAR','SPHERICAL','EXPONENTIAL'
  MVARIOGRAM [MODELTYPE=#Mod; PRINT=model,summary,estimates;\
             WEIGHTING=counts]  LogKvar; COUNTS=Kcounts; DISTANCES=Midpoints
ENDFOR
Updated on March 7, 2019

Was this article helpful?