Informe geoestadístico de temperatura en árboles de aguacate

library(readxl)
library(dplyr)
library(lubridate)
library(ggplot2)
library(sf)
library(sp)
library(gstat)
library(automap)
library(viridis)
library(knitr)
library(tidyr)
library(kableExtra)

1 .Introducción

El presente informe desarrolla un análisis geoestadístico aplicado a datos de medición de árboles de producción de aguacate. El objetivo es evaluar si la variable ambiental Temperature presenta autocorrelación espacial en las observaciones correspondientes al día 01/10/2020.

Para ello, se realiza primero un análisis exploratorio de los datos y de su distribución espacial. Posteriormente, se evalúa si existe una tendencia espacial previa en la variable, debido a que la presencia de tendencia puede afectar el supuesto de estacionariedad requerido para el análisis geoestadístico tradicional.

Luego, se calcula el semivariograma experimental, se ajustan tres modelos teóricos —esférico, exponencial y gaussiano— y se selecciona el mejor modelo con base en el menor error de ajuste. Finalmente, se realiza una predicción espacial mediante Kriging ordinario y se construye el mapa de predicción de temperatura.

#Carga de datos
ruta <- "C:/Users/Erick Caicedo Ruiz/Desktop/Maestría Ciencia de Datos/Datos_Completos_Aguacate.xlsx"

datos <- read_excel(ruta)

glimpse(datos)

#filtro de datos observaciones del 01/10/2020

datos_filtrados <- datos %>%
  mutate(
    FORMATTED_DATE_TIME = as.Date(FORMATTED_DATE_TIME, format = "%d/%m/%Y")
  ) %>%
  filter(FORMATTED_DATE_TIME == as.Date("2020-10-01")) %>%
  filter(
    !is.na(Temperature),
    !is.na(Latitude),
    !is.na(Longitude)
  )

nrow(datos_filtrados)

head(datos_filtrados)

2 .Análisis exploratorio de la variable Temperature

El análisis exploratorio de la variable Temperature se realizó con un total de 534 observaciones correspondientes al 01/10/2020. Los valores de temperatura se encuentran entre 22.20 y 29.70, con una media de 25.83 y una mediana de 25.80. La cercanía entre la media y la mediana indica que la distribución de la temperatura es relativamente equilibrada y no presenta una asimetría marcada.

El primer cuartil se ubica en 24.50 y el tercer cuartil en 27.17, lo que significa que el 50% central de las observaciones se concentra aproximadamente entre estos dos valores. Esto muestra que la mayor parte de las mediciones se agrupa en un rango moderado de temperatura, sin variaciones extremas muy pronunciadas. La desviación estándar de 1.77 también indica una dispersión relativamente baja, lo que sugiere que las temperaturas observadas son bastante consistentes entre los puntos medidos (ver tabla 1).

tabla_resumen_temperature <- datos_filtrados %>%
  summarise(
    `Número de observaciones` = n(),
    `Valor mínimo` = min(Temperature, na.rm = TRUE),
    `Primer cuartil` = quantile(Temperature, 0.25, na.rm = TRUE),
    `Mediana` = median(Temperature, na.rm = TRUE),
    `Media` = mean(Temperature, na.rm = TRUE),
    `Tercer cuartil` = quantile(Temperature, 0.75, na.rm = TRUE),
    `Valor máximo` = max(Temperature, na.rm = TRUE),
    `Desviación estándar` = sd(Temperature, na.rm = TRUE)
  ) %>%
  pivot_longer(
    cols = everything(),
    names_to = "Estadístico",
    values_to = "Valor"
  ) %>%
  mutate(
    Valor = round(Valor, 2)
  )

kable(
  tabla_resumen_temperature,
  caption = "Resumen estadístico de la variable Temperature",
  align = c("l", "c")
) %>%
  kable_styling(
    full_width = FALSE,
    bootstrap_options = c("striped", "hover", "condensed")
  )

Resumen estadístico de la variable Temperature

Estadístico	Valor
Número de observaciones	534.00
Valor mínimo	22.20
Primer cuartil	24.50
Mediana	25.80
Media	25.83
Tercer cuartil	27.17
Valor máximo	29.70
Desviación estándar	1.77

El histograma muestra que la mayor frecuencia de observaciones se encuentra alrededor de valores cercanos a 25.5 y 26.0. Aunque existen registros hacia valores más bajos, cercanos a 22, y más altos, cercanos a 30, estos no dominan la distribución. En general, la forma del histograma permite observar una distribución con concentración en valores medios y sin presencia evidente de valores atípicos extremos.

ggplot(datos_filtrados, aes(x = Temperature)) +
  geom_histogram(bins = 30, color = "white") +
  labs(
    title = "Distribución de la variable Temperature",
    x = "Temperature",
    y = "Frecuencia"
  ) +
  theme_minimal()

