Kriging_Puntual

Kriging

Kriging es una técnica de interpolación espacial que permite estimar valores desconocidos en una superficie a partir de los valores conocidos en ubicaciones cercanas. Es especialmente útil para datos de alta densidad y para datos con estructuras espaciales complejas. En este cuaderno, vamos a explorar los diferentes tipos de Kriging y cómo implementarlos en R.

Antes de comenzar, asegúrate de tener instalados los siguientes paquetes de R: “gstat”, “sp”, “raster” y “ggplot2”. Si no los tienes, puedes instalarlos mediante el siguiente código:

# Instalar paquetes si es necesario

library(gstat)
library(sp)
library(raster)
library(ggplot2)
library(phylin)

Adjuntando el paquete: 'phylin'

The following object is masked from 'package:gstat':

    idw

# Generación de datos simulados para los ejemplos
set.seed(1029)
coords <- data.frame(x = runif(10, 0, 100), y = runif(10, 0, 100))
values <- data.frame(value = rnorm(10))

# Combinar 'coords' y 'values' en un solo data frame
data <- cbind(coords, values)

# Convertir 'data' en un objeto espacial
coordinates(data) <- ~ x + y
proj4string(data) <- CRS("+proj=utm +zone=33 +datum=WGS84")

Introducción Kriging

Kriging es una técnica de interpolación espacial que se basa en la teoría de la estadística espacial. La idea fundamental detrás de Kriging es que la estimación del valor desconocido en una ubicación dada se hace como una combinación lineal ponderada de los valores conocidos en las ubicaciones vecinas. Los pesos de esta combinación se determinan utilizando la información de covarianza entre los puntos en el espacio. Esta técnica se utiliza comúnmente en estadística espacial, geología, minería, geofísica, oceanografía y otros campos que involucran la toma de medidas en lugares dispersos en el espacio.

El nombre “Kriging” proviene del geólogo Danie G. Krige, quien fue el primero en aplicar esta técnica para la exploración de minerales en Sudáfrica en la década de 1950

Kriging Simple

En esta sección se explica el método de kriging simple, que es una técnica de interpolación espacial que utiliza una función de tendencia conocida y un modelo de variograma ajustado a los datos. El objetivo es estimar el valor de una variable Z en una ubicación s0 a partir de los valores observados en n ubicaciones si cercanas. El kriging simple asume que la variable Z se puede descomponer como:

\[Z(s) = μ(s) + ε(s),\]

donde \(μ(s)\) es la función de tendencia conocida y \(ε(s)\) es un proceso espacial aleatorio con media cero y covarianza \(C(s_i,s_j) = Cov(ε(s_i),ε(s_j))\). El predictor del kriging simple es una combinación lineal de los valores observados:

\[p(Z,s_0) = μ(s_0) + ∑_{i=1}^n λ_i (Z(s_i) - μ(s_i)),\]

donde los pesos \(λ_i\) se obtienen resolviendo el sistema de ecuaciones:

\[∑_{j=1}^n λ_j C(s_i,s_j) - C(s_i,s_0) = 0, \ i = 1,...,n.\]

El kriging simple minimiza el error cuadrático medio entre el valor estimado y el verdadero, que se denomina varianza kriging y se calcula como:

\[σ^2 KS(s_0) = C(s_0,s_0) - ∑_{i=1}^n ∑_{j=1}^n λ_i λ_j C(s_i,s_j).\]

Para implementar el kriging simple en R se puede utilizar la función krig del paquete phylin, que requiere como argumentos los valores observados, las coordenadas de las ubicaciones muestreadas y las ubicaciones a interpolar, el modelo de variograma ajustado y el valor conocido de la media (o NA para usar kriging ordinario). A continuación se muestra un ejemplo con datos simulados:

# Calcular el variograma experimental
vgram <- variogram(value ~ 1, data)

# Ajustar un modelo de variograma
model <- vgm(psill = 1, model = "Sph", range = 50, nugget = 0)
vfit <- fit.variogram(vgram, model)

Warning in fit.variogram(vgram, model): linear model has singular covariance
matrix

Warning in fit.variogram(vgram, model): singular model in variogram fit

Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
fit.ranges, : No convergence after 200 iterations: try different initial
values?

# Crear una cuadrícula de predicción
grd <- expand.grid(x = seq(0, 100, by = 5), y = seq(0, 100, by = 5))
coordinates(grd) <- ~ x + y
gridded(grd) <- TRUE
proj4string(grd) <- CRS("+proj=utm +zone=33 +datum=WGS84")

# Kriging Simple
kriged_values <- krige(value ~ 1, data, grd, model = vfit, nmax = 5)

[using ordinary kriging]

# Visualización
spplot(kriged_values["var1.pred"], main = "Kriging Simple Prediction")

Kriging Ordinario

El kriging ordinario es el tipo más general y más utilizado de kriging. Presupone que el valor medio constante es desconocido y lo estima a partir de los datos disponibles. Esta es una suposición razonable a menos que haya una razón científica para rechazarla.

