Caso 1. Estudio de Variabilidad Espacial de Variables Climáticas

Se recopilaron datos climaticos en una finca en cierta zona de Colombia. En esta finca se tomaron 394 registros de diversas variables climáticas como la Humedad Relativa y la Temperatura (“datos_geo.RData”).

Se realiza un estudio de la variabilidad espacial de estas dos caracteristicas (Humedad Relativa y Temperatura) por medio de la metodología de Geoestadística

Carga de datos y paquetes

pacman::p_load(sf,mapview,spData,ggplot2,terra,sp,sf,malariaAtlas,spatstat,
               sparr,geoFourierFDA,rnaturalearth,epiflows,stars,ncdf4,
               rnaturalearth,geodata,raster,geoR,patchwork,viridis,leaflet,tmap,
               htmltools,leafsync,htmlwidgets,webshot,RColorBrewer,plotly,
               flowmapblue,chirps,elevatr,osmdata,wbstats,spocc,rdhs,gstat,
               malariaAtlas,openair,ggmap,opencage,spdep,gridExtra,dplyr)



load("caso 1/datos_geo.RData")
zona <- st_read("caso 1/zona/area.shp")

## Reading layer `area' from data source 
##   `C:\Users\PC\Desktop\Univalle\8VO SEMESTRE\ANÁLISIS DE DATOS ESPACIALES\caso 1\zona\area.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 1 feature and 1 field
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -76.61417 ymin: 2.380766 xmax: -76.61251 ymax: 2.382726
## Geodetic CRS:  WGS 84

Análisis Exploratorio de las Variables

Estadísticas descriptivas

medias <- apply(datos_geo[, 1:2], 2, mean)
desviaciones <- apply(datos_geo[, 1:2], 2, sd)
coef_var <- apply(datos_geo[, 1:2], 2, function(x) sd(x) / mean(x))
rangos <- apply(datos_geo[, 1:2], 2, range)

descr <- data.frame(
  min = rangos[1, ],
  media = medias,
  desvs = desviaciones,
  coef_var = coef_var,
  max = rangos[2, ]
)
print(descr)

##                    min    media     desvs  coef_var  max
## Temperature       16.0 22.62030  4.288311 0.1895780 34.6
## Relative_Humidity 28.5 58.08503 14.014297 0.2412721 91.4

En el análisis exploratorio de las variables, se observa que la Temperatura presenta un rango de valores entre 16.0 °C y 33.0 °C, con una media de 22.62 °C y una desviación estándar de 4.29 °C. El coeficiente de variación es del 18.96%, lo que indica una variación moderada en los datos, sugiriendo que la temperatura es relativamente consistente en la región estudiada.

Por otro lado, la Humedad Relativa muestra un rango más amplio, oscilando entre el 28.5% y el 91.4%, con una media de 58.09% y una desviación estándar de 14.01%. Su coeficiente de variación es del 24.17%, reflejando una mayor dispersión en comparación con la temperatura. Esto sugiere que la humedad relativa tiene fluctuaciones más significativas, posiblemente influenciadas por factores locales como la vegetación o la proximidad a fuentes de agua. Los resultados descriptivos reflejan una mayor estabilidad en la variable de temperatura frente a una mayor variabilidad en la humedad relativa.

Histogramas

h1 <- ggplot(datos_geo, aes(x = Temperature)) +
  geom_histogram(binwidth = 1, fill = "red", color = "black", alpha = 0.7) +
  labs(title = "", x = "Temperatura (°C)", y = "Frecuencia") +
  theme_minimal()

h2 <- ggplot(datos_geo, aes(x = Relative_Humidity)) +
  geom_histogram(binwidth = 5, fill = "skyblue", color = "black", alpha = 0.7) +
  labs(title = "", x = "Humedad Relativa (%)", y = "Frecuencia") +
  theme_minimal()

grid.arrange(h1, h2, nrow = 1)

Figura 1: Histograma: Temperatura y Humedad Relativa

Como se aprecia en la Figura 1, el histograma de la Temperatura (color rojo) muestra un único pico significativo alrededor de los 20 °C, lo que indica que esta es la temperatura más frecuente en el conjunto de datos. La frecuencia disminuye consistentemente a medida que los valores se acercan a los 25 °C, con un ligero incremento luego de esta temperatura, sin embargo decrece muy rapido al llegar a los 30°C. Por lo tanto, se evidencia una ligera bimodalidad, además de cierta sesgadez hacia los valores más bajos de temperatura.