El diagrama de caja confirma esta lectura. La caja se extiende aproximadamente entre 24.5 y 27.2, con una mediana cercana a 25.8. Los bigotes cubren el rango general de la variable, desde valores cercanos a 22.2 hasta 29.7, sin mostrar puntos atípicos visibles. Esto indica que las variaciones de temperatura se mantienen dentro de un rango esperado y que no hay registros extremos que puedan distorsionar de manera importante el análisis posterior.

En conjunto, estos resultados muestran que Temperature es una variable adecuada para continuar con el análisis geoestadístico. Presenta suficiente variabilidad para evaluar patrones espaciales, pero no evidencia valores extremos fuertes que puedan afectar de forma considerable el cálculo del semivariograma o la predicción espacial mediante Kriging.

ggplot(datos_filtrados, aes(y = Temperature)) +
  geom_boxplot() +
  labs(
    title = "Diagrama de caja de la variable Temperature",
    y = "Temperature"
  ) +
  theme_minimal()

3 .Distribución espacial de los puntos observados

La distribución espacial de los puntos observados muestra que las mediciones de Temperature se concentran en un área relativamente delimitada, con una forma alargada y una cobertura densa de puntos. Esto es positivo para el análisis geoestadístico, porque existe una cantidad suficiente de observaciones distribuidas en el espacio para evaluar posibles patrones de autocorrelación espacial.

En términos visuales, se observa que la temperatura no está distribuida completamente al azar. Hay zonas con valores más altos, representadas por tonos amarillos y verdes claros, especialmente hacia algunos bordes del área de estudio. También se identifican zonas con temperaturas más bajas, representadas por tonos morados y azules, principalmente en sectores internos y hacia la parte central-inferior del conjunto de puntos.

Este patrón sugiere que existen agrupamientos espaciales de valores similares. Es decir, algunos puntos cercanos tienden a compartir temperaturas parecidas, mientras que otros sectores presentan contrastes más marcados. Esta primera lectura visual es importante porque anticipa la posible existencia de autocorrelación espacial en la variable.

El mapa también permite observar que no todos los sectores tienen el mismo comportamiento térmico. Mientras algunas zonas presentan temperaturas cercanas o superiores a 28, otras se ubican alrededor de valores cercanos a 23 o 24. Esta variabilidad espacial justifica continuar con el cálculo del semivariograma experimental, ya que este permitirá evaluar formalmente si la similitud entre temperaturas disminuye a medida que aumenta la distancia entre los puntos.

En conjunto, la distribución espacial observada indica que Temperature presenta diferencias territoriales dentro del área de medición y que existen condiciones adecuadas para avanzar hacia el análisis geoestadístico mediante semivariograma, ajuste de modelos teóricos y predicción espacial con Kriging.

ggplot(datos_filtrados, aes(x = Longitude, y = Latitude, color = Temperature)) +
  geom_point(size = 3, alpha = 0.8) +
  scale_color_viridis_c() +
  labs(
    title = "Distribución espacial de Temperature observada",
    x = "Longitud",
    y = "Latitud",
    color = "Temperature"
  ) +
  theme_minimal()

#Conversión de datos a objeto espacial: Para el análisis geoestadístico es necesario trabajar con un sistema de coordenadas proyectadas en unidades métricas. Por esta razón, las coordenadas originales en latitud y longitud se convierten primero a un objeto espacial y luego se proyectan a un sistema UTM.

datos_sf <- st_as_sf(
  datos_filtrados,
  coords = c("Longitude", "Latitude"),
  crs = 4326
)

datos_sf_utm <- st_transform(datos_sf, crs = 32618)

datos_sp <- as(datos_sf_utm, "Spatial")

4 .Análisis Exploratorio Espacial de Tendencia

Antes de calcular el semivariograma, se evalúa si la variable Temperature presenta una tendencia espacial asociada a las coordenadas. Este paso es importante porque el semivariograma clásico y el Kriging ordinario asumen que la media de la variable es aproximadamente constante en el espacio.

Si existe una tendencia espacial fuerte, es recomendable analizar los residuos de un modelo de tendencia o aplicar Kriging universal. Si no se identifica una tendencia fuerte, se puede continuar con el semivariograma de la variable original y con Kriging ordinario.

datos_coord <- datos_sf_utm %>%
  mutate(
    x = st_coordinates(.)[, 1],
    y = st_coordinates(.)[, 2]
  ) %>%
  st_drop_geometry()

tabla_coordenadas <- datos_coord %>%
  select(
    FORMATTED_DATE_TIME,
    Temperature,
    x,
    y
  ) %>%
  head(10) %>%
  mutate(
    Temperature = round(Temperature, 2),
    x = round(x, 2),
    y = round(y, 2)
  ) %>%
  rename(
    `Fecha de observación` = FORMATTED_DATE_TIME,
    `Temperatura` = Temperature,
    `Coordenada Este` = x,
    `Coordenada Norte` = y
  )

