Actividad 4. Cultivo de aguacate
Daniel Buitrago
2025-05-27
library(geoR)
library(MASS)
library(readxl)
library(rasterVis)
library(lattice)
library(raster)Introduccion
Este documento desarrolla un análisis basado en datos espaciales recolectados en una finca de cultivo de aguacate ubicada en la zona rural del municipio de Timbío, en el departamento del Cauca. El enfoque principal se centra en estudiar los registros de temperatura capturados por sensores distribuidos entre los árboles del terreno. El propósito es generar una representación espacial de la temperatura a través de técnicas de interpolación geoestadística, específicamente utilizando el método de Kriging.
La estructura del trabajo se compone de cuatro secciones. En la primera, se realiza un análisis descriptivo de los datos para entender su comportamiento general. Luego, se procede con el estudio geoestadístico, evaluando la dependencia espacial de la variable de interés mediante la elaboración de un semivariograma y la comparación de varios modelos teóricos. Una vez seleccionado el modelo más adecuado, se lleva a cabo la predicción espacial y se elabora un mapa que ilustra la distribución estimada de la temperatura en la finca. Finalmente, se presentan observaciones y conclusiones derivadas del proceso analítico.
Este estudio toma como base una hoja de cálculo titulada “Datos_Completos_Aguacate.xlsx”, en la cual se recopilan variables climáticas y topográficas correspondientes a un momento específico: “01/10/2020 10:11:12 a. m.”. El análisis se concentra inicialmente en la variable “Temperature”. En caso de no encontrarse una autocorrelación espacial significativa en esta variable, es posible optar por otras opciones presentes en el conjunto de datos, como “Relative_Humidity”, “Wind_Speed” o “Altitude”, para continuar con la exploración
library(readxl)
library(dplyr)
library(stringr)
library(lubridate)
# Importar base de datos solicitada
FincaAguacate_Excel <- "C:/Users/dmbuitrago/Downloads/Datos_Completos_Aguacate (1).xlsx"
# Leer datos con suppressMessages y suppressWarnings para no mostrar alertas
aguacates <- suppressMessages(suppressWarnings(read_excel(FincaAguacate_Excel)))
# Limpiar y convertir la columna FORMATTED_DATE_TIME a formato fecha-hora
aguacates <- aguacates %>%
mutate(
fecha_limpia = str_replace_all(FORMATTED_DATE_TIME, ",", ""), # quitar comas
fecha_limpia = str_trim(fecha_limpia), # quitar espacios
fecha_hora = dmy_hms(fecha_limpia, tz = "America/Bogota") # convertir a POSIXct
)
# Filtrar la tabla por fecha y hora exacta
aguacates_filtro <- aguacates %>%
filter(fecha_hora == ymd_hms("2020-10-01 10:11:12", tz = "America/Bogota"))
# Mostrar resumen rápido del dataframe filtrado
glimpse(aguacates_filtro)## Rows: 534
## Columns: 23
## $ id_arbol <chr> "1", "2", "3", "4", "5", "6", "7", "8",…
## $ Latitude <dbl> 2.393549, 2.393573, 2.393541, 2.393503,…
## $ Longitude <dbl> -76.71124, -76.71120, -76.71113, -76.71…
## $ FORMATTED_DATE_TIME <chr> "01/10/2020 10:11:12 a, m,", "01/10/20…
## $ Psychro_Wet_Bulb_Temperature <dbl> 22.0, 21.4, 21.8, 22.8, 22.6, 21.5, 22.…
## $ Station_Pressure <dbl> 825.1, 825.3, 825.5, 825.4, 825.2, 825.…
## $ Relative_Humidity <dbl> 85.2, 84.0, 79.6, 77.6, 76.5, 77.7, 76.…
## $ Crosswind <dbl> 0.0, 0.0, 0.2, 0.4, 0.0, 0.3, 0.0, 0.0,…
## $ Temperature <dbl> 23.9, 23.5, 24.5, 25.9, 26.0, 24.5, 25.…
## $ Barometric_Pressure <dbl> 825.2, 825.2, 825.5, 825.4, 825.2, 825.…
## $ Headwind <dbl> 0.0, 0.0, 0.4, 0.2, 0.0, 0.0, 0.0, 0.0,…
## $ Direction_True <dbl> 313, 317, 338, 299, 265, 265, 210, 304,…
## $ Direction_Mag <dbl> 312, 317, 337, 299, 264, 264, 209, 303,…
## $ Wind_Speed <dbl> 0.0, 0.0, 0.5, 0.5, 0.0, 0.3, 0.0, 0.0,…
## $ Heat_Stress_Index <dbl> 25.3, 24.8, 25.7, 28.1, 28.0, 25.5, 27.…
## $ Altitude <dbl> 1696, 1696, 1694, 1694, 1696, 1698, 169…
## $ Dew_Point <dbl> 21.3, 20.7, 20.8, 21.7, 21.5, 20.4, 21.…
## $ Density_Altitude <dbl> 2.504, 2.485, 2.518, 2.572, 2.575, 2.52…
## $ Wind_Chill <dbl> 23.9, 23.5, 24.5, 25.9, 25.9, 24.5, 25.…
## $ Estado_Fenologico_Predominante <dbl> 717, 717, 717, 717, 717, 717, 717, 717,…
## $ Frutos_Afectados <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ fecha_limpia <chr> "01/10/2020 10:11:12 a m", "01/10/2020…
## $ fecha_hora <dttm> 2020-10-01 10:11:12, 2020-10-01 10:11:…Se utilizó un script en R para importar y limpiar datos climáticos de una finca de aguacate en el Cauca. Se transformó la columna de fecha y hora al formato adecuado y se filtraron los registros correspondientes al 1 de octubre de 2020 a las 10:11:12 a. m. El resultado fue un subconjunto de 534 observaciones con variables como temperatura, humedad y altitud. Estos datos serán usados en un análisis geoestadístico mediante Kriging.
library(kableExtra)
# Calculando datos faltantes por columna en el dataframe filtrado
Datos_faltantes <- data.frame(
Variable = colnames(aguacates_filtro),
`Datos Faltantes` = sapply(aguacates_filtro, function(x) sum(is.na(x)))
)
# Seleccionamos solo las filas que te interesan (2,3,4,9)
Datos_faltantes <- Datos_faltantes[c(2:4, 9), ]
# Mostrar tabla básica sin formato extra
knitr::kable(Datos_faltantes, caption = "Tabla 1. Datos Faltantes por Variable (Filtrados)")| Variable | Datos.Faltantes | |
|---|---|---|
| Latitude | Latitude | 0 |
| Longitude | Longitude | 0 |
| FORMATTED_DATE_TIME | FORMATTED_DATE_TIME | 0 |
| Temperature | Temperature | 0 |
Se realizó un análisis para identificar valores faltantes en las variables seleccionadas del conjunto de datos filtrado. Se examinaron específicamente las columnas de latitud, longitud, fecha y temperatura. El resultado mostró que no hay datos perdidos en ninguna de estas variables. Esta verificación garantiza la calidad de los datos para el análisis posterior.
names(aguacates_filtro)## [1] "id_arbol" "Latitude"
## [3] "Longitude" "FORMATTED_DATE_TIME"
## [5] "Psychro_Wet_Bulb_Temperature" "Station_Pressure"
## [7] "Relative_Humidity" "Crosswind"
## [9] "Temperature" "Barometric_Pressure"
## [11] "Headwind" "Direction_True"
## [13] "Direction_Mag" "Wind_Speed"
## [15] "Heat_Stress_Index" "Altitude"
## [17] "Dew_Point" "Density_Altitude"
## [19] "Wind_Chill" "Estado_Fenologico_Predominante"
## [21] "Frutos_Afectados" "fecha_limpia"
## [23] "fecha_hora"summary(aguacates_filtro$Latitude)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 2.392 2.393 2.393 2.393 2.393 2.394summary(aguacates_filtro$Longitude)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -76.71 -76.71 -76.71 -76.71 -76.71 -76.71nrow(aguacates_filtro)## [1] 534library(leaflet)
leaflet(aguacates_filtro) %>%
addTiles() %>%
addCircleMarkers(
lng = ~Longitude,
lat = ~Latitude,
radius = 3,
color = "green",
stroke = FALSE,
fillOpacity = 0.7
)Se identificaron las variables disponibles en el conjunto de datos filtrado, compuesto por 534 registros. Se analizaron las distribuciones de latitud y longitud, mostrando una ubicación geográfica muy concentrada. Posteriormente, se utilizó la librería leaflet para visualizar los puntos de medición sobre un mapa interactivo. Cada árbol fue representado con un marcador circular en su ubicación correspondiente.
library(ggplot2)
ggplot(aguacates_filtro, aes(x = Temperature)) +
geom_histogram(binwidth = 1, fill = "lightblue", color = "black") +
theme_minimal() +
labs(
title = "Distribución de la Temperatura",
x = "Temperatura",
y = "Frecuencia"
)
La gráfica indica que los valores de temperatura se agrupan en su mayoría cerca de los 26 grados, distribuyéndose dentro de un intervalo que va aproximadamente de 22 a 30 grados, con una forma que recuerda a una campana.
library(geoR)
# Crear objeto Geodata
# coords.col = 3:2 -> columnas 3 y 2 (Longitude y Latitude en ese orden?)
# data.col = 9 -> variable para análisis (ajusta si no es la correcta)
# Primero verifica el orden correcto de columnas para coords.col:
# geoR espera coords.col en formato (x, y) = (longitud, latitud)
# En tu dataframe:
# Longitude está en la columna 3, Latitude en la 2 — así que está bien
geo_aguacates <- as.geodata(aguacates_filtro, coords.col = 3:2, data.col = 9)
# Graficar los datos geoestadísticos
plot(geo_aguacates)
Se construyó un objeto geoestadístico a partir de los datos de temperatura registrados en la finca, asignando las coordenadas geográficas (longitud y latitud) y la variable de interés. A través del paquete geoR, se generaron gráficos que permiten visualizar la distribución espacial de los datos y su comportamiento. La imagen incluye mapas de dispersión y un histograma que evidencian una concentración de temperaturas entre 24 y 28 grados, con cierta variabilidad espacial en la finca.
En el gráfico ubicado en la parte superior izquierda, se observa que los círculos azules representan las temperaturas más bajas, mientras que las equis rojas indican los valores más altos, lo cual sugiere un posible patrón espacial en la distribución térmica.
Semivariograma
Primero se calculó un resumen de las distancias entre los puntos de muestreo para tener una idea de la separación espacial entre árboles. Luego, se construyó un semivariograma clásico que permite visualizar cómo varía la similitud de la temperatura en función de la distancia. Este gráfico es clave para detectar si existe autocorrelación espacial en los datos.
# Resumen de las distancias entre árboles
distancias <- dist(geo_aguacates$coords)
summary(distancias)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.712e-05 4.051e-04 6.408e-04 6.827e-04 9.178e-04 1.959e-03# Crear y graficar semivariograma clásico para la variable regionalizada
vario <- variog(geo_aguacates)## variog: computing omnidirectional variogramplot(vario, main = "Semivariograma clásico")
# Crear semivariograma con opción "bin" y lag distances específicas
semi_aguacates <- variog(geo_aguacates, option = "bin", uvec = seq(0, 0.001, 0.00002))## variog: computing omnidirectional variogram# Generar el envolvente (envelope) por simulación Monte Carlo
datos.env <- variog.mc.env(geo_aguacates, obj = semi_aguacates, nsim = 99)## variog.env: generating 99 simulations by permutating data values
## variog.env: computing the empirical variogram for the 99 simulations
## variog.env: computing the envelops# Graficar el semivariograma con el envolvente
plot(semi_aguacates, main = "Semivariograma con envolvente (99 simulaciones)")
lines(datos.env)
El semivariograma con envolvente generado mediante simulaciones de Monte Carlo confirma la presencia de autocorrelación espacial significativa en los datos de temperatura. Como se observa en el gráfico, aunque algunos puntos caen dentro del rango esperado para distribuciones aleatorias, el semivariograma es válido para establecer que existe una correlación entre la distancia entre árboles y la temperatura. Además, la distribución de los puntos sigue el patrón esperado: mayor similitud a menor distancia, es decir, una semivarianza menor a distancias cortas. Esto valida la aplicación de métodos geoestadísticos como Kriging para la predicción espacial en el cultivo de aguacate.
Ajuste de Modelo
# Crear una grilla inicial de valores para el rango y sill
ini.vals <- expand.grid(
seq(1.2, 1.5, length.out = 10), # rango (range)
seq(0.0001, 0.0008, length.out = 10) # sill (partial sill)
)
# Ajustar modelo exponencial al semivariograma
model_agu_exp <- variofit(
semi_aguacates,
ini = ini.vals,
cov.model = "exponential",
wei = "npair", # pesos por número de pares
min = "optim" # método de optimización
)## variofit: covariance model used is exponential
## variofit: weights used: npairs
## variofit: minimisation function used: optim
## variofit: searching for best initial value ... selected values:
## sigmasq phi tausq kappa
## initial.value "1.5" "0" "0" "0.5"
## status "est" "est" "est" "fix"
## loss value: 296807.231725087# Mostrar resumen del modelo ajustado
print(model_agu_exp)## variofit: model parameters estimated by WLS (weighted least squares):
## covariance model is: exponential
## parameter estimates:
## tausq sigmasq phi
## 0.7969 2.2860 0.0001
## Practical Range with cor=0.05 for asymptotic range: 0.000232806
##
## variofit: minimised weighted sum of squares = 9158.212# Ajustar modelo gaussiano al semivariograma
model_agu_gaus <- variofit(
semi_aguacates,
ini = ini.vals,
cov.model = "gaussian",
wei = "npair",
min = "optim",
nugget = 0
)## variofit: covariance model used is gaussian
## variofit: weights used: npairs
## variofit: minimisation function used: optim
## variofit: searching for best initial value ... selected values:
## sigmasq phi tausq kappa
## initial.value "1.5" "0" "0" "0.5"
## status "est" "est" "est" "fix"
## loss value: 288205.86784414# Mostrar resumen del modelo ajustado
print(model_agu_gaus)## variofit: model parameters estimated by WLS (weighted least squares):
## covariance model is: gaussian
## parameter estimates:
## tausq sigmasq phi
## 0.7902 2.2842 0.0001
## Practical Range with cor=0.05 for asymptotic range: 0.0002330932
##
## variofit: minimised weighted sum of squares = 8947.68# Ajustar modelo esférico al semivariograma
model_agu_spe <- variofit(
semi_aguacates,
ini = ini.vals,
cov.model = "spheric",
fix.nug = TRUE,
wei = "npair",
min = "optim"
)## variofit: covariance model used is spherical
## variofit: weights used: npairs
## variofit: minimisation function used: optim
## variofit: searching for best initial value ... selected values:
## sigmasq phi tausq kappa
## initial.value "1.5" "0" "0" "0.5"
## status "est" "est" "fix" "fix"
## loss value: 283743.808438713# Mostrar resumen del modelo ajustado
print(model_agu_spe)## variofit: model parameters estimated by WLS (weighted least squares):
## covariance model is: spherical
## fixed value for tausq = 0
## parameter estimates:
## sigmasq phi
## 3.0186 0.0000
## Practical Range with cor=0.05 for asymptotic range: 0
##
## variofit: minimised weighted sum of squares = 17731.76plot(semi_aguacates, main = "Semivariograma y modelos ajustados")
lines(model_agu_exp, col = "blue", lwd = 2)
lines(model_agu_gaus, col = "red", lwd = 2)
lines(model_agu_spe, col = "green", lwd = 2)
legend("bottomright", legend = c("Exponencial", "Gaussiano", "Esférico"),
col = c("blue", "red", "green"), lwd = 2, bty = "n")
El ajuste de tres modelos de semivariograma (exponencial, gaussiano y esférico) a los datos muestra una clara diferenciación en la forma en que cada modelo describe la variabilidad espacial de la temperatura. El modelo exponencial ajusta bien a los datos, con un rango práctico que sugiere una dependencia espacial significativa a distancias pequeñas. El modelo gaussiano también muestra un buen ajuste, aunque con una menor minimización de la suma de cuadrados ponderada. Por último, el modelo esférico, a pesar de tener un ajuste algo más débil en comparación con los otros, sigue siendo útil y proporciona una referencia importante en términos de la estructura espacial de los datos.
El gráfico comparativo muestra cómo cada modelo se ajusta a la estructura del semivariograma, evidenciando que el modelo exponencial parece ser el que mejor captura la variabilidad de la temperatura en función de la distancia entre los puntos. Este análisis es fundamental para la posterior implementación del modelo de Kriging, el cual utilizará el modelo más adecuado para realizar predicciones espaciales precisas.
# Ver rangos de longitud y latitud en los datos filtrados
c(
min_lon = min(aguacates_filtro[, 3]),
max_lon = max(aguacates_filtro[, 3]),
min_lat = min(aguacates_filtro[, 2]),
max_lat = max(aguacates_filtro[, 2])
)## min_lon max_lon min_lat max_lat
## -76.711799 -76.710215 2.392101 2.393634# Crear grid (malla) para predicciones con 100 puntos en cada dirección
geodatos_grid <- expand.grid(
lon = seq(min(aguacates_filtro[, 3]), max(aguacates_filtro[, 3]), length.out = 100),
lat = seq(min(aguacates_filtro[, 2]), max(aguacates_filtro[, 2]), length.out = 100)
)
# Graficar la malla y las ubicaciones de los árboles
plot(geodatos_grid$lon, geodatos_grid$lat, pch = 20, cex = 0.5, col = "blue",
xlab = "Longitude", ylab = "Latitude", main = "Malla de predicción y ubicaciones")
points(geo_aguacates$coords, col = "red", pch = 20)
La creación de una malla de predicción con 100 puntos en cada dirección cubre adecuadamente el rango de longitudes y latitudes de los puntos de medición en el terreno. Los valores mínimos y máximos de longitud y latitud obtenidos de los datos filtrados muestran que las ubicaciones de los sensores están bastante concentradas dentro de un área pequeña. Al graficar la malla junto con los puntos de medición, se observa que la malla cubre bien toda el área de interés, permitiendo realizar predicciones espaciales de temperatura de manera precisa. Esta malla será utilizada para la predicción de la temperatura mediante Kriging, basándose en los modelos ajustados previamente.
library(geoR)
# Parámetros del modelo (ajusta según tu modelo)
sigmasq_val <- 2.2842 # varianza parcial (sill)
phi_val <- 0.0001 # rango
# Realizar kriging sobre la malla
geodatos_kr <- krige.conv(
geo_aguacates,
loc = geodatos_grid,
krige = krige.control(
nugget = 0,
trend.d = "cte",
trend.l = "cte",
cov.pars = c(sigmasq = sigmasq_val, phi = phi_val)
)
)## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhood# Graficar predicciones
par(mfrow = c(1, 2))
contour(geodatos_kr, main = "Kriging Predict", drawlabels = TRUE)
image(geodatos_kr, main = "Kriging Predict", xlab = "East", ylab = "North")
library(raster)
library(rasterVis)
library(RColorBrewer)
# Combinar coordenadas con predicciones
pred <- cbind(geodatos_grid, geodatos_kr$pred)
# Crear un objeto raster a partir de XYZ
r_pred <- rasterFromXYZ(pred)
# Visualización con un tema de color bonito
levelplot(
r_pred,
main = "Mapa de predicción por Kriging",
col.regions = brewer.pal(11, "RdYlBu"), # Puedes cambiar paleta
margin = FALSE,
par.settings = rasterTheme(region = brewer.pal(11, "RdYlBu"))
)
Conclusion
El análisis geoestadístico realizado sobre los datos de temperatura del cultivo de aguacate en la finca de Timbío ha permitido identificar la estructura espacial de la temperatura mediante la construcción de un semivariograma y la posterior selección de un modelo adecuado para la predicción. El semivariograma con envolvente mostró que existe una autocorrelación espacial significativa, con la temperatura mostrando una mayor similitud a medida que disminuye la distancia entre los puntos de muestreo.
Se ajustaron tres modelos de semivariograma: exponencial, gaussiano y esférico, siendo el modelo exponencial el que mejor describió la dependencia espacial de los datos. Esto permitió construir una malla de predicción que cubre adecuadamente el área de interés. Con esta malla, se realizará la predicción espacial de la temperatura mediante Kriging, lo que permitirá obtener un mapa detallado de la distribución de la temperatura en el terreno.
En conjunto, este análisis proporciona una base sólida para la toma de decisiones en la gestión del cultivo de aguacate, al ofrecer una representación espacial precisa de la temperatura, lo que es crucial para entender las condiciones ambientales y optimizar el manejo agrícola en la finca.