El histograma de la Humedad Relativa (color azul claro) presenta un pico principal entre 60% y 70% aproximadamente, que corresponde a la mayor concentración de datos en esta región. También hay un aumento muy notable de frecuencia significativo en el 40%, pero este no alcanza la misma magnitud del pico principal. Esto sugiere una clara bimodalidad, los valores más altos de humedad relativa(%) dominan claramente la distribución.

Boxplots

bp1 <- ggplot(datos_geo, aes(y = Temperature)) +
  geom_boxplot(fill = "red", color = "black", alpha = 0.7) +
  labs(title = "", y = "Temperatura (°C)") +
  theme_minimal()

bp2 <- ggplot(datos_geo, aes(y = Relative_Humidity)) +
  geom_boxplot(fill = "skyblue", color = "black", alpha = 0.7) +
  labs(title = "", y = "Humedad Relativa (%)") +
  theme_minimal()

grid.arrange(bp1, bp2, ncol = 2)

Figura 2: Boxplot: Temperatura - Humedad Relativa

Como se observa en la Figura 2, El diagrama de caja de la Temperatura refleja que la mayor parte de los datos se encuentran concentrados entre 20 °C y 25 °C. La mediana, ubicada cerca de los 22 °C, muestra que los datos están distribuidos ligeramente hacia valores más bajos. Los bigotes indican una variación moderada, abarcando el rango general de las temperaturas, sin presencia de valores atípicos destacados.

En el caso de la Humedad Relativa, la mayor concentración de valores está entre el 50% y el 65%, con una mediana cercana al 58%. Esto indica una distribución más amplia en comparación con la temperatura.os bigotes indican una variación moderada, abarcando el rango general de los porcentajes de humedad relativa, es decir, tampoco se observan valores atípicos.

La Temperatura tiene una distribución más compacta y centrada, mientras que la Humedad Relativa presenta una mayor dispersión, esto refuerza lo observado en los respectivos histogramas.

Relación entre Humedad Relativa y Temperatura

gdis <- ggplot(datos_geo, aes(x = Temperature, y = Relative_Humidity)) +
  geom_point(alpha = 0.7, color = "black", size = 2) +
  geom_smooth(method = "lm", color = "red", linetype = "dashed") +
  labs(
    title = "",
    x = "Temperatura (°C)",
    y = "Humedad Relativa (%)"
  ) +
  theme_minimal()

gdis

Figura 3: Relacion entre Humedad relativa y temperatura

La Figura 3 evidencia una relación negativa entre la Humedad Relativa y la Temperatura, donde a medida que la temperatura aumenta, los niveles de humedad relativa tienden a disminuir. Esta tendencia es claramente capturada por la línea de regresión ajustada, lo que sugiere una correlación inversa consistente entre ambas variables. Aunque los datos muestran cierta dispersión alrededor de la línea, especialmente en los extremos de la temperatura, el patrón general refleja cómo temperaturas más altas están asociadas con menores niveles de humedad, una dinámica común en condiciones meteorológicas.

Análisis Exploratorio Espacial

g1 = ggplot(zona) + geom_sf() +
  geom_point(data = datos_geo,aes(x = Longitude, y = Latitude, size = Temperature,
                                  color = Temperature)) + 
  theme_minimal() + scale_color_gradient(low = "blue3", high = "red") + 
  ggtitle("Temperatura") +
  xlab("Longitud") + ylab("Latitud")

g2 = ggplot(zona) + geom_sf() +
  geom_point(data = datos_geo, aes(x = Longitude, y = Latitude, size = Relative_Humidity,
                                   color = Relative_Humidity)) +
  theme_minimal() + scale_color_gradient(low = "blue", high = "red") +
  ggtitle("Humedad Relativa") +
  xlab("Longitud") + ylab("Latitud")

grid.arrange(g1,g2,ncol=2)

Figura 4: Relacion Espacial entre Humedad relativa y temperatura

Los mapas ilustrados en la Figura 4 muestran la distribución espacial de la Temperatura y la Humedad Relativa, revelando un patrón inverso entre ambas variables. Las temperaturas más altas, concentradas en el sur, coinciden con niveles más bajos de humedad relativa, mientras que las temperaturas más bajas, predominantes en el norte, están asociadas con una mayor humedad. Este comportamiento refuerza la relación negativa observada previamente, consistente con patrones climáticos típicos donde temperaturas elevadas reducen la humedad relativa debido a la mayor capacidad del aire para retener vapor de agua.

Análisis Estructural

Semivariogramas Experimentales

coordinates(datos_geo) <- ~Longitude + Latitude

# Semivariograma para Temperatura
vgm_temp <- variogram(Temperature ~ 1, datos_geo)
plot(vgm_temp, main = "Semivariograma Experimental - Temperatura",
     xlab = "Distancia", ylab = "Semivarianza", pch = 16)