kable(
  tabla_coordenadas,
  caption = "Primeras observaciones con coordenadas proyectadas para el análisis espacial",
  align = c("c", "c", "c", "c")
) %>%
  kable_styling(
    full_width = FALSE,
    bootstrap_options = c("striped", "hover", "condensed")
  )

4.1 .Relación entre Temperature y coordenada Este

La relación entre Temperature y la coordenada Este muestra una tendencia negativa leve. Esto significa que, a medida que los puntos se desplazan hacia el este dentro del área de estudio, la temperatura tiende a disminuir ligeramente.

La línea de tendencia presenta una pendiente descendente, lo cual sugiere que existe cierto gradiente espacial en dirección este-oeste. Sin embargo, la nube de puntos muestra una dispersión amplia, con temperaturas altas y bajas distribuidas a lo largo de casi todo el rango de la coordenada Este. Por esta razón, la relación no parece ser fuerte, sino más bien moderada o débil.

Este resultado indica que la ubicación en el eje Este puede explicar una parte del comportamiento espacial de la temperatura, pero no toda su variabilidad. Es decir, aunque existe una tendencia general de disminución hacia el este, también hay variaciones locales importantes que deben ser analizadas mediante el semivariograma.

En términos geoestadísticos, este gráfico es útil porque permite identificar una posible tendencia espacial previa al cálculo de la autocorrelación. Si el modelo exploratorio de tendencia confirma que esta relación es estadísticamente significativa y tiene un poder explicativo relevante, podría considerarse el uso de residuos o Kriging universal. Si la tendencia es débil, se puede continuar con el semivariograma de la variable original y aplicar Kriging ordinario.

ggplot(datos_coord, aes(x = x, y = Temperature)) +
  geom_point(alpha = 0.7) +
  geom_smooth(method = "lm", se = TRUE) +
  labs(
    title = "Relación entre Temperature y coordenada Este",
    x = "Coordenada Este",
    y = "Temperature"
  ) +
  theme_minimal()

4.2 .Relación entre Temperature y coordenada Norte

La relación entre Temperature y la coordenada Norte también muestra una tendencia negativa leve. Esto indica que, a medida que los puntos se ubican más hacia el norte dentro del área de estudio, la temperatura tiende a disminuir ligeramente.

La línea de tendencia tiene una pendiente descendente, similar a lo observado con la coordenada Este. Sin embargo, la nube de puntos presenta una dispersión amplia, con valores de temperatura altos, medios y bajos distribuidos a lo largo de casi todo el rango de la coordenada Norte. Esto sugiere que la relación entre la temperatura y la ubicación norte-sur existe, pero no es suficientemente fuerte para explicar por sí sola toda la variabilidad de la variable.

En términos exploratorios, este resultado permite identificar una posible tendencia espacial general de la temperatura. No obstante, la dispersión de los puntos muestra que también hay variaciones locales importantes que no dependen únicamente de la posición norte-sur. Por esta razón, es necesario continuar con el análisis del semivariograma, ya que este permite evaluar formalmente si los valores de temperatura presentan autocorrelación espacial.

En conjunto con el gráfico anterior de la coordenada Este, se observa que la temperatura presenta un comportamiento espacial con una ligera disminución en ambas direcciones proyectadas. Esto puede indicar la existencia de un gradiente espacial suave dentro del área de estudio. Sin embargo, dado que la tendencia no parece ser muy marcada visualmente, el análisis geoestadístico mediante semivariograma y Kriging ordinario sigue siendo pertinente.

ggplot(datos_coord, aes(x = y, y = Temperature)) +
  geom_point(alpha = 0.7) +
  geom_smooth(method = "lm", se = TRUE) +
  labs(
    title = "Relación entre Temperature y coordenada Norte",
    x = "Coordenada Norte",
    y = "Temperature"
  ) +
  theme_minimal()

4.3 .Superficie exploratoria de tendencia espacial

La superficie exploratoria de tendencia espacial permite observar la distribución de Temperature usando las coordenadas proyectadas Este y Norte. En el mapa se evidencia que la temperatura no se distribuye de manera completamente homogénea dentro del área de estudio, sino que existen sectores con valores relativamente altos y otros con valores más bajos.

Los valores más altos de temperatura, representados por tonos amarillos y verdes claros, se concentran principalmente hacia algunos bordes del área de observación, especialmente en sectores del costado inferior izquierdo y en algunos puntos del costado derecho. Por el contrario, los valores más bajos, representados por tonos morados y azules oscuros, aparecen principalmente en zonas internas y hacia el sector central-inferior del área.

Este comportamiento sugiere la presencia de una estructura espacial en la variable, ya que se observan agrupamientos de puntos con temperaturas similares. Es decir, algunos sectores tienden a presentar temperaturas más altas o más bajas de forma localizada, lo cual indica que la ubicación espacial influye en el comportamiento de la variable.

Sin embargo, también se observa una mezcla de valores en distintas partes del área, por lo que la tendencia no parece ser completamente uniforme ni estrictamente lineal. Esto significa que, además de una posible tendencia espacial general, existen variaciones locales que deben ser evaluadas mediante el semivariograma experimental.

