Kriging
This section describes the estimators used in the Kriging solver. Most users don't want to use estimators directly because they lack features such as neighborhood search and change of support. If that is your case, please refer to the solver documentation.
A Kriging estimator has the form:
\[\hat{Z}(\x_0) = \lambda_1 Z(\x_1) + \lambda_2 Z(\x_2) + \cdots + \lambda_n Z(\x_n),\quad \x_i \in \R^m, \lambda_i \in \R\]
with $Z\colon \R^m \times \Omega \to \R$ a random field.
This package implements the following Kriging variants:
- Simple Kriging
- Ordinary Kriging
- Universal Kriging
- External Drift Kriging
All these variants follow the same interface: an estimator object is first created with a given set of parameters, it is then combined with the data to obtain predictions at new locations.
The fit
function takes care of building the Kriging system and factorizing the LHS with an appropriate decomposition (e.g. Cholesky, LU):
StatsAPI.fit
— Methodfit(estimator, data)
Build Kriging system from data
and return a fitted estimator.
and the predict
function performs the estimation for a given variable and location:
GeoStatsBase.predict
— Methodpredict(estimator, var, uₒ)
Compute mean and variance of variable var
using the estimator
at point or geometry uₒ
.
Alternative constructors are provided for convenience that will immediately fit the Kriging parameters to the data. In this case, the data is passed as the first argument.
For advanced users, the Kriging weights and Lagrange multipliers at a given location can be accessed with the weights
method. This method returns a KrigingWeights
object containing a field λ
for the weights and a field ν
for the Lagrange multipliers:
KrigingEstimators.weights
— Functionweights(estimator, uₒ)
Compute the weights λ (and Lagrange multipliers ν) for the estimator
at point or geometry uₒ
.
These weights can be combined with a vector of observations using the combine
function:
KrigingEstimators.combine
— Functioncombine(estimator, weights, z)
Combine weights
with values z
to produce mean and variance using the appropriate formulas for the estimator
.
This functionality can be useful in real-time applications when the locations of the observations are fixed and an online stream of data is available.
Finally, all estimators work with general Hilbert spaces, meaning that one can interpolate any data type that implements scalar multiplication, vector addition and inner product.
Simple Kriging
In Simple Kriging, the mean $\mu$ of the random field is assumed to be constant and known. The resulting linear system is:
\[\begin{bmatrix} cov(\x_1,\x_1) & cov(\x_1,\x_2) & \cdots & cov(\x_1,\x_n) \\ cov(\x_2,\x_1) & cov(\x_2,\x_2) & \cdots & cov(\x_2,\x_n) \\ \vdots & \vdots & \ddots & \vdots \\ cov(\x_n,\x_1) & cov(\x_n,\x_2) & \cdots & cov(\x_n,\x_n) \end{bmatrix} \begin{bmatrix} \lambda_1 \\ \lambda_2 \\ \vdots \\ \lambda_n \end{bmatrix} = \begin{bmatrix} cov(\x_1,\x_0) \\ cov(\x_2,\x_0) \\ \vdots \\ cov(\x_n,\x_0) \end{bmatrix}\]
or in matricial form $\C\l = \c$. We subtract the given mean from the observations $\boldsymbol{y} = \z - \mu \1$ and compute the mean and variance at location $\x_0$:
\[\mu(\x_0) = \mu + \boldsymbol{y}^\top \l\]
\[\sigma^2(\x_0) = cov(0) - \c^\top \l\]
KrigingEstimators.SimpleKriging
— TypeSimpleKriging(γ, μ)
SimpleKriging(data, γ, μ)
Simple Kriging with variogram model γ
and constant mean μ
.
Optionally, pass the geospatial data
to the fit
function.
Notes
- Simple Kriging requires stationary variograms
Ordinary Kriging
In Ordinary Kriging the mean of the random field is assumed to be constant and unknown. The resulting linear system is:
\[\begin{bmatrix} \G & \1 \\ \1^\top & 0 \end{bmatrix} \begin{bmatrix} \l \\ \nu \end{bmatrix} = \begin{bmatrix} \g \\ 1 \end{bmatrix}\]
with $\nu$ the Lagrange multiplier associated with the constraint $\1^\top \l = 1$. The mean and variance at location $\x_0$ are given by:
\[\mu(\x_0) = \z^\top \l\]
\[\sigma^2(\x_0) = \begin{bmatrix} \g \\ 1 \end{bmatrix}^\top \begin{bmatrix} \l \\ \nu \end{bmatrix}\]
KrigingEstimators.OrdinaryKriging
— TypeOrdinaryKriging(γ)
OrdinaryKriging(data, γ)
Ordinary Kriging with variogram model γ
.
Optionally, pass the geospatial data
to the fit
function.
Universal Kriging
In Universal Kriging, the mean of the random field is assumed to be a polynomial of the spatial coordinates:
\[\mu(\x) = \sum_{k=1}^{N_d} \beta_k f_k(\x)\]
with $N_d$ monomials $f_k$ of degree up to $d$. For example, in 2D there are $6$ monomials of degree up to $2$:
\[\mu(x_1,x_2) = \beta_1 1 + \beta_2 x_1 + \beta_3 x_2 + \beta_4 x_1 x_2 + \beta_5 x_1^2 + \beta_6 x_2^2\]
The choice of the degree $d$ determines the size of the polynomial matrix
\[\F = \begin{bmatrix} f_1(\x_1) & f_2(\x_1) & \cdots & f_{N_d}(\x_1) \\ f_1(\x_2) & f_2(\x_2) & \cdots & f_{N_d}(\x_2) \\ \vdots & \vdots & \ddots & \vdots \\ f_1(\x_n) & f_2(\x_n) & \cdots & f_{N_d}(\x_n) \end{bmatrix}\]
and polynomial vector $\f = \begin{bmatrix} f_1(\x_0) & f_2(\x_0) & \cdots & f_{N_d}(\x_0) \end{bmatrix}^\top$.
The variogram determines the variogram matrix:
\[\G = \begin{bmatrix} \gamma(\x_1,\x_1) & \gamma(\x_1,\x_2) & \cdots & \gamma(\x_1,\x_n) \\ \gamma(\x_2,\x_1) & \gamma(\x_2,\x_2) & \cdots & \gamma(\x_2,\x_n) \\ \vdots & \vdots & \ddots & \vdots \\ \gamma(\x_n,\x_1) & \gamma(\x_n,\x_2) & \cdots & \gamma(\x_n,\x_n) \end{bmatrix}\]
and the variogram vector $\g = \begin{bmatrix} \gamma(\x_1,\x_0) & \gamma(\x_2,\x_0) & \cdots & \gamma(\x_n,\x_0) \end{bmatrix}^\top$.
The resulting linear system is:
\[\begin{bmatrix} \G & \F \\ \F^\top & \boldsymbol{0} \end{bmatrix} \begin{bmatrix} \l \\ \boldsymbol{\nu} \end{bmatrix} = \begin{bmatrix} \g \\ \f \end{bmatrix}\]
with $\boldsymbol{\nu}$ the Lagrange multipliers associated with the universal constraints. The mean and variance at location $\x_0$ are given by:
\[\mu(\x_0) = \z^\top \l\]
\[\sigma^2(\x_0) = \begin{bmatrix}\g \\ \f\end{bmatrix}^\top \begin{bmatrix}\l \\ \boldsymbol{\nu}\end{bmatrix}\]
KrigingEstimators.UniversalKriging
— TypeUniversalKriging(γ, degree, dim)
UniversalKriging(data, γ, degree)
Universal Kriging with variogram model γ
and polynomial degree
on a spatial domain of dimension dim
.
Optionally, pass the geospatial data
to the fit
function.
Notes
OrdinaryKriging
is recovered for 0th degree polynomial- For non-polynomial mean, see
ExternalDriftKriging
External Drift Kriging
In External Drift Kriging, the mean of the random field is assumed to be a combination of known smooth functions:
\[\mu(\x) = \sum_k \beta_k m_k(\x)\]
Differently than Universal Kriging, the functions $m_k$ are not necessarily polynomials of the spatial coordinates. In practice, they represent a list of variables that is strongly correlated (and co-located) with the variable being estimated.
External drifts are known to cause numerical instability. Give preference to other Kriging variants if possible.
KrigingEstimators.ExternalDriftKriging
— TypeExternalDriftKriging(γ, drifts)
ExternalDriftKriging(data, γ, drifts)
External Drift Kriging with variogram model γ
and external drifts
functions.
Optionally, pass the geospatial data
to the fit
function.
Notes
- External drift functions should be smooth
- Kriging system with external drift is often unstable
- Include a constant drift (e.g.
x->1
) for unbiased estimation OrdinaryKriging
is recovered fordrifts = [x->1]
- For polynomial mean, see
UniversalKriging