Theoretical variograms

Models

In an intrinsic isotropic model, the variogram is only a function of the distance between any two points $\x_1,\x_2 \in \R^m$:

\[\gamma(\x_1,\x_2) = \gamma(||\x_1 - \x_2||) = \gamma(h)\]

Under the additional assumption of 2nd-order stationarity, the well-known covariance is directly related via $\gamma(h) = \cov(0) - \cov(h)$. This package implements a few commonly used as well as other more exotic variogram models. Most of these models share a set of default parameters (e.g. sill=1.0, range=1.0), which can be set with keyword arguments.

Functions are provided to query properties of variogram models programatically:

Meshes.isisotropicMethod
isisotropic(γ)

Tells whether or not variogram γ is isotropic.

Gaussian

\[\gamma(h) = (s - n) \left[1 - \exp\left(-3\left(\frac{h}{r}\right)^2\right)\right] + n \cdot \1_{(0,\infty)}(h)\]

Variography.GaussianVariogramType
GaussianVariogram(range=r, sill=s, nugget=n)
GaussianVariogram(ball; sill=s, nugget=n)

A Gaussian variogram with range r, sill s and nugget n. Optionally, use a custom metric ball.

plot(GaussianVariogram())

Exponential

\[\gamma(h) = (s - n) \left[1 - \exp\left(-3\left(\frac{h}{r}\right)\right)\right] + n \cdot \1_{(0,\infty)}(h) \]

Variography.ExponentialVariogramType
ExponentialVariogram(range=r, sill=s, nugget=n)
ExponentialVariogram(ball; sill=s, nugget=n)

An exponential variogram with range r, sill s and nugget n. Optionally, use a custom metric ball.

plot(ExponentialVariogram())

Matern

\[\gamma(h) = (s - n) \left[1 - \frac{2^{1-\nu}}{\Gamma(\nu)} \left(\sqrt{2\nu}\frac{h}{r}\right)^\nu K_\nu\left(\sqrt{2\nu}\frac{h}{r}\right)\right] + n \cdot \1_{(0,\infty)}(h)\]

Variography.MaternVariogramType
MaternVariogram(range=r, sill=s, nugget=n, order=ν)
MaternVariogram(ball; sill=s, nugget=n, order=ν)

A Matérn variogram with range r, sill s and nugget n. The parameter ν is the order of the Bessel function. Optionally, use a custom metric ball.

plot(MaternVariogram())

Spherical

\[\gamma(h) = (s - n) \left[\left(\frac{3}{2}\left(\frac{h}{r}\right) + \frac{1}{2}\left(\frac{h}{r}\right)^3\right) \cdot \1_{(0,r)}(h) + \1_{[r,\infty)}(h)\right] + n \cdot \1_{(0,\infty)}(h)\]

Variography.SphericalVariogramType
SphericalVariogram(range=r, sill=s, nugget=n)
SphericalVariogram(ball; sill=s, nugget=n)

A spherical variogram with range r, sill s and nugget n. Optionally, use a custom metric ball.

plot(SphericalVariogram())

Cubic

\[\gamma(h) = (s - n) \left[\left(7\left(\frac{h}{r}\right)^2 - \frac{35}{4}\left(\frac{h}{r}\right)^3 + \frac{7}{2}\left(\frac{h}{r}\right)^5 - \frac{3}{4}\left(\frac{h}{r}\right)^7\right) \cdot \1_{(0,r)}(h) + \1_{[r,\infty)}(h)\right] + n \cdot \1_{(0,\infty)}(h)\]

Variography.CubicVariogramType
CubicVariogram(range=r, sill=s, nugget=n)
CubicVariogram(ball; sill=s, nugget=n)

A cubic variogram with range r, sill s and nugget n. Optionally, use a custom metric ball.

plot(CubicVariogram())

Pentaspherical

\[\gamma(h) = (s - n) \left[\left(\frac{15}{8}\left(\frac{h}{r}\right) - \frac{5}{4}\left(\frac{h}{r}\right)^3 + \frac{3}{8}\left(\frac{h}{r}\right)^5\right) \cdot \1_{(0,r)}(h) + \1_{[r,\infty)}(h)\right] + n \cdot \1_{(0,\infty)}(h)\]

Variography.PentasphericalVariogramType
PentasphericalVariogram(range=r, sill=s, nugget=n)
PentasphericalVariogram(ball; sill=s, nugget=n)

A pentaspherical variogram with range r, sill s and nugget n. Optionally, use a custom metric ball.

plot(PentasphericalVariogram())

Power

\[\gamma(h) = sh^a + n \cdot \1_{(0,\infty)}(h)\]

Variography.PowerVariogramType
PowerVariogram(scaling=s, exponent=a, nugget=n)

A power variogram with scaling s, exponent a and nugget n.

plot(PowerVariogram())

Sine hole

\[\gamma(h) = (s - n) \left[1 - \frac{\sin(\pi h / r)}{\pi h / r}\right] + n \cdot \1_{(0,\infty)}(h)\]

Variography.SineHoleVariogramType
SineHoleVariogram(range=r, sill=s, nugget=n)
SineHoleVariogram(ball; sill=s, nugget=n)

A sine hole variogram with range r, sill s and nugget n. Optionally, use a custom metric ball.

plot(SineHoleVariogram())

Nugget

\[\gamma(h) = n \cdot \1_{(0,\infty)}(h)\]

plot(NuggetEffect(1.0))

Anisotropy

Anisotropic models are easily obtained by defining an ellipsoid metric in place of the default Euclidean metric as shown in the following example. First, we create a custom metric ball, which specifies the ranges and angles of rotation of the ellipsoid:

ellipsoid = MetricBall((3.0, 2.0, 1.0), EulerAngles(0.0, 0.0, 0.0))
MetricBall((3.0, 2.0, 1.0), Mahalanobis)

We then pass this ball as the first argument to the variogram model instead of specifying a single range with a keyword argument:

GaussianVariogram(ellipsoid, sill=2.0)
GaussianVariogram (anisotropic)
  └─sill ⇨ 2.0
  └─nugget ⇨ 0.0
  └─ranges ⇨ (3.0, 2.0, 1.0)
  └─metric ⇨ Mahalanobis

To illustrate the concept, consider the following 2D data set:

dim, nobs = 2, 50
table = (Z=rand(nobs),)
coord = 100rand(dim, nobs)

𝒮 = georef(table, coord)

plot(𝒮)

and the corresponding estimation problem on a Cartesian grid:

problem = EstimationProblem(𝒮, CartesianGrid(100,100), :Z)
2D EstimationProblem
  data:      50 MeshData{2,Float64}
  domain:    100×100 CartesianGrid{2,Float64}
  variables: Z (Float64)

We solve the problem with different ellipsoids by varying the angle of rotation from $0$ to $2\pi$ clockwise:

anim = @animate for θ in range(0, stop=2π, length=10)
  # ellipsoid rotated clockwise by angle θ
  e = MetricBall((20.,5.), ClockwiseAngle(θ))

  # anisotropic variogram model
  γ = GaussianVariogram(e)

  # solve the problem with Kriging
  sol = solve(problem, Kriging(:Z => (variogram=γ,)))

  # plot current frame
  plot(sol)
end

# generate gif from list of frames
gif(anim, "anisotropy.gif", fps=1)