En conjunto, la figura respalda la necesidad de continuar con el análisis geoestadístico, ya que la variable Temperature presenta diferencias espaciales visibles y posibles agrupamientos locales. Esta exploración previa permite justificar el cálculo del semivariograma y el posterior ajuste de modelos teóricos para aplicar Kriging.

ggplot(datos_coord, aes(x = x, y = y, color = Temperature)) +
  geom_point(size = 3, alpha = 0.85) +
  scale_color_viridis_c() +
  labs(
    title = "Exploración espacial de Temperature",
    subtitle = "Distribución espacial de la temperatura observada",
    x = "Coordenada Este",
    y = "Coordenada Norte",
    color = "Temperature"
  ) +
  theme_minimal()

## .Modelo exploratorio de tendencia
modelo_tendencia <- lm(Temperature ~ x + y, data = datos_coord)

summary(modelo_tendencia)

r2_tendencia <- summary(modelo_tendencia)$r.squared

p_modelo <- pf(
  summary(modelo_tendencia)$fstatistic[1],
  summary(modelo_tendencia)$fstatistic[2],
  summary(modelo_tendencia)$fstatistic[3],
  lower.tail = FALSE
)

r2_tendencia
p_modelo

if (p_modelo < 0.05 & r2_tendencia > 0.20) {
  mensaje_tendencia <- "Se identifica evidencia de tendencia espacial en Temperature. Se recomienda analizar residuos o aplicar Kriging universal."
} else {
  mensaje_tendencia <- "No se identifica una tendencia espacial fuerte en Temperature. Se puede continuar con semivariograma y Kriging ordinario."
}

mensaje_tendencia

5 .Evaluación de autocorrelación espacial mediante semivariograma

El semivariograma experimental permite analizar si existe dependencia espacial entre las observaciones. En términos generales, si los puntos cercanos tienen valores de temperatura más parecidos que los puntos lejanos, se esperaría observar una estructura espacial en el semivariograma.

variograma_exp <- variogram(Temperature ~ 1, datos_sp)

plot(
  variograma_exp,
  main = "Semivariograma experimental de Temperature"
)

El semivariograma experimental de la variable Temperature evidencia una estructura de autocorrelación espacial. Se observa que la semivarianza aumenta conforme crece la distancia entre los puntos de medición, lo que indica que las observaciones cercanas tienden a presentar valores de temperatura más similares entre sí que aquellas ubicadas a mayor distancia.

El incremento de la semivarianza es más marcado en las primeras distancias, lo que sugiere una dependencia espacial fuerte a escala local. Posteriormente, la semivarianza tiende a estabilizarse alrededor de valores cercanos a 3.0 y 3.2, aproximadamente a partir de una distancia entre 55 y 65 unidades. Esta estabilización puede interpretarse como la presencia de una meseta o sill, indicando que a partir de esa distancia la autocorrelación espacial disminuye considerablemente.

El gráfico también sugiere un posible efecto pepita, dado que la semivarianza inicial no parte de cero. Esto puede estar asociado a variabilidad local no explicada, errores de medición o cambios ambientales a escalas muy pequeñas.

En conjunto, el semivariograma permite concluir que la variable Temperature presenta dependencia espacial, por lo cual resulta pertinente continuar con el ajuste de modelos teóricos de semivariograma y aplicar una interpolación espacial mediante Kriging.

6 .Ajuste de modelos teóricos de semivariograma

Se ajustan tres modelos teóricos comúnmente utilizados en geoestadística:

-Modelo esférico. -Modelo exponencial. -Modelo gaussiano.

La tabla de ajuste de modelos teóricos muestra los parámetros estimados para los modelos esférico, exponencial y gaussiano aplicados al semivariograma experimental de Temperature.

En los tres modelos, el efecto pepita es igual a 0. Esto indica que, bajo los modelos ajustados, no se identifica una variabilidad inicial importante atribuible a errores de medición o a microvariabilidad espacial no capturada por la escala de muestreo. En otras palabras, el comportamiento espacial de la temperatura parece estar explicado principalmente por la estructura espacial modelada.

El modelo esférico presenta una meseta total de 2.8054 y un rango de 21.6902. Esto significa que, según este modelo, la dependencia espacial de la temperatura se mantendría aproximadamente hasta una distancia de 21.69 unidades. Después de esa distancia, la autocorrelación espacial tendería a estabilizarse. Su error de ajuste SSErr es de 4.9853.

El modelo exponencial presenta una meseta total de 2.9769 y un rango de 10.4440. Aunque su rango es menor que el del modelo esférico, este modelo alcanza el mejor ajuste numérico, con un SSErr de 1.6714. Esto indica que el modelo exponencial es el que mejor representa la estructura espacial observada en el semivariograma experimental, según el criterio de menor error.