Figura 5: Semivariogramas Experimentales Temperatura y Humedad Relativa

# Semivariograma para Humedad Relativa
vgm_hum <- variogram(Relative_Humidity ~ 1, datos_geo)
plot(vgm_hum, main = "Semivariograma Experimental - Humedad Relativa",
     xlab = "Distancia", ylab = "Semivarianza", pch = 16)

Figura 5: Semivariogramas Experimentales Temperatura y Humedad Relativa

De la Figura 5, el semivariograma experimental para la Temperatura muestra cómo la semivarianza aumenta conforme incrementa la distancia entre los puntos de muestreo. Esto indica una relación espacial clara, donde los valores de temperatura están más correlacionados en distancias más cortas y la correlación disminuye al aumentar la distancia. El patrón creciente sugiere la presencia de una estructura espacial que podría modelarse adecuadamente para entender mejor la variabilidad de la temperatura en la región estudiada.

El semivariograma experimental para la Humedad Relativa evidencia un aumento progresivo de la semivarianza con la distancia entre los puntos de muestreo. Este patrón indica que los valores de humedad relativa están más correlacionados en distancias cortas y que dicha correlación disminuye conforme la separación espacial aumenta. Este comportamiento sugiere una estructura espacial bien definida, que podría ser modelada para analizar y predecir la distribución espacial de la humedad relativa en la región estudiada.

Evaluación con Monte Carlo

temp <- as.geodata(datos_geo, data.col = 1)

# temp_res = semi(temp, tex = "Temperatura", c(50, 0.0008))

hume <- as.geodata(datos_geo, data.col = 2)

# hume_res <- semi(hume, tex = "Humedad Relativa", c(150, 0.0008))

## Bandas de no autocorrelación por Monte Carlo ----

vg_temp <- variog(temp)

## variog: computing omnidirectional variogram

vg_hume <- variog(hume)

## variog: computing omnidirectional variogram

# Generar simulaciones de Monte Carlo
mc_temp <- variog.mc.env(temp, obj.variog = vg_temp, nsim = 100)

## variog.env: generating 100 simulations by permutating data values
## variog.env: computing the empirical variogram for the 100 simulations
## variog.env: computing the envelops

mc_hume <- variog.mc.env(hume, obj.variog = vg_hume, nsim = 100)

## variog.env: generating 100 simulations by permutating data values
## variog.env: computing the empirical variogram for the 100 simulations
## variog.env: computing the envelops

# Graficar el semivariograma con las bandas de no correlación
plot(vg_temp, envelope = mc_temp, xlab = "Distancia", ylab = "Semivarianza",
     main = "Semivariograma con Bandas de Monte Carlo - Temperatura")

Figura 6: Semivariograma con evaluación montecarlo

plot(vg_hume, envelope = mc_hume, xlab = "Distancia", ylab = "Semivarianza",
     main = "Semivariograma con Bandas de Monte Carlo - Humedad rel.")

Figura 6: Semivariograma con evaluación montecarlo

En la Figura 6, el semivariograma con bandas de Monte Carlo para la Temperatura muestra que la semivarianza aumenta con la distancia, indicando una fuerte autocorrelación espacial positiva en distancias cortas. Los puntos de semivarianza experimental se encuentran en su mayoría por encima de la banda inferior de Monte Carlo, lo que confirma que los valores de temperatura en ubicaciones cercanas están más correlacionados que aquellos en ubicaciones más distantes, reflejando un patrón espacial característico.

El semivariograma con bandas de Monte Carlo para la Humedad Relativa revela un comportamiento similar al de la Temperatura, donde la semivarianza aumenta con la distancia, indicando una fuerte autocorrelación espacial en distancias cortas. Sin embargo, los valores de semivarianza para la humedad son significativamente más altos, reflejando una mayor variabilidad en esta variable en comparación con la temperatura. En ambos casos, la dependencia espacial disminuye a medida que aumenta la distancia, pero la estabilización de la semivarianza en la humedad ocurre a valores más altos.

Ajuste de Modelo Teórico

# Función para ajustar modelos

