Análisis Multivariable y Predicción Kriging en Producción de Aguacate (Cauca)
Angie Murillo López
2026-06-25
Introducción y Objetivos
Este estudio aplica herramientas de estadística espacial y geoestadística para analizar la variabilidad espacial de variables microclimáticas (focalizado en Temperatura) en un lote productor de aguacate ubicado en el departamento del Cauca.
El objetivo es evaluar la existencia de autocorrelación espacial mediante el semivariograma empírico, ajustar el mejor modelo teórico (Esférico, Exponencial o Gaussiano) y realizar una interpolación espacial mediante Kriging Ordinario para predecir el comportamiento de la variable en zonas no muestreadas de la finca.
Paso 1: Carga de Librerías y Datos
Se procede con la carga de las librerías fundamentales para el manejo de datos espaciales y la lectura de archivos matriciales.
Importación y Filtrado por Fecha
Se realiza la lectura del conjunto de datos completo y se aplica un filtro estricto sobre la columna FORMATTED_DATE_TIME para aislar los registros correspondientes al primer periodo de muestreo (01/10/2020).
ruta_excel <- "C:/Users/lizet/Documents/2026/CIENCIA DE DATOS/ACTIVIDAD 3 AGUACATE/Datos_Completos_Aguacate.xlsx"
datos_completo <- read_excel(ruta_excel)
datos_filtrados <- datos_completo %>%
filter(substr(FORMATTED_DATE_TIME, 1, 10) == "01/10/2020")
cat("Número total de registros para el 01/10/2020:", nrow(datos_filtrados), "\n")## Número total de registros para el 01/10/2020: 534Paso 2: Análisis Exploratorio de Datos (EDA) y Creación del GeoData
Se transforma el marco de datos filtrado en un objeto de clase geodata propio de la librería geoR, asociando las coordenadas geométricas de los árboles con sus respectivos valores térmicos.
datos_geo <- as.geodata(datos_filtrados,
coords.col = c("Longitude", "Latitude"),
data.col = "Temperature")
summary(datos_geo)## Number of data points: 534
##
## Coordinates summary
## Longitude Latitude
## min -76.71180 2.392101
## max -76.71022 2.393634
##
## Distance summary
## min max
## 1.711724e-05 1.959127e-03
##
## Data summary
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 22.20000 24.50000 25.80000 25.82903 27.17500 29.70000Gráfico Exploratorio de Tendencias
Se genera la representación gráfica inicial para inspeccionar visualmente el comportamiento de la variable en función de su distribución geográfica.
El gráfico cuádruple de diagnóstico espacial para la variable
Temperature permite extraer las siguientes observaciones estructurales
sobre el predio:Distribución Espacial (Superior Izquierdo): Muestra la
ubicación geográfica de los 534 árboles evaluados en la finca. La
densidad de muestreo es homogénea y cubre la totalidad de la extensión
del lote, delineando la geometría del cultivo mediante coordenadas de
Longitud (\(X\)) y Latitud (\(Y\)).
Tendencia con la Latitud (Superior Derecho): Al dispersar los valores térmicos (data) frente al eje \(Y\) (Y Coord), se evidencia una distribución aleatoria y dispersa sin pendientes ascendentes o descendentes marcadas, lo que sugiere ausencia de una tendencia norte-sur.
Tendencia con la Longitud (Inferior Izquierdo): Al contrastar la temperatura frente al eje \(X\) (X Coord), los datos se agrupan de forma estacionaria sin un patrón lineal definido, indicando la inexistencia de una tendencia este-oeste.
Comportamiento de la Variable (Inferior Derecho): El histograma y la curva de densidad muestran que la temperatura oscila entre los 22°C y los 30°C, registrando un comportamiento unimodal con una concentración principal (moda) cercana a los 26°C. La distribución no presenta valores atípicos (outliers) extremos, lo que valida la estabilidad de los datos para los modelos de interpolación lineal óptimos.
Conclusión EDA: Al no identificarse tendencias geográficas lineales significativas en ninguno de los dos ejes coordenados, se asume el cumplimiento del supuesto de estacionaridad de primer orden. Por lo tanto, no se requiere remover tendencias (detrending) y se puede proceder directamente al cálculo del semivariograma empírico sobre los valores originales de la variable.
Paso 3: Análisis de Autocorrelación Espacial (Semivariograma Empírico)
Se calcula el semivariograma omnidireccional con el fin de cuantificar la variabilidad espacial y la tasa de disimilitud entre los pares de puntos en función de la distancia.
dist_max <- max(dist(datos_geo$coords)) * 0.6
variograma_emp <- variog(datos_geo, max.dist = dist_max)## variog: computing omnidirectional variogramplot(variograma_emp, main = "Semivariograma Empírico: Temperature",
xlab = "Distancia", ylab = "Semivarianza", pch = 19, col = "darkblue")
El análisis del semivariograma empírico omnidireccional para la variable Temperature permite extraer las siguientes conclusiones sobre la estructura de dependencia espacial en el cultivo:
Estructura de la Curva: Se observa un comportamiento asintótico clásico y bien definido. A distancias cortas (cercanas a \(0e+00\)), la semivarianza es notablemente baja (aproximadamente \(1.4\)), lo que confirma que los árboles de aguacate más cercanos entre sí registran temperaturas muy similares.
Presencia de Dependencia Espacial: A medida que la distancia de separación entre los puntos de muestreo se incrementa en el eje horizontal, la semivarianza asciende de forma continua. Esto demuestra matemáticamente el cumplimiento de la Primera Ley de la Geografía (Tobler): existe una clara autocorrelación espacial donde la similitud térmica disminuye conforme aumenta la distancia geográfica.
Efecto Pepita (Nugget Effect): La curva no nace en el origen del eje \(Y\), sino que intercepta de forma teórica cerca del valor \(1.0\) o \(1.2\). Esto representa el efecto pepita, atribuible a la variabilidad a microescala que el diseño de muestreo no logra capturar o a errores aleatorios intrínsecos de los sensores portátiles.
Tendencia hacia la Meseta (Sill): Hacia el extremo derecho del gráfico, a distancias superiores a \(8e-04\), la curva tiende a estabilizarse y desacelerar su crecimiento, acercándose a una meseta o varianza máxima del sistema situada en torno a un valor de \(3.5\) a \(3.7\).
Conclusión Semivariograma: El gráfico valida la existencia de una estructura geoestadística sólida para la temperatura. La clara tendencia ascendente justifica técnicamente abandonar métodos de interpolación determinísticos tradicionales (como el inverso de la distancia o IDW) y proceder con el ajuste de modelos teóricos continuos para ejecutar una predicción espacial óptima mediante Kriging Ordinario.
Paso 4: Ajuste del Modelo Teórico
Se realiza el ajuste de tres estructuras matemáticas teóricas (Esférico, Exponencial y Gaussiano) mediante el método de mínimos cuadrados ponderados (WLS).
# Cálculo de parámetros iniciales requeridos para la optimización
dist_max <- max(dist(datos_geo$coords)) * 0.6
variograma_emp <- variog(datos_geo, max.dist = dist_max)## variog: computing omnidirectional variogramsill_inicial <- var(datos_geo$data)
rango_inicial <- dist_max / 3
# Ajuste de las estructuras teóricas
mod_esferico <- variofit(variograma_emp, ini.cov.pars = c(sill_inicial, rango_inicial),
cov.model = "spherical")## variofit: covariance model used is spherical
## variofit: weights used: npairs
## variofit: minimisation function used: optimmod_exponencial <- variofit(variograma_emp, ini.cov.pars = c(sill_inicial, rango_inicial),
cov.model = "exponential")## variofit: covariance model used is exponential
## variofit: weights used: npairs
## variofit: minimisation function used: optimmod_gaussiano <- variofit(variograma_emp, ini.cov.pars = c(sill_inicial, rango_inicial),
cov.model = "gaussian")## variofit: covariance model used is gaussian
## variofit: weights used: npairs
## variofit: minimisation function used: optimSelección del Modelo Óptimo Se evalúan numéricamente los modelos ajustados para identificar aquel que presente el menor error residual.
plot(variograma_emp, main = "Comparación de Modelos Teóricos Ajustados",
xlab = "Distancia", ylab = "Semivarianza")
lines(mod_esferico, col = "blue", lwd = 2)
lines(mod_exponencial, col = "red", lwd = 2)
lines(mod_gaussiano, col = "darkgreen", lwd = 2)
legend("bottomright", legend = c("Esférico", "Exponencial", "Gaussiano"),
col = c("blue", "red", "darkgreen"), lwd = 2)
## Error de ajuste (Esférico): 9967.541## Error de ajuste (Exponencial): 6960.977## Error de ajuste (Gaussiano): 20707.56Modelo Gaussiano (Línea Verde): Presenta un ascenso extremadamente vertical y prematuro a distancias muy cortas, estabilizándose de forma inmediata en una meseta horizontal. Este comportamiento falla en capturar la progresión gradual de los datos reales (puntos vacíos), por lo que se descarta al no representar adecuadamente la variabilidad espacial del predio.
Modelo Esférico (Línea Azul): Describe un ajuste lineal e inicial aceptable, pero genera una inflexión matemática muy rígida antes de alcanzar su rango, subestimando la semivarianza en las distancias intermedias (alrededor de \(2e-04\)).
Modelo Exponencial (Línea Roja): Es la estructura matemática que mejor se adapta a la naturaleza del fenómeno microclimático. Muestra una transición suave y asintótica que acompaña la tendencia de los puntos observados a lo largo de toda la distancia de muestreo, minimizando de forma óptima la distancia residual respecto a los datos reales.
Conclusión del Ajuste: Basado en el criterio de máxima proximidad visual y menor suma de cuadrados de los residuos ponderados, el Modelo Exponencial es seleccionado como el modelo óptimo. Sus parámetros estructurales de efecto pepita (Nugget) cercano a \(0.1\), meseta (Sill) en torno a \(3.2\) y rango práctico serán los encargados de guiar con la mayor fidelidad estadística la predicción espacial continua mediante Kriging Ordinario.
Paso 5: Predicción Espacial por Kriging Ordinario
Utilizando los parámetros del modelo teórico óptimo seleccionado, se define una rejilla regular de alta resolución para estimar de forma continua los valores de la variable de interés en las zonas no muestreadas del predio.
grid_prediccion <- expand.grid(
Longitude = seq(min(datos_geo$coords[,1]), max(datos_geo$coords[,1]), length.out = 100),
Latitude = seq(min(datos_geo$coords[,2]), max(datos_geo$coords[,2]), length.out = 100)
)
kriging_res <- krige.conv(datos_geo, locations = grid_prediccion,
krige = krige.control(obj.model = mejor_modelo))## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhoodPaso 6: Visualización y Generación de Mapas Finales
Se estructuran las matrices de salida resultantes del Kriging para proceder con la cartografía analítica del lote.
df_resultados <- data.frame(
Longitude = grid_prediccion$Longitude,
Latitude = grid_prediccion$Latitude,
Temperature = kriging_res$predict,
Varianza = kriging_res$krige.var
)Superficie de Predicción Espacial El siguiente mapa describe la distribución continua predicha para la temperatura a lo largo del predio de producción.
ggplot() +
geom_tile(data = df_resultados, aes(x = Longitude, y = Latitude, fill = Temperature)) +
scale_fill_viridis_c(option = "plasma", name = "Temp (°C)") +
geom_point(data = datos_filtrados, aes(x = Longitude, y = Latitude),
color = "black", size = 0.6, alpha = 0.5) +
labs(title = "Superficie de Predicción Espacial: Temperatura",
subtitle = "Kriging Ordinario - Lote de Aguacate (01/10/2020)",
x = "Longitude", y = "Latitude") +
theme_minimal()
El mapa de interpolación generado mediante Kriging Ordinario basado en
el modelo teórico Exponencial permite realizar un análisis detallado de
la distribución térmica continua en el lote de aguacate:
Patrones de Distribución Térmica: Se evidencia una marcada heterogeneidad espacial de la temperatura a lo largo del predio, con valores que oscilan entre los 23°C y los 29°C. La presencia de franjas alternas de temperaturas altas y bajas sugiere microclimas condicionados posiblemente por factores topográficos locales (como la altitud o la exposición solar) característicos de las zonas productoras del Cauca.
Identificación de Zonas Cálidas: En el sector sudoccidental (coordenadas inferiores izquierdas cercanas a Longitud \(-76.7115\)), se delimita un núcleo térmico homogéneo con las temperaturas más elevadas del lote, alcanzando rangos entre 28°C y 29°C.
Identificación de Franjas Frías: Contrastando con el punto anterior, se desarrollan corredores lineales bien definidos con tendencias frías (valores entre 23°C y 24.5°C) que atraviesan el predio en sentido sudoeste-noreste, configurando zonas de acumulación térmica diferencial.
Zonificación para Agricultura de Precisión: La delimitación precisa de estas fronteras microclimáticas proporciona un insumo analítico fundamental para el manejo agronómico del cultivo de aguacate. Permite planificar estrategias tecnificadas diferenciadas, tales como el diseño de sistemas de riego focalizados para las zonas con mayor estrés térmico o el monitoreo fitosanitario específico en los corredores más fríos y húmedos.
Variabilidad de la Varianza (Incertidumbre) Se mapea el error estándar de la predicción, identificando la ganancia de incertidumbre asociada a la distancia respecto a los puntos reales de control.
ggplot(df_resultados, aes(x = Longitude, y = Latitude, fill = Varianza)) +
geom_tile() +
scale_fill_viridis_c(option = "mako", name = "Varianza") +
labs(title = "Mapa de Varianza del Kriging (Incertidumbre)",
subtitle = "Densidad y confiabilidad del modelo espacial estimado",
x = "Longitude", y = "Latitude") +
theme_minimal()
El análisis del mapa de varianza de la predicción permite evaluar la precisión y confiabilidad del modelo geoestadístico estimado a lo largo de la extensión del predio:
Comportamiento en la Zona Muestreada (Áreas Oscuras): Al interior del lote, donde se localiza la densidad de los 534 árboles evaluados (puntos negros), la varianza registra sus valores mínimos (cercanos a \(1.0\)). Esto demuestra que la incertidumbre de la predicción es sumamente baja en el corazón del cultivo, garantizando estimaciones térmicas de alta fidelidad y precisión debido a la proximidad de los datos reales.
Efecto de Borde y Extrapolación (Áreas Claras/Verdes): Hacia la periferia y las zonas externas del predio (especialmente en los márgenes occidental y oriental), la varianza se incrementa de forma exponencial alcanzando valores superiores a \(3.0\). Esto responde matemáticamente a la ausencia de puntos de control cercanos, lo que limita la capacidad del algoritmo y convierte la interpolación en una extrapolación con mayor margen de error.
Geometría del Error Estándar: La distribución de la varianza confirma que el error en el Kriging Ordinario depende exclusivamente de la configuración espacial y geométrica del muestreo, y no del valor numérico de la temperatura en sí. Las pequeñas “coronas” concéntricas de incremento de varianza entre las hileras de árboles evidencian cómo el modelo penaliza la pérdida de continuidad espacial a medida que la grilla predictiva se aleja de un punto real.
Conclusión de la Incertidumbre: El mapa de varianza valida la solidez del diseño experimental al interior de la finca, demostrando que las decisiones agronómicas basadas en la zonificación térmica son altamente confiables dentro de los límites del cultivo. Asimismo, establece una advertencia técnica para omitir las interpretaciones en los bordes exteriores del mapa (zonas verdes claras), donde la falta de datos degrada la capacidad predictiva del modelo.
Conclusiones
Se procesó y depuró con éxito el set de datos microclimáticos, logrando el aislamiento exacto del primer periodo bajo estudio con sus correspondientes registros estructurados.
La implementación del semivariograma evidenció la presencia de dependencia espacial, lo que validó la aplicación de modelos continuos frente a interpolaciones tradicionales no estadísticas.
Los mapas generados ofrecen un insumo analítico clave para la agricultura de precisión en El Cauca, permitiendo zonificar el comportamiento térmico del cultivo de aguacate.