El modelo gaussiano presenta una meseta total de 2.5928 y un rango de 7.4306, pero tiene el mayor error de ajuste, con un SSErr de 12.3110. Esto sugiere que, entre los tres modelos evaluados, el gaussiano es el que menos se ajusta al comportamiento espacial observado de la temperatura.

modelo_esferico <- fit.variogram(
  variograma_exp,
  model = vgm(model = "Sph")
)

modelo_exponencial <- fit.variogram(
  variograma_exp,
  model = vgm(model = "Exp")
)

modelo_gaussiano <- fit.variogram(
  variograma_exp,
  model = vgm(model = "Gau")
)

extraer_parametros_modelo <- function(modelo, nombre_modelo) {
  
  modelo_df <- as.data.frame(modelo)
  
  efecto_pepita <- modelo_df$psill[modelo_df$model == "Nug"]
  meseta_parcial <- modelo_df$psill[modelo_df$model != "Nug"]
  rango <- modelo_df$range[modelo_df$model != "Nug"]
  
  if (length(efecto_pepita) == 0) {
    efecto_pepita <- 0
  }
  
  if (length(meseta_parcial) == 0) {
    meseta_parcial <- NA
  }
  
  if (length(rango) == 0) {
    rango <- NA
  }
  
  data.frame(
    Modelo = nombre_modelo,
    `Efecto pepita` = efecto_pepita,
    `Meseta parcial` = meseta_parcial,
    `Meseta total` = efecto_pepita + meseta_parcial,
    Rango = rango,
    `Error de ajuste SSErr` = attr(modelo, "SSErr")
  )
}

tabla_modelos <- bind_rows(
  extraer_parametros_modelo(modelo_esferico, "Esférico"),
  extraer_parametros_modelo(modelo_exponencial, "Exponencial"),
  extraer_parametros_modelo(modelo_gaussiano, "Gaussiano")
) %>%
  mutate(
    across(
      where(is.numeric),
      ~ round(.x, 4)
    )
  )

kable(
  tabla_modelos,
  caption = "Parámetros de los modelos teóricos ajustados al semivariograma experimental",
  align = c("l", "c", "c", "c", "c", "c")
) %>%
  kable_styling(
    full_width = FALSE,
    bootstrap_options = c("striped", "hover", "condensed")
  )

Parámetros de los modelos teóricos ajustados al semivariograma experimental

Modelo	Efecto.pepita	Meseta.parcial	Meseta.total	Rango	Error.de.ajuste.SSErr
Esférico	0	2.8054	2.8054	21.6902	4.9853
Exponencial	0	2.9769	2.9769	10.4440	1.6714
Gaussiano	0	2.5928	2.5928	7.4306	12.3110

comparacion_modelos <- data.frame(
  Modelo = c("Esférico", "Exponencial", "Gaussiano"),
  SSErr = c(
    attr(modelo_esferico, "SSErr"),
    attr(modelo_exponencial, "SSErr"),
    attr(modelo_gaussiano, "SSErr")
  )
)

kable(
  comparacion_modelos,
  caption = "Comparación de modelos teóricos de semivariograma"
)

Comparación de modelos teóricos de semivariograma

Modelo	SSErr
Esférico	4.985254
Exponencial	1.671427
Gaussiano	12.310972

mejor_modelo_nombre <- comparacion_modelos$Modelo[which.min(comparacion_modelos$SSErr)]

mejor_modelo <- switch(
  mejor_modelo_nombre,
  "Esférico" = modelo_esferico,
  "Exponencial" = modelo_exponencial,
  "Gaussiano" = modelo_gaussiano
)

En conjunto, los resultados permiten seleccionar el modelo exponencial como el modelo teórico más adecuado para continuar con la predicción espacial mediante Kriging ordinario. Esta selección se justifica porque presenta el menor valor de SSErr, lo que indica una mejor correspondencia entre el semivariograma experimental y la curva teórica ajustada.

La comparación gráfica de los modelos teóricos ajustados al semivariograma experimental muestra que el modelo exponencial presenta el mejor ajuste según el menor valor de SSErr. Este modelo reproduce adecuadamente el aumento progresivo de la semivarianza conforme se incrementa la distancia entre los puntos de medición.

El modelo esférico también representa parcialmente la estructura espacial, especialmente en distancias cortas e intermedias, pero alcanza la meseta de manera más rápida, lo que puede limitar su capacidad para representar la continuidad espacial de la temperatura en distancias mayores. Por su parte, el modelo gaussiano presenta una curva más suave, pero tiende a estabilizarse en valores de semivarianza inferiores a varios puntos experimentales observados en distancias medias y largas.

En general, el comportamiento del semivariograma confirma la presencia de autocorrelación espacial en la variable Temperature, ya que la semivarianza aumenta con la distancia y posteriormente tiende a estabilizarse. Esto indica que los puntos cercanos presentan temperaturas más similares entre sí, mientras que la similitud disminuye a medida que aumenta la distancia. Por esta razón, resulta adecuado utilizar el modelo exponencial para realizar la interpolación espacial mediante Kriging ordinario.