El kriging ordinario se puede expresar como:

\[Z^*(x_0) = \sum_{i=1}^n \lambda_i Z(x_i)\]

Donde \(Z^*(x_0)\) es el valor estimado en el punto \(x_0\), \(Z(x_i)\) son los valores medidos en los puntos \(x_i\), \(n\) es el número de puntos medidos y \(\lambda_i\) son los pesos óptimos que minimizan la varianza del error de estimación.

Los pesos óptimos se obtienen resolviendo un sistema de ecuaciones lineales conocido como sistema krigiano:

\[ \begin{bmatrix} \gamma_{11} & \gamma_{12} & \cdots & \gamma_{1n} & 1 \\ \gamma_{21} & \gamma_{22} & \cdots & \gamma_{2n} & 1 \\ \vdots & \vdots & \ddots & \vdots & 1 \\ \gamma_{n1} & \gamma_{n2} & \cdots & \gamma_{nn} & 1 \\ 1 & 1 & \cdots & 1 & 0 \end{bmatrix} \begin{bmatrix} \lambda_1 \\ \lambda_2 \\ \vdots \\ \lambda_n \\ m \end{bmatrix} = \begin{bmatrix} \gamma(x_0,x_1) \\ \gamma(x_0,x_2) \\ \vdots \\ \gamma(x_0,x_n) \\ -1 \end{bmatrix} \]

Donde \(\gamma_{ij}\) son los semivariogramas entre los puntos medidos, \(\gamma(x_0,x_i)\) son los semivariogramas entre el punto a estimar y los puntos medidos y \(m\) es un multiplicador lagrangeano que representa la estimación del valor medio constante.

Para implementar el kriging ordinario en R, se pueden utilizar las librerías “gstat” o “kriging”. A continuación se muestra un ejemplo usando la librería “gstat”.

# Ajustar modelo de variograma para kriging ordinario
vgram <- variogram(value ~ 1, data)
vfit <- fit.variogram(vgram, model)

Warning in fit.variogram(vgram, model): linear model has singular covariance
matrix

Warning in fit.variogram(vgram, model): singular model in variogram fit

Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
fit.ranges, : No convergence after 200 iterations: try different initial
values?

# Kriging ordinario
kriged_values_ordinary <- krige(value ~ 1, data, grd, model = vfit, nmax = 5)

[using ordinary kriging]

# Visualización
spplot(kriged_values_ordinary["var1.pred"], main = "Kriging Ordinario Prediction")

Kriging Universal

El kriging universal es una variante del kriging que permite incorporar variables auxiliares que influyen en la variable de interés. Por ejemplo, si queremos interpolar la temperatura en una región, podemos usar como variables auxiliares la altitud, la latitud o la distancia al mar. El kriging universal asume que la variable de interés se puede expresar como una combinación lineal de las variables auxiliares más un error espacialmente correlacionado. Es decir:

\[ Z(s) = p ∑_{j=0} X_j(s)β_j +ε(s), \]

donde \(Z(s)\) es la variable de interés en el punto \(s\), \(X_j(s)\) son las variables auxiliares (incluyendo una constante), βj son los coeficientes desconocidos y ε(s) es el error espacial con variograma conocido.

Para estimar \(Z(s_0)\) en un punto no observado \(s_0\), se usa el predictor del kriging universal:

\[ p(Z,s_0) = ∑_{i=1}^n λ_iZ(s_i), \]

donde \(λ_i\) son los pesos que minimizan el error cuadrático medio de predicción sujeto a las restricciones:

\[λ^⊤X = x^⊤ 0,\]

donde \(X\) es una matriz con las variables auxiliares en los puntos observados y \(x_0\) es un vector con las variables auxiliares en el punto \(s_0\).

Para implementar el kriging universal en R, podemos usar el paquete gstat y la función krige. Esta función requiere una fórmula que defina la variable dependiente como un modelo lineal de las variables independientes, un objeto espacial con las coordenadas y los datos observados, un objeto espacial con las coordenadas de los puntos a predecir y un modelo de variograma para el error espacial.

A continuación se muestra un ejemplo de cómo usar krige para realizar kriging universal con dos variables auxiliares: \(x\) e \(y\).

# Añadir variables auxiliares al conjunto de datos
values$x <- coordinates(data)[,1]
values$y <- coordinates(data)[,2]

# Ajustar variograma
vgram_universal <- variogram(value ~ x + y, data)
vfit_universal <- fit.variogram(vgram_universal, vgm(model = "Sph", psill = 1, range = 50, nugget = 0))

Warning in fit.variogram(vgram_universal, vgm(model = "Sph", psill = 1, :
singular model in variogram fit

Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
fit.ranges, : No convergence after 200 iterations: try different initial
values?

# Kriging Universal
kriged_values_universal <- krige(value ~ x + y, data, grd, model = vfit_universal)

[using universal kriging]

# Visualización
spplot(kriged_values_universal["var1.pred"], main = "Kriging Universal Prediction")