semis <- function(j, tex = '', initial = c(0, 1)) {
  # Carga de librerías necesarias
  if (!requireNamespace("geoR")) stop("Por favor instala el paquete 'geoR'")
  if (!requireNamespace("viridis")) stop("Por favor instala el paquete 'viridis'")
  if (!requireNamespace("raster")) stop("Por favor instala el paquete 'raster'")
  if (!requireNamespace("leaflet")) stop("Por favor instala el paquete 'leaflet'")
  
  # Modelos de covarianza disponibles
  models <- c("matern", "exponential", "gaussian", "spherical", 
              "circular", "cubic", "wave", "power", "powered.exponential")
  
  # Función para ajustar modelos y graficar líneas
  nmodels <- function(x, y, semi) {
    mod_ <- variofit(vario = semi, ini.cov.pars = y, cov.model = x)
    lines(mod_, col = viridis::viridis(length(models))[which(models == x)], lwd = 2)
    return(mod_$value)
  }
  
  
  
  # Exploración de las coordenadas
  # plot(j)
  options(scipen = 999)
  print(dist(j$coords))
  
  # Calcular las distancias y crear el semivariograma experimental
  s <- summary(dist(j$coords))
  distancias <- seq(s[1], s[5], length.out = 20)
  semi <- variog(geodata = j, uvec = distancias)
  
  # Graficar el semivariograma experimental
  plot(semi, xlab = 'Distancia', ylab = 'Semivarianza', panel.first = grid(),
       pch = 16, main = paste('Semivariograma de la', tex))
  
  # Ajustar modelos teóricos y seleccionar el mejor
  k <- lapply(models, nmodels, y = initial, semi = semi)
  best_model <- models[which.min(unlist(k))]
  print(paste("Mejor modelo:", best_model))
  
  # Ajustar el modelo teórico al semivariograma experimental
  mod_vario <- variofit(vario = semi, ini.cov.pars = initial, cov.model = best_model)
  
  # Graficar el semivariograma ajustado
  plot(semi, main = paste("Semivariograma de", tex), panel.first = grid(),
       ylab = 'Semivarianza', xlab = 'Distancia', pch = 16)
  lines(mod_vario, col = "blue", lwd = 2)
  
  # Predicción espacial con kriging
  x1 <- min(j$coords[, 1])
  x2 <- max(j$coords[, 1])
  y1 <- min(j$coords[, 2])
  y2 <- max(j$coords[, 2])
  
  puntos <- expand.grid(seq(x1, x2, length.out = 100),
                        seq(y1, y2, length.out = 100))
  
  plot(j$coords, pch = 16, xlab = 'Longitud', ylab = 'Latitud', panel.first = grid())
  points(puntos, col = "red")
  
  print(mod_vario$cov.pars)
  
  pred <- krige.conv(geodata = j, loc = puntos,
                     krige = krige.control(cov.pars = mod_vario$cov.pars))
  
  # Crear raster con predicciones
  res <- data.frame(puntos, pred$predict)
  im <- rasterFromXYZ(res)
  crs(im) <- "+init=EPSG:4326"
  
  # Graficar el raster
  plot(im, main = paste("Predicción espacial -", tex))
  
  # Aplicar máscara con shapefile (si se incluye)
  if (exists("shape")) {
    plot(mask(crop(im, shape), shape), main = "Predicción recortada")
  }
  
  # Visualización interactiva con Leaflet
  g <- leaflet() %>%
    addProviderTiles(providers$Esri.WorldImagery) %>%
    addRasterImage(mask(crop(im, shape), shape))
  
  print(g)
  
  # Devolver resultados
  return(list(res = res, mod_vario = mod_vario))
}


# semi_teo_temp <- semis(temp, "Temperatura", c(779792.2807, 1.4942))
# semi_teo_hume <- semis(hume, "Humedad Relativa", c(387473.1730, 1.0596))

variog: computing omnidirectional variogram
variofit: covariance model used is matern 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is exponential 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is gaussian 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is spherical 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is circular 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is cubic 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is wave 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is power 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is powered.exponential 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
[1] "Mejor modelo: gaussian"
variofit: covariance model used is gaussian 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
[1] 779792.2806999      0.1698245
krige.conv: model with constant mean
krige.conv: Kriging performed using global neighbourhood 

Figura 7: Modelos teoricos- Temperatura

Figura 8: Modelo ajustado gaussian- Temperatura

Interpretación de los Resultados

La salida del ajuste para la Temperatura indica que se evaluaron múltiples modelos de covarianza teóricos (Matern, Exponential, Gaussian, Spherical, Circular, Cubic, Wave, Power, y Powered Exponential), seleccionándose el modelo gaussiano como el que mejor representa la estructura de autocorrelación espacial de los datos, con base en la minimización de la función objetivo.

Los parámetros ajustados del modelo gaussiano fueron:

• Varianza parcial (sill): 779792.28

• Parámetro de rango (range): 0.1698

Esto sugiere que la dependencia espacial es significativa dentro de una distancia aproximada definida por el rango, y disminuye rápidamente fuera de este. La predicción espacial realizada mediante kriging empleó este modelo, lo que permitió estimar la variable de temperatura en ubicaciones no muestreadas. Este enfoque refuerza la utilidad del modelo gaussiano en la representación y predicción de patrones espaciales observados en la temperatura.