# Crear líneas teóricas para cada modelo ajustado
linea_esferico <- variogramLine(
  modelo_esferico,
  maxdist = max(variograma_exp$dist),
  n = 100
) %>%
  mutate(Modelo = "Esférico")

linea_exponencial <- variogramLine(
  modelo_exponencial,
  maxdist = max(variograma_exp$dist),
  n = 100
) %>%
  mutate(Modelo = "Exponencial")

linea_gaussiano <- variogramLine(
  modelo_gaussiano,
  maxdist = max(variograma_exp$dist),
  n = 100
) %>%
  mutate(Modelo = "Gaussiano")

lineas_modelos <- bind_rows(
  linea_esferico,
  linea_exponencial,
  linea_gaussiano
)

# Graficar semivariograma experimental y modelos teóricos
ggplot() +
  geom_point(
    data = variograma_exp,
    aes(x = dist, y = gamma),
    size = 2
  ) +
  geom_line(
    data = lineas_modelos,
    aes(x = dist, y = gamma, linetype = Modelo),
    linewidth = 1
  ) +
  labs(
    title = "Comparación de modelos teóricos del semivariograma",
    subtitle = paste("Mejor modelo según menor SSErr:", mejor_modelo_nombre),
    x = "Distancia",
    y = "Semivarianza",
    linetype = "Modelo teórico"
  ) +
  theme_minimal()

7 .Validación cruzada del modelo

La validación cruzada del modelo permite evaluar qué tan bien el Kriging predice los valores observados de Temperature. Para ello, se comparan los valores reales con los valores estimados por el modelo y se analizan los residuos, entendidos como la diferencia entre el valor observado y el valor predicho.

Los residuos presentan una media de -0.0104, valor muy cercano a cero. Esto indica que el modelo no tiene un sesgo sistemático importante; es decir, no tiende a sobreestimar ni a subestimar de forma marcada la temperatura. La mediana también se encuentra muy cerca de cero, lo que refuerza la idea de que los errores de predicción están relativamente equilibrados.

El rango de los residuos va aproximadamente desde -3.73 hasta 2.81. Esto muestra que, aunque la mayoría de los errores son moderados, existen algunos casos donde la diferencia entre la temperatura observada y la predicha es más alta. Sin embargo, estos valores parecen ser casos puntuales y no dominan el comportamiento general del modelo.

validacion <- krige.cv(
  Temperature ~ 1,
  datos_sp,
  model = mejor_modelo
)

tabla_resumen_residuos <- data.frame(
  Estadístico = c(
    "Número de observaciones",
    "Valor mínimo",
    "Primer cuartil",
    "Mediana",
    "Media",
    "Tercer cuartil",
    "Valor máximo",
    "Desviación estándar",
    "Error medio absoluto",
    "Raíz del error cuadrático medio"
  ),
  Valor = c(
    length(validacion$residual),
    min(validacion$residual, na.rm = TRUE),
    quantile(validacion$residual, 0.25, na.rm = TRUE),
    median(validacion$residual, na.rm = TRUE),
    mean(validacion$residual, na.rm = TRUE),
    quantile(validacion$residual, 0.75, na.rm = TRUE),
    max(validacion$residual, na.rm = TRUE),
    sd(validacion$residual, na.rm = TRUE),
    mean(abs(validacion$residual), na.rm = TRUE),
    sqrt(mean(validacion$residual^2, na.rm = TRUE))
  )
) %>%
  mutate(
    Valor = round(Valor, 4)
  )

kable(
  tabla_resumen_residuos,
  caption = "Resumen estadístico de los residuos de la validación cruzada del Kriging",
  align = c("l", "c")
) %>%
  kable_styling(
    full_width = FALSE,
    bootstrap_options = c("striped", "hover", "condensed")
  )

Resumen estadístico de los residuos de la validación cruzada del Kriging

Estadístico	Valor
Número de observaciones	534.0000
Valor mínimo	-3.7305
Primer cuartil	-0.6379
Mediana	0.0009
Media	-0.0104
Tercer cuartil	0.6041
Valor máximo	2.8148
Desviación estándar	1.0163
Error medio absoluto	0.7875
Raíz del error cuadrático medio	1.0154

El histograma de residuos muestra una distribución concentrada alrededor de cero, con forma aproximadamente simétrica. Esto es un resultado favorable, ya que sugiere que la mayor parte de las predicciones presentan errores pequeños y que el modelo tiene un comportamiento razonablemente estable.

validacion_df <- as.data.frame(validacion)

ggplot(validacion_df, aes(x = residual)) +
  geom_histogram(bins = 30, color = "white") +
  labs(
    title = "Distribución de residuos de validación cruzada",
    x = "Residuo",
    y = "Frecuencia"
  ) +
  theme_minimal()

Por su parte, el gráfico de valores observados frente a valores predichos muestra una relación positiva clara. La mayoría de los puntos se distribuyen cerca de la línea diagonal, lo que indica que el modelo logra reproducir adecuadamente el patrón general de la temperatura observada. Sin embargo, también se observa cierta dispersión alrededor de la línea, especialmente en los valores medios y altos, lo que indica que el modelo no predice perfectamente todos los casos.

En conjunto, la validación cruzada indica que el modelo de Kriging tiene un desempeño aceptable para la predicción espacial de Temperature. La media de los residuos cercana a cero, la concentración de errores alrededor de cero y la relación positiva entre valores observados y predichos respaldan el uso del modelo seleccionado para generar la superficie de predicción espacial.

ggplot(validacion_df, aes(x = observed, y = var1.pred)) +
  geom_point(alpha = 0.7) +
  geom_abline(intercept = 0, slope = 1, linetype = "dashed") +
  labs(
    title = "Valores observados vs valores predichos",
    x = "Temperature observada",
    y = "Temperature predicha"
  ) +
  theme_minimal()

#Creación de grilla de predicción
#Para realizar la interpolación espacial mediante Kriging, se construye una grilla regular sobre el área donde se encuentran los puntos observados.

bbox <- st_bbox(datos_sf_utm)

grilla <- expand.grid(
  x = seq(bbox["xmin"], bbox["xmax"], length.out = 100),
  y = seq(bbox["ymin"], bbox["ymax"], length.out = 100)
)

coordinates(grilla) <- ~ x + y
proj4string(grilla) <- proj4string(datos_sp)
gridded(grilla) <- TRUE

kriging_resultado <- krige(
  Temperature ~ 1,
  datos_sp,
  grilla,
  model = mejor_modelo
)

kriging_df <- as.data.frame(kriging_resultado)

head(kriging_df)

El mapa de predicción espacial mediante Kriging ordinario muestra la distribución estimada de Temperature en el área de estudio a partir del modelo exponencial, que fue seleccionado como el mejor ajuste del semivariograma experimental.

La superficie interpolada evidencia una variación espacial clara de la temperatura. Las zonas representadas con tonos amarillos y verdes claros corresponden a sectores con temperaturas predichas más altas, cercanas o superiores a 28. Estas áreas se observan principalmente hacia el sector inferior izquierdo y en algunas franjas del costado derecho del área de análisis.

Por el contrario, los tonos morados y azules oscuros representan zonas con temperaturas predichas más bajas, cercanas a valores entre 23 y 25. Estas áreas aparecen principalmente en sectores internos del área de estudio, especialmente hacia la zona central e inferior-central.

El patrón espacial obtenido muestra que la temperatura no se distribuye de manera uniforme, sino que presenta franjas y agrupamientos de valores similares. Esto es coherente con los resultados del semivariograma, donde se identificó autocorrelación espacial: los puntos cercanos tienden a presentar temperaturas parecidas y la similitud disminuye con la distancia.

Los puntos negros representan las observaciones originales utilizadas para construir la predicción. Se observa que la superficie de Kriging responde a la distribución de estos puntos, generando una estimación continua de temperatura para zonas donde no se tienen mediciones directas.

En términos generales, el mapa permite identificar áreas con mayor y menor temperatura dentro del cultivo o zona de medición. Esta información puede ser útil para reconocer patrones ambientales internos, posibles microclimas o sectores que podrían requerir análisis complementarios en relación con la producción de aguacate.

puntos_df <- as.data.frame(coordinates(datos_sp))
names(puntos_df) <- c("x", "y")

puntos_df <- puntos_df %>%
  bind_cols(as.data.frame(datos_sp))

ggplot() +
  geom_raster(
    data = kriging_df,
    aes(x = x, y = y, fill = var1.pred)
  ) +
geom_point(
  data = puntos_df,
  aes(x = x, y = y),
  size = 1,
  color = "black"
) +
  scale_fill_viridis_c() +
  labs(
    title = paste("Predicción espacial de Temperature mediante Kriging - Modelo", mejor_modelo_nombre),
    x = "Coordenada Este",
    y = "Coordenada Norte",
    fill = "Temperature predicha"
  ) +
  theme_minimal()

El mapa de varianza de predicción del Kriging muestra el nivel de incertidumbre asociado a la estimación espacial de Temperature. A diferencia del mapa de predicción, que muestra los valores estimados de temperatura, este mapa permite identificar qué tan confiable es la predicción en cada zona del área analizada.

Las zonas con tonos morados y azules corresponden a áreas con menor varianza de predicción, es decir, donde el modelo tiene menor incertidumbre. Estas zonas coinciden principalmente con los sectores donde existe una mayor concentración de puntos observados. Esto ocurre porque el Kriging tiene más información cercana para estimar la temperatura, lo que mejora la precisión de la predicción.

Por el contrario, las zonas con tonos verdes claros y amarillos presentan mayor varianza de predicción. Estas áreas se ubican especialmente en los bordes y esquinas del mapa, donde hay menor presencia o ausencia de puntos de medición. En estos sectores, el modelo debe estimar la temperatura con menos información cercana, por lo que la incertidumbre aumenta.

