Understanding Variogram in Geospatial Analysis

A Note on Spatial Relationships

The First Law of Geography, proposed by Waldo Tobler, states:

“Everything is related to everything else, but near things are more related than distant things.”

This principle forms the backbone of variogram analysis and Kriging.

Key Components of a Variogram

1. Lag Distance (h): The separation distance between data points.
2. Semivariance (γ(h)): A measure of the dissimilarity between data points separated by distance h.
3. Sill: The semivariance value where the variogram levels off, representing the data’s overall variance.
4. Range: The distance at which the variogram reaches the sill, indicating the limit of spatial correlation.
5. Nugget: The y-intercept of the variogram, accounting for measurement errors or microscale variations.

Variogram Models

Variogram models are mathematical functions used to fit experimental variograms. These models capture the spatial dependence of the data and are essential for Kriging. Commonly used models include:

1. Spherical Model:
   • Formula: $$\gamma(h) = \begin{cases} C_0 + C \left[\frac{3h}{2a} - \frac{h^3}{2a^3}\right], & 0 \leq h \leq a \\ C_0 + C, & h > a \end{cases}$$
   • Characteristics:
     • Commonly used in geosciences.
     • Levels off at the range (a).
     • Represents a gradual transition of spatial correlation.
2. Exponential Model:
   • Formula: $$\gamma(h) = C_0 + C \left(1 - e^{-h/a}\right)$$
   • Characteristics:
     • Suitable for phenomena with continuous but less abrupt changes.
     • Correlation decreases exponentially with distance.
3. Gaussian Model:
   • Formula: $$\gamma(h) = C_0 + C \left(1 - e^{-\left(\frac{h}{a}\right)^2}\right)$$
   • Characteristics:
     • Smooth transitions, often used in environmental studies.
     • Effective for data with strong local correlations.
4. Linear Model:
   • Formula: $$\gamma(h) = C_0 + C h$$
   • Characteristics:
     • Represents a direct, proportional increase in semivariance with distance.
     • Simplistic but useful for certain datasets.
5. Power Model:
   • Formula: $$\gamma(h) = C_0 + C h^k, \quad 0 < k < 2$$
   • Characteristics:
     • Flexible model for datasets with non-linear spatial relationships.

Each model serves specific purposes based on the dataset and the spatial behavior being analyzed. Selecting the appropriate variogram model is crucial for accurate Kriging.

The summary of selecting a variogram model presented below:

Experimental Variogram

if (!require(sp)) { install.packages('sp'); library(sp)}
if (!require(gstat)) { install.packages('gstat'); library(gstat)}

# Convert to spatial object
coordinates(location) <- ~x + y

# Compute experimental variogram
variogram_exp <- variogram(water_level ~ 1, location)
head(variogram_exp)copy

##   np      dist     gamma dir.hor dir.ver   id
## 1  2  7.912548 12.129825       0       0 var1
## 2  2 11.366795  9.534399       0       0 var1
## 3  2 14.046246 17.299652       0       0 var1
## 4  6 15.922439 23.613955       0       0 var1
## 5  4 19.786952 31.321265       0       0 var1
## 6  5 22.195585 18.261615       0       0 var1copy

The formula used above have format independent_variable ~ 1 or independent_variable ~ x+y. X or/and Y, any independent variables or covariates (if applicable). If the formula is ~ 1, it means there are no covariates, and the variogram will be calculated solely based on spatial relationships. The output of variogram_exp is has field such as np, dist, gamma etc…

• np: the number of point pairs for this estimate; in case of a variogramCloud see below
• dist: the average distance of all point pairs considered for this estimate
• gamma: the actual sample variogram estimate
• dir.hor: the horizontal direction
• dir.ver: the vertical direction
• id: the combined id pair

for more information, read the help of gstat::variogram help page in R.

Semivariance Calculation

Experimental variogram will be created with Semivariance against distance. Process is explained is, first, the semivariance $\gamma(h)$ is computed for each distance bin $h$. The formula is:

$$\gamma(h) = \frac{1}{2N(h)} \sum_{i=1}^{N(h)} \left( z(x_i) - z(x_i + h) \right)^2$$

Where:

• $N(h)$: Number of point pairs in the distance bin $h$.
• $z(x_i)$: Observed value at location $x_i$.
• $z(x_i + h)$: Observed value at location $x_i + h$.

This gives the average squared difference of observed values separated by $h$.

The formula shown is for computing semivariance for a given lag or distance ℎ, but it is not directly tied to a specific variogram model such as spherical, exponential, or Gaussian. This formula is a general expression used to calculate experimental semivariance, which is then plotted as an experimental variogram Once the experimental variogram is computed, you can fit it with a specific theoretical variogram model, such as Spherical, Exponential, Gaussian, etc. To identify the specific model, you would calculate the semivariances using the given formula, plot them against distance h, and then choose the model that best fits the curve.

There is a very nice youtube video by YT@MiningGeologist explaing the computation semivariance using the above equation.

plot(variogram_exp, main = "Experimental Variogram")copy

Fitting a Theoretical Model

This involves selecting a mathematical function (e.g., spherical, exponential, Gaussian) to approximate the experimental variogram. The goal is to create a theoretical variogram that can be used in kriging or other geostatistical methods. The output is a smooth, continuous variogram model that represents the spatial variability of the data.

# Fit a variogram model
variogram_model <- fit.variogram(variogram_exp, model = vgm("Sph", range = 50, nugget = 2))

# Plot the fitted model
plot(variogram_exp, variogram_model, main = "Fitted Variogram Model (Spherical)")copy

Trying Another Model

# Fit a variogram model
variogram_model_exp <- fit.variogram(variogram_exp, model = vgm("Exp", nugget = 2))

# Plot the fitted model
plot(variogram_exp, variogram_model_exp, main = "Fitted Variogram Model (Exponential)")copy