variog: computing omnidirectional variogram
variofit: covariance model used is matern 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is exponential 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is gaussian 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is spherical 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is circular 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is cubic 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is wave 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is power 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
variofit: covariance model used is powered.exponential 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
[1] "Mejor modelo: cubic"
variofit: covariance model used is cubic 
variofit: weights used: npairs 
variofit: minimisation function used: optim 
[1] 387473.17309048      0.09708241
krige.conv: model with constant mean
krige.conv: Kriging performed using global neighbourhood

Interpretación de los Resultados para Humedad Relativa

La salida del ajuste para la Humedad Relativa indica que se evaluaron diversos modelos de covarianza teóricos (Matern, Exponential, Gaussian, Spherical, Circular, Cubic, Wave, Power, y Powered Exponential), siendo seleccionado el modelo cúbico como el mejor para describir la estructura de autocorrelación espacial, según la minimización de la función objetivo.

Los parámetros ajustados del modelo cúbico fueron:

• Varianza parcial (sill): 387473.17

• Parámetro de rango (range): 0.0971

Esto sugiere que la dependencia espacial de la humedad relativa es significativa dentro de una distancia más corta, definida por el rango, en comparación con la temperatura. La predicción espacial, realizada mediante kriging utilizando este modelo, permitió estimar valores de humedad relativa en ubicaciones no muestreadas. El modelo cúbico, al capturar esta variabilidad espacial, resalta la influencia de factores locales sobre la distribución de la humedad en el área de estudio, diferenciándola de la temperatura, cuyo patrón espacial es más amplio.

Figura 9: Modelos teóricos - Humedad Relativa

Figura 10: Modelo ajustado cubic - Humedad Relativa

Predicción Espacial

# Definir los límites espaciales
x_range <- seq(min(temp$coords[, 1]), max(temp$coords[, 1]), length.out = 100)
y_range <- seq(min(temp$coords[, 2]), max(temp$coords[, 2]), length.out = 100)

grid <- expand.grid(x = x_range, y = y_range)

# Predicción espacial para Temperatura
pred_temp <- krige.conv(
  geodata = temp,
  loc = grid,
  krige = krige.control(
    cov.pars = c(779792.2807, 0.1698),  # Parámetros del modelo Gaussian
    nugget = 0  # Ajustar según el modelo
  )
)

## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhood

# Graficar predicción para Temperatura
res_temp <- data.frame(grid, pred_temp = pred_temp$predict)
im_temp <- rasterFromXYZ(res_temp)
plot(im_temp, main = "Predicción espacial de Temperatura")

Figura 11: Predicciones Espaciales

# Predicción espacial para Humedad Relativa
pred_hume <- krige.conv(
  geodata = hume,
  loc = grid,
  krige = krige.control(
    cov.pars = c(387473.1731, 0.0971),  # Parámetros del modelo Cubic
    nugget = 0  # Ajustar según el modelo
  )
)

## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhood

# Graficar predicción para Humedad Relativa
res_hume <- data.frame(grid, pred_hume = pred_hume$predict)
im_hume <- rasterFromXYZ(res_hume)
plot(im_hume, main = "Predicción espacial de Humedad Relativa")

Figura 11: Predicciones Espaciales

De la Figura 11, el mapa de predicción espacial de la Temperatura muestra un gradiente claro, donde las temperaturas más bajas (20-25 °C) se concentran en el noroeste, representadas por tonos rosados, mientras que las temperaturas más altas (30 °C o más) predominan en el sureste, con tonos verdes. Las áreas centrales presentan temperaturas moderadas, entre 25 y 30 °C, representadas en amarillo. Este patrón espacial sugiere una transición gradual de temperaturas, posiblemente influenciada por factores geográficos o ambientales, como altitud o exposición solar, proporcionando un panorama continuo de las variaciones térmicas en el área estudiada.

El mapa de predicción espacial de la Humedad Relativa muestra un gradiente donde los niveles más altos (70-80%) se concentran en el noroeste, representados por tonos verdes, mientras que los niveles más bajos (30-50%) predominan hacia el sureste, con tonos rosados. En las áreas centrales, los niveles de humedad relativa se encuentran en un rango intermedio (50-70%), marcados en amarillo. Este patrón espacial es inverso al observado en la temperatura, lo que refuerza la relación negativa entre ambas variables. La distribución refleja cómo factores locales, como la topografía o la vegetación, podrían estar influyendo en la variabilidad de la humedad en la región estudiada.