Este resultado es esperable en un análisis de Kriging, ya que la precisión de la predicción depende en gran medida de la cercanía y densidad de las observaciones disponibles. Las áreas más alejadas de los puntos medidos tienden a presentar mayor incertidumbre, mientras que las áreas con mayor densidad de datos presentan predicciones más confiables.

En conjunto, el mapa confirma que la predicción espacial es más robusta dentro del área cubierta por los puntos observados y menos confiable hacia los bordes externos de la grilla de interpolación. Por esta razón, la interpretación del mapa de temperatura predicha debe hacerse con mayor cautela en las zonas amarillas, especialmente en los extremos donde no existen observaciones cercanas.

ggplot() +
  geom_raster(
    data = kriging_df,
    aes(x = x, y = y, fill = var1.var)
  ) +
  geom_point(
    data = puntos_df,
    aes(x = coords.x1, y = coords.x2),
    size = 1,
    color = "black"
  ) +
  scale_fill_viridis_c() +
  labs(
    title = "Varianza de predicción del Kriging",
    x = "Coordenada Este",
    y = "Coordenada Norte",
    fill = "Varianza"
  ) +
  theme_minimal()

8 .Conclusiones generales

El análisis geoestadístico realizado permitió evaluar el comportamiento espacial de la variable Temperature en las observaciones de árboles de producción de aguacate correspondientes al 01/10/2020. A partir de 534 registros válidos, se identificó que la temperatura presentó valores entre 22.20 y 29.70, con una media de 25.83 y una mediana de 25.80. Esto indica una distribución relativamente equilibrada, sin presencia evidente de valores atípicos extremos que pudieran afectar de forma considerable el análisis espacial.

La exploración inicial mostró que la temperatura no se distribuye de manera completamente homogénea en el área de estudio. Los mapas de puntos evidenciaron sectores con temperaturas más altas y otros con valores más bajos, lo que sugiere la existencia de agrupamientos espaciales. Esta primera lectura permitió justificar la aplicación de herramientas geoestadísticas para evaluar formalmente la autocorrelación espacial.

El análisis exploratorio de tendencia mostró una disminución leve de la temperatura en relación con las coordenadas Este y Norte. Aunque esta tendencia no parece ser suficientemente fuerte para explicar por sí sola toda la variabilidad de la variable, sí indica la presencia de un gradiente espacial suave. Sin embargo, la dispersión observada en los gráficos también mostró que existen variaciones locales importantes, por lo que fue pertinente continuar con el cálculo del semivariograma experimental.

El semivariograma experimental evidenció una estructura clara de autocorrelación espacial. La semivarianza aumentó conforme se incrementó la distancia entre los puntos, lo que indica que las observaciones cercanas tienden a tener temperaturas más similares entre sí que aquellas ubicadas a mayor distancia. Este resultado confirma que la variable Temperature presenta dependencia espacial y que, por tanto, es adecuado aplicar técnicas de interpolación espacial como el Kriging.

Al comparar los modelos teóricos de semivariograma, el modelo exponencial presentó el mejor ajuste, con el menor valor de error SSErr frente a los modelos esférico y gaussiano. Por esta razón, se seleccionó como el modelo más adecuado para representar la estructura espacial de la temperatura y para realizar la predicción espacial mediante Kriging ordinario.

La validación cruzada mostró que el modelo tiene un desempeño aceptable. La media de los residuos fue cercana a cero, lo que indica que no existe un sesgo sistemático importante en las predicciones. Además, el histograma de residuos mostró una distribución concentrada alrededor de cero, y el gráfico de valores observados frente a valores predichos evidenció una relación positiva clara. Estos resultados respaldan el uso del modelo seleccionado para generar la superficie de predicción espacial.

El mapa de Kriging permitió identificar zonas con temperaturas predichas más altas, especialmente hacia algunos bordes del área de estudio, y zonas con temperaturas más bajas en sectores internos y centrales. Este resultado confirma que la temperatura presenta variaciones espaciales relevantes dentro del área analizada, lo cual puede ser útil para reconocer posibles microclimas o diferencias ambientales internas asociadas al cultivo de aguacate.

El mapa de varianza de predicción mostró que la incertidumbre es menor en las zonas con mayor concentración de puntos observados y mayor hacia los bordes de la grilla de interpolación. Esto indica que las predicciones son más confiables dentro del área cubierta por las observaciones y deben interpretarse con mayor cautela en los extremos donde hay menor densidad de datos.

En conclusión, la variable Temperature presenta una estructura espacial definida, autocorrelación espacial y condiciones adecuadas para ser modelada mediante semivariograma y Kriging ordinario. El modelo exponencial fue el más apropiado para representar dicha estructura, y la interpolación permitió construir una superficie continua de temperatura que facilita la interpretación espacial del comportamiento ambiental dentro del área de estudio.