Código
#install.packages(c('tidyverse','sf','terra','stars','tmap','remotes','geodata','rstac','viridis','gt','quarto','tinytex', 'knitr', 'rmarkdown','gstat','gridExtra','spdep','automap'))Taller 3: Análisis de Autocorrelación Espacial
1 Introducción
La autocorrelación espacial mide cómo los valores de una variable están relacionados espacialmente, es decir, cuánto se desvían de un patrón aleatorio debido a su distribución geográfica. Se espera que valores similares, ya sean altos o bajos, tiendan a agruparse espacialmente.
Para verificar si existe autocorrelación espacial significativa en los datos de temperatura utilizados, se aplicaron dos enfoques estadísticos. Primero, se calculó el índice de Moran, tanto en su versión global como local, y posteriormente, se empleó el variograma para evaluar la dependencia espacial entre las distintas estaciones con los datos recopilados.
Finalmente, tras ajustar los variogramas experimentales generados, se consideraron los resultados del taller anterior para desarrollar un modelo de residuos. Este modelo fue sometido a análisis, ajustes e interpretación, tomando en cuenta los meses de trabajo previos en el taller n°2.
2 Instalación de Paquetes y Carga de Librerías
Las siguientes librerías son fundamentales para la manipulación de datos espaciales y la visualización de resultados en mapas.
Código
library(remotes)
library(tidyverse)
library(sf)
library(terra)
library(tmap)
library(geodata)
library(rstac)
library(viridis)
library(gt)
library(quarto)
library(tinytex)
library(knitr)
library(rmarkdown)
library(gstat)
library(stars)
library(gridExtra)
library(spdep)
library(automap)3 1. Preparación de Datos
La preparación de los datos consistio en cargar los datos previamente generados, luego se visualizó la ubicación de las estaciones de medición de temperatura, tanto las reales como las simuladas. Finalmente, se preparó una grilla que sirviera como base para realizar las interpolaciones necesarias.
Código
region <- read_sf("data/procesada/region_biobio.gpkg") |> st_transform(crs =32718)
data_temp <- read_rds('data/procesada/data_estas_temp_con_predictores.rds') # Cargar los datos de temperatura con predictores
preds <- rast('data/procesada/predictores/todos_predictores.tif') # Cargar los datos raster de predictores (MDE, orientación, pendiente, etc.)
# Agregar una columna que clasifica las estaciones en ficticias o reales según el station_id
data_temp <- data_temp |>
mutate(estacion_ficticia = ifelse(station_id >= 1000, "Ficticia", "Real"))
# 1.1 Visualizar la ubicación de las estaciones de temperatura (reales y ficticias)
tmap_mode('view') # Configurar el modo interactivo para visualización en tmap
tm_shape(region) + # Mostrar la región
tm_borders() + # Dibujar los bordes de la región
tm_shape(data_temp) + # Añadir las estaciones de temperatura
tm_dots(col = "estacion_ficticia", size = 0.05,
palette = c("blue", "red"), # Colores para estaciones reales y ficticias
title = "Tipo de Estación") + # Leyenda que indica el tipo de estación
tm_layout(legend.outside = TRUE, # Leyenda colocada fuera del mapa
legend.title.size = 0.8, # Tamaño del título de la leyenda
legend.text.size = 0.6) # Tamaño del texto de la leyendaCódigo
#1.2 Preparar una grilla para la interpolación
bb <- st_bbox(region) # Obtener los límites geográficos de la región
grilla <- rast(xmin = bb$xmin, xmax = bb$xmax, ymin = bb$ymin, ymax = bb$ymax,
res = 500, crs = "EPSG:32718") # Crear una grilla con resolución de 500 metros y sistema de coordenadas UTM
grilla_st <- grilla |> st_as_stars() # Convertir la grilla en formato 'stars'4 2. Índice de Moran Global para analizar la temperatura de los meses trabajados
El Índice de Moran es una herramienta útil para evaluar cómo se distribuye una variable en el espacio, ayudando a identificar si dicha distribución es agrupada, aleatoria o dispersa. En particular, el Índice de Moran Global permite medir la autocorrelación espacial de un conjunto de datos en su totalidad, proporcionando una visión general de los patrones espaciales. En este contexto, se utiliza para analizar si la temperatura en la Región Metropolitana presenta algún tipo de estructura espacial, es decir, si existen tendencias claras o agrupamientos en los valores de temperatura a nivel regional.
Código
# Inicializar listas para almacenar resultados
moran_global_results <- list()
moran_mc_results <- list()
moran_plots <- list()
# Meses a analizar
meses <- c(3, 6, 9, 12) # (Marzo, Junio, Septiembre y Diciembre)
# Realizar el análisis para cada mes
for (mes in meses) {
# Filtrar los datos de temperatura para el mes correspondiente
data_temp_filtered <- data_temp %>% filter(mes == mes)
# Eliminar duplicados espaciales
data_temp_filtered <- data_temp_filtered %>%
distinct(geom, .keep_all = TRUE)
# Verificar la existencia de estaciones duplicadas (mismo punto espacial)
if (any(duplicated(data_temp_filtered$geom))) {
warning(paste("Existen geometrías duplicadas en los datos filtrados para el mes", mes)) # Advertencia si hay duplicados
}
# Calcular la matriz de distancias entre las estaciones
dist_matrix <- as.matrix(st_distance(data_temp_filtered)) # Generar la matriz de distancias entre estaciones
w <- 1 / dist_matrix # Calcular los pesos inversos para las distancias
diag(w) <- 0 # Poner a cero los elementos diagonales (distancia consigo mismo)
w <- units::drop_units(w) # Quitar las unidades de la matriz
w <- mat2listw(w, style = "W") # Convertir a formato de lista de vecinos para análisis espacial y especificar el estilo
# Calcular el índice de Moran Global para detectar autocorrelación espacial en la temperatura
moran_global <- moran(data_temp_filtered$temp, w, n = length(w$neighbours), S0 = Szero(w))
moran_global_results[[as.character(mes)]] <- moran_global # Guardar el resultado del índice de Moran Global
# Aplicar la prueba de Monte Carlo para evaluar la significancia de la autocorrelación espacial
moran_mc_result <- moran.mc(data_temp_filtered$temp, w, nsim = 10000) # Test Montecarlo con 10,000 simulaciones
# Extraer los resultados relevantes del test Montecarlo
moran_mc_summary <- data.frame(
statistic = moran_mc_result$statistic,
p.value = moran_mc_result$p.value,
month = mes
)
moran_mc_results[[as.character(mes)]] <- moran_mc_summary # Guardar los resultados de la prueba
# Graficar el Moran Plot (visualización de la relación espacial)
plot <- moran.plot(data_temp_filtered$temp, w, main = paste("Moran Plot para el mes de", mes), xlab = "Temperatura Promedio", ylab = "Espacio")
moran_plots[[as.character(mes)]] <- plot # Almacenar gráfico
}
Código
# 2.1 Convertir los resultados en una tabla y guardarla en un archivo RDS
moran_global_table <- bind_rows(lapply(moran_global_results, as.data.frame), .id = "mes")
print(moran_global_table) mes I K
1 3 0.3555352 1.605292
2 6 0.3555352 1.605292
3 9 0.3555352 1.605292
4 12 0.3555352 1.605292Código
write_rds(moran_global_table, "data/procesada/taller3/moran_global_results.rds") # Guardar en carpeta procesada
# Convertir los resultados del método de Montecarlo en una tabla y guardarla en un archivo RDS
moran_mc_table <- bind_rows(moran_mc_results) # No es necesario usar lapply aquí, ya que ya es un data.frame
print(moran_mc_table) statistic p.value month
statistic...1 0.3555352 9.999e-05 3
statistic...2 0.3555352 9.999e-05 6
statistic...3 0.3555352 9.999e-05 9
statistic...4 0.3555352 9.999e-05 12Código
write_rds(moran_mc_table, "data/procesada/taller3/moran_mc_results.rds") # Guardar en carpeta procesadaLos resultados del análisis del Índice de Moran Global para los meses 3, 6, 9 y 12 muestran un valor consistente de 0.3555, lo que indica una autocorrelación espacial positiva en la variable analizada, sugiriendo que las estaciones con valores similares tienden a agruparse. Además, los p-values extremadamente bajos (9.999e-05) confirman que esta autocorrelación es altamente significativa. En resumen, existe un patrón espacial claro en los datos, con una distribución no aleatoria a lo largo de los meses evaluados.
5 3. Índice de Moran Local para analizar la temperatura de los meses trabajados
El Índice de Moran Local se diferencia del índice global en que permite identificar patrones espaciales a nivel local dentro de una región, destacando zonas específicas donde existe autocorrelación espacial. A través de este análisis, es posible identificar agrupamientos de valores similares o detectar áreas con valores atípicos, proporcionando un nivel de detalle más granular en la evaluación espacial.
Además, el Índice de Moran Local también permite detectar outliers o valores anómalos. Estos surgen cuando el índice local se aproxima a 0, pero con p-values pequeños, lo que indica que los valores de la estación son significativamente diferentes a los de las estaciones vecinas, señalando posibles anomalías en la distribución espacial de la variable analizada, como las temperaturas.
Código
# Inicializar lista para almacenar resultados del índice de Moran local
local_moran_results <- list()
# Realizar el análisis para cada mes
for (mes in meses) {
# Filtrar los datos de temperatura para el mes correspondiente
data_temp_filtered <- data_temp %>% filter(mes == mes)
# Calcular la matriz de distancias entre las estaciones
dist_matrix <- as.matrix(st_distance(data_temp_filtered))
# Quitar las unidades de la matriz de distancias antes de operar
dist_matrix <- units::drop_units(dist_matrix)
# Reemplazar distancias cero con un valor pequeño para evitar infinitos
dist_matrix[dist_matrix == 0] <- min(dist_matrix[dist_matrix > 0]) / 10
# Calcular los pesos inversos para las distancias
w <- 1 / dist_matrix
# Verificar si hay NA o infinitos en la matriz de pesos
w[is.na(w) | is.infinite(w)] <- 0
diag(w) <- 0 # Poner a cero los elementos diagonales (distancia consigo mismo)
w <- mat2listw(w, style = "W") # Convertir a formato de lista de vecinos para análisis espacial
# Calcular el índice de Moran Local
lm <- localmoran(data_temp_filtered$temp, w) # Cálculo del índice de Moran Local
quad <- attr(lm, 'quadr')
P_val <- lm[, 5]
data_temp_filtered <- data_temp_filtered |>
bind_cols(quad = quad$median) |> # Agregar los cuadrantes a los datos filtrados
mutate(quad = case_when(
P_val < 0.05 ~ quad,
.default = 'No Significativo'
))
# Guardar resultados en la lista
local_moran_results[[as.character(mes)]] <- data_temp_filtered # Guardar resultados para el mes actual
tmap_mode('plot')
local_moran_map <- tm_shape(region) + # Añadir capa de fondo (límites o fondo geográfico)
tm_borders() + # Dibujar las fronteras de la capa de fondo
tm_shape(data_temp_filtered) + # Añadir la capa de puntos
tm_dots(col = "quad", style = "cat",
palette = c('High-High' = "red", 'Low-Low' = "blue", 'High-Low' = "green", 'Low-High' = "yellow", 'No Significativo' = "black"),
title = "Cuadrantes Índice de M.L.",
size = 0.4) + # Aumentar el tamaño de los puntos
tm_layout(title = paste("Índice de Moran Local - Temperatura Promedio -", mes), # Título del mapa
title.size = 2.5, # Aumentar el tamaño del título del mapa
legend.outside = TRUE, # Colocar la leyenda fuera del mapa
legend.title.size = 0.9, # Tamaño del título de la leyenda
legend.text.size = 0.9) # Tamaño del texto de la leyenda
# Mostrar el mapa
print(local_moran_map)
}
Código
# Convertir resultados a un data frame para guardar
local_moran_df <- bind_rows(local_moran_results, .id = "mes") # Asegurarse de que todos los resultados se guardan
# Guardar los resultados en un archivo RDS
write_rds(local_moran_df, "data/procesada/taller3/local_moran_results.rds") # Ajusta la ruta según tu estructura de carpetas6 4.Análisis de autocorrelación espacial mediante el uso de variogramas para todos los meses trabajados
El análisis de autocorrelación espacial mediante variogramas permite evaluar la dependencia espacial de la temperatura para cada mes analizado. A través del variograma, se observan tres componentes clave: el nugget, que refleja la variabilidad a pequeñas distancias; el sill, que indica el punto donde la autocorrelación se estabiliza; y el rango, que señala la distancia máxima a la que los valores de temperatura están correlacionados. Al calcular los variogramas para los meses 3, 6, 9 y 12, se puede determinar si los patrones espaciales de la temperatura son consistentes o varían estacionalmente. Este análisis facilita la identificación de cambios en la estructura espacial a lo largo del año, proporcionando una visión más clara de cómo la temperatura se distribuye en el espacio en distintos momentos.
Código
#4.1 Análisis variográfico para el mes de marzo
# Filtrar los datos para el mes de marzo
data_temp_filtered_marzo <- data_temp[data_temp$mes == 3, ] # Filtrar datos para marzo
# Crear el variograma experimental de las temperaturas
v_marzo <- variogram(temp ~ 1, data_temp_filtered_marzo,
cutoff = 160000, width = 6000) # Variograma experimental
# Graficar el variograma experimental
plot(v_marzo, main = "Variograma Experimental - Mes: Marzo",
xlab = "Distancia", ylab = "Semivarianza") # Visualizar el variograma
Código
# Ajustar el variograma usando autoajuste
auto_fit_marzo <- autofitVariogram(temp ~ 1, data_temp_filtered_marzo) # Autoajuste del variograma
#| results: hide print(auto_fit_marzo) # Mostrar los resultados del ajuste
plot(auto_fit_marzo) # Graficar el modelo ajustado
Código
# Definir manualmente los parámetros iniciales del modelo Gaussiano
nugget_gau_marzo <- 1.1 # Estimación inicial del nugget
sill_gau_marzo <- 110 # Estimación inicial del sill
psill_gau_marzo <- sill_gau_marzo - nugget_gau_marzo # Cálculo del psill
range_gau_marzo <- 160000 # Estimación del rango
# Crear el modelo de variograma Gaussiano
model_vgm_gau_marzo <- vgm(psill = psill_gau_marzo, model = "Gau",range = range_gau_marzo, nugget = nugget_gau_marzo) # Modelo Gaussiano
# Ajustar el modelo al variograma experimental
fit_v_gau_marzo <- fit.variogram(v_marzo, model_vgm_gau_marzo)
print(fit_v_gau_marzo) # Ver resultados del ajuste model psill range
1 Nug 0.06306149 0.00
2 Gau 5.12962270 90220.57Código
# Graficar el variograma ajustado
plot(v_marzo, fit_v_gau_marzo, main = "Variograma Ajustado - Mes: Marzo",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_gau_marzo, "data/procesada/taller3/fit_variograma_gau_marzo.rds")
# Análisis y descripción variográfica anisotrópica para las cuatro direcciones principales del mes de marzo
# Crear el variograma experimental anisotrópico en diferentes direcciones
variograma_anisotropico_marzo <- variogram(temp ~ 1, data_temp_filtered_marzo,
alpha = c(0, 45, 90, 135), # Direcciones a analizar
cutoff = 160000, width = 6000,
tol.hor = 22.5) # Tolerancia horizontal en grados
# Graficar el variograma experimental anisotrópico
plot(variograma_anisotropico_marzo, main = "Variograma Anisotrópico - Mes: Marzo",
xlab = "Distancia", ylab = "Semivarianza") # Graficar variograma anisotrópico
Código
# Definir modelo anisotrópico
model_vgm_aniso_marzo <- vgm(psill = psill_gau_marzo, model = "Gau",
range = range_gau_marzo, nugget = nugget_gau_marzo,
anis = c(45, 0.5)) # Ajuste de anisotropía
# Ajustar el modelo anisotrópico al variograma experimental
fit_v_aniso_marzo <- fit.variogram(variograma_anisotropico_marzo, model_vgm_aniso_marzo)
print(fit_v_aniso_marzo) # Mostrar el modelo ajustado model psill range ang1 anis1
1 Nug 0.05333708 0.0 0 1.0
2 Gau 5.08705316 140967.8 45 0.5Código
# Graficar el variograma ajustado anisotrópico
plot(variograma_anisotropico_marzo, fit_v_aniso_marzo,
main = "Variograma Ajustado Anisotrópico - Mes: Marzo",
xlab = "Distancia", ylab = "Semivarianza") # Graficar
Código
# Guardar el modelo ajustado anisotrópico
saveRDS(fit_v_aniso_marzo, "data/procesada/taller3/fit_variograma_aniso_marzo.rds") # Guardar en un archivo .rds
#4.2 Análisis variográfico para el mes de junio
# Filtrar los datos para el mes de junio
data_temp_filtered_junio <- data_temp[data_temp$mes == 6, ] # Filtrar datos para junio
# Crear el variograma experimental de las temperaturas
v_junio <- variogram(temp ~ 1, data_temp_filtered_junio,
cutoff = 160000, width = 6000) # Variograma experimental
# Graficar el variograma experimental
plot(v_junio, main = "Variograma Experimental - Mes: Junio",
xlab = "Distancia", ylab = "Semivarianza") # Visualizar el variograma
Código
# Ajustar el variograma usando autoajuste
auto_fit_junio <- autofitVariogram(temp ~ 1, data_temp_filtered_junio) # Autoajuste del variograma
#| results: hide print(auto_fit_junio) # Mostrar los resultados del ajuste
plot(auto_fit_junio) # Graficar el modelo ajustado
Código
# Definir manualmente los parámetros iniciales del modelo Ste (Stein)
nugget_ste_junio <- 0.5 # Estimación inicial del nugget
sill_ste_junio <- 200 # Estimación inicial del sill
psill_ste_junio <- sill_ste_junio - nugget_ste_junio # Cálculo del psill
range_ste_junio <- 160000 # Estimación del rango
# Crear el modelo de variograma Stein (Ste)
model_vgm_ste_junio <- vgm(psill = psill_ste_junio, model = "Ste",range = range_ste_junio, nugget = nugget_ste_junio, kappa = 1.3)
# Ajustar el modelo al variograma experimental
fit_v_ste_junio <- fit.variogram(v_junio, model_vgm_ste_junio)
print(fit_v_ste_junio) # Ver resultados del ajuste model psill range kappa
1 Nug 0.00000 0.0 0.0
2 Ste 18.19581 190863.4 1.3Código
# Graficar el variograma ajustado
plot(v_junio, fit_v_ste_junio, main = "Variograma Ajustado - Mes: Junio",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_ste_junio, "data/procesada/taller3/fit_variograma_ste_junio.rds")
# Crear el variograma experimental anisotrópico en diferentes direcciones
variograma_anisotropico_junio <- variogram(temp ~ 1, data_temp_filtered_junio,
alpha = c(0, 45, 90, 135), # Direcciones a analizar
cutoff = 160000, width = 6000,
tol.hor = 22.5) # Tolerancia horizontal en grados
# Graficar el variograma experimental anisotrópico
plot(variograma_anisotropico_junio, main = "Variograma Anisotrópico - Mes: Junio",
xlab = "Distancia", ylab = "Semivarianza") # Graficar variograma anisotrópico
Código
# Definir modelo anisotrópico
model_vgm_aniso_junio <- vgm(psill = psill_ste_junio, model = "Ste",range = range_ste_junio, nugget = nugget_ste_junio, kappa = 1.3, anis = c(45, 0.5)) # Ajuste de anisotropía
# Ajustar el modelo anisotrópico al variograma experimental
fit_v_aniso_junio <- fit.variogram(variograma_anisotropico_junio, model_vgm_aniso_junio)
print(fit_v_aniso_junio) # Mostrar el modelo ajustado model psill range kappa ang1 anis1
1 Nug 0.0000 0.0 0.0 0 1.0
2 Ste 18.1279 300804.4 1.3 45 0.5Código
# Graficar el variograma ajustado anisotrópico
plot(variograma_anisotropico_junio, fit_v_aniso_junio,
main = "Variograma Ajustado Anisotrópico - Mes: Junio",
xlab = "Distancia", ylab = "Semivarianza") # Graficar
Código
# Guardar el modelo ajustado anisotrópico
saveRDS(fit_v_aniso_junio, "data/procesada/taller3/fit_variograma_aniso_junio.rds") # Guardar en un archivo .rds
#4.3Análisis variográfico para el mes de septiembre----
# Filtrar los datos para el mes de septiembre
data_temp_filtered_septiembre <- data_temp[data_temp$mes == 9, ] # Filtrar datos para septiembre
# Crear el variograma experimental de las temperaturas
v_septiembre <- variogram(temp ~ 1, data_temp_filtered_septiembre,
cutoff = 160000, width = 6000) # Variograma experimental
# Graficar el variograma experimental
plot(v_septiembre, main = "Variograma Experimental - Mes: Septiembre",
xlab = "Distancia", ylab = "Semivarianza") # Visualizar el variograma
Código
# Ajustar el variograma usando autoajuste
auto_fit_septiembre <- autofitVariogram(temp ~ 1, data_temp_filtered_septiembre) # Autoajuste del variograma
#| results: hide print(auto_fit_septiembre) # Mostrar los resultados del ajuste
plot(auto_fit_septiembre) # Graficar el modelo ajustado
Código
# Definir manualmente los parámetros iniciales del modelo Ste (Stein)
nugget_ste_septiembre <- 0.43 # Estimación inicial del nugget
sill_ste_septiembre <- 225 # Estimación inicial del sill
psill_ste_septiembre <- sill_ste_septiembre - nugget_ste_septiembre # Cálculo del psill
range_ste_septiembre <- 160000 # Estimación del rango
# Crear el modelo de variograma Stein (Ste)
model_vgm_ste_septiembre <- vgm(psill = psill_ste_septiembre, model = "Ste",range = range_ste_septiembre, nugget = nugget_ste_septiembre, kappa = 1.2)
# Ajustar el modelo al variograma experimental
fit_v_ste_septiembre <- fit.variogram(v_septiembre, model_vgm_ste_septiembre)
print(fit_v_ste_septiembre) # Ver resultados del ajuste model psill range kappa
1 Nug 0.00000 0.0 0.0
2 Ste 16.50387 147648.7 1.2Código
# Graficar el variograma ajustado
plot(v_septiembre, fit_v_ste_septiembre, main = "Variograma Ajustado - Mes: Septiembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_ste_septiembre, "data/procesada/taller3/fit_variograma_ste_septiembre.rds")
# Crear el variograma experimental anisotrópico en diferentes direcciones
variograma_anisotropico_septiembre <- variogram(temp ~ 1, data_temp_filtered_septiembre,
alpha = c(0, 45, 90, 135), # Direcciones a analizar
cutoff = 160000, width = 6000,
tol.hor = 22.5) # Tolerancia horizontal en grados
# Graficar el variograma experimental anisotrópico
plot(variograma_anisotropico_septiembre, main = "Variograma Anisotrópico - Mes: Septiembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar variograma anisotrópico
Código
# Definir modelo anisotrópico
model_vgm_aniso_septiembre <- vgm(psill = psill_ste_septiembre, model = "Ste",range = range_ste_septiembre, nugget = nugget_ste_septiembre, kappa = 1.2, anis = c(45, 0.5)) # Ajuste de anisotropía
# Ajustar el modelo anisotrópico al variograma experimental
fit_v_aniso_septiembre <- fit.variogram(variograma_anisotropico_septiembre, model_vgm_aniso_septiembre)
print(fit_v_aniso_septiembre) # Mostrar el modelo ajustado model psill range kappa ang1 anis1
1 Nug 0.00000 0.0 0.0 0 1.0
2 Ste 16.46434 232864.9 1.2 45 0.5Código
# Graficar el variograma ajustado anisotrópico
plot(variograma_anisotropico_septiembre, fit_v_aniso_septiembre,
main = "Variograma Ajustado Anisotrópico - Mes: Septiembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar
Código
# Guardar el modelo ajustado anisotrópico
saveRDS(fit_v_aniso_septiembre, "data/procesada/taller3/fit_variograma_aniso_septiembre.rds") # Guardar en un archivo .rds
#4.4 Análisis variográfico para el mes de diciembre----
# Filtrar los datos para el mes de diciembre
data_temp_filtered_diciembre <- data_temp[data_temp$mes == 12, ] # Filtrar datos para diciembre
# Crear el variograma experimental de las temperaturas
v_diciembre <- variogram(temp ~ 1, data_temp_filtered_diciembre,
cutoff = 160000, width = 6000) # Variograma experimental
# Graficar el variograma experimental
plot(v_diciembre, main = "Variograma Experimental - Mes: Diciembre",
xlab = "Distancia", ylab = "Semivarianza") # Visualizar el variograma
Código
# Ajustar el variograma usando autoajuste
auto_fit_diciembre <- autofitVariogram(temp ~ 1, data_temp_filtered_diciembre) # Autoajuste del variograma
#| results: hide print(auto_fit_diciembre) # Mostrar los resultados del ajuste
plot(auto_fit_diciembre) # Graficar el modelo ajustado
Código
# Definir manualmente los parámetros iniciales del modelo Gaussiano para diciembre
nugget_gau_diciembre <- 0.58 # Estimación inicial del nugget
sill_gau_diciembre <- 149 # Estimación inicial del sill
psill_gau_diciembre <- sill_gau_diciembre - nugget_gau_diciembre # Cálculo del psill
range_gau_diciembre <- 160000 # Estimación del rango
# Crear el modelo de variograma Gaussiano
model_vgm_gau_diciembre <- vgm(psill = psill_gau_diciembre, model = "Gau",
range = range_gau_diciembre, nugget = nugget_gau_diciembre) # Modelo Gaussiano
# Ajustar el modelo al variograma experimental
fit_v_gau_diciembre <- fit.variogram(v_diciembre, model_vgm_gau_diciembre)
print(fit_v_gau_diciembre) # Ver resultados del ajuste model psill range
1 Nug 0.04923307 0.00
2 Gau 8.76967089 86117.67Código
# Graficar el variograma ajustado
plot(v_diciembre, fit_v_gau_diciembre, main = "Variograma Ajustado - Mes: Diciembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_gau_diciembre, "data/procesada/taller3/fit_variograma_gau_diciembre.rds")
# Crear el variograma experimental anisotrópico en diferentes direcciones
variograma_anisotropico_diciembre <- variogram(temp ~ 1, data_temp_filtered_diciembre,
alpha = c(0, 45, 90, 135), # Direcciones a analizar
cutoff = 160000, width = 6000,
tol.hor = 22.5) # Tolerancia horizontal en grados
# Graficar el variograma experimental anisotrópico
plot(variograma_anisotropico_diciembre, main = "Variograma Anisotrópico - Mes: Diciembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar variograma anisotrópico
Código
# Definir modelo anisotrópico
model_vgm_aniso_diciembre <- vgm(psill = psill_gau_diciembre, model = "Gau",
range = range_gau_diciembre, nugget = nugget_gau_diciembre,
anis = c(45, 0.5)) # Ajuste de anisotropía
# Ajustar el modelo anisotrópico al variograma experimental
fit_v_aniso_diciembre <- fit.variogram(variograma_anisotropico_diciembre, model_vgm_aniso_diciembre)
print(fit_v_aniso_diciembre) # Mostrar el modelo ajustado model psill range ang1 anis1
1 Nug 0.04051726 0.0 0 1.0
2 Gau 8.73214456 135302.7 45 0.5Código
# Graficar el variograma ajustado anisotrópico
plot(variograma_anisotropico_diciembre, fit_v_aniso_diciembre,
main = "Variograma Ajustado Anisotrópico - Mes: Diciembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar
Código
# Guardar el modelo ajustado anisotrópico
saveRDS(fit_v_aniso_diciembre, "data/procesada/taller3/fit_variograma_aniso_diciembre.rds") # Guardar en un archivo .rds7 5.Análisis variográfico de los residuos
El análisis variográfico de los residuos consiste en estudiar la autocorrelación espacial de los errores o residuos de un modelo, verificando si existe algún patrón espacial no capturado por el modelo ajustado. Mediante el cálculo del variograma de los residuos, se puede evaluar si la variabilidad de estos depende de la distancia entre las observaciones. Si el variograma muestra autocorrelación, esto sugiere que el modelo no ha explicado completamente la estructura espacial de los datos. Por el contrario, si el variograma no presenta autocorrelación, los residuos están distribuidos de forma aleatoria, lo que indica un buen ajuste del modelo en términos espaciales.
Código
# Filtrar los datos para el mes de marzo
data_temp_filtered_marzo <- data_temp[data_temp$mes == 3, ] # Filtrar datos para MARZO
# Crear el variograma experimental para los residuos del modelo con predictores
vresiduos_marzo <- variogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_marzo + lst_marzo, data_temp_filtered_marzo, cutoff = 80000, width = 6000)
# Graficar el variograma experimental de los residuos OK
plot(vresiduos_marzo)
Código
# Ajustar el modelo de variograma para los residuos
auto_fit_residuos_marzo <- autofitVariogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_diciembre + lst_diciembre, data_temp_filtered_marzo)
#| results: hide print(auto_fit_residuos_marzo)
plot(auto_fit_residuos_marzo)
Código
# Definir manualmente los parámetros iniciales del modelo Ste (Stein)
nugget_ste_marzo <- 2.18 # Estimación inicial del nugget
sill_ste_marzo <- 6 # Estimación inicial del sill
psill_ste_marzo <- sill_gau_marzo - nugget_gau_marzo # Cálculo del psill
range_ste_marzo <- 70000 # Estimación del rango
# Crear el modelo de variograma Stein (Ste)
model_vgm_ste_marzo <- vgm(psill = psill_ste_marzo, model = "Ste", range = range_ste_marzo, nugget = nugget_ste_marzo, kappa = 0.3)
# Ajustar el modelo al variograma experimental
fit_v_ste_res_marzo <- fit.variogram(vresiduos_marzo, model_vgm_ste_marzo)
print(fit_v_ste_res_marzo) # Ver resultados del ajuste model psill range kappa
1 Nug 0.1089285 0.00 0.0
2 Ste 0.6224526 29879.87 0.3Código
# Graficar el variograma ajustado residuos
plot(vresiduos_marzo, fit_v_ste_res_marzo, main = "Variograma Ajustado Residuos- Mes: Marzo",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_ste_res_marzo, "data/procesada/taller3/fit_variograma_ste_res_marzo.rds")
#5.2 Residuos del mes de junio
# Filtrar los datos para el mes de junio
data_temp_filtered_junio <- data_temp[data_temp$mes == 6, ] # Filtrar datos para JUNIO
# Crear el variograma experimental para los residuos del modelo con predictores
vresiduos_junio <- variogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_junio + lst_junio, data_temp_filtered_junio, cutoff = 80000, width = 6000)
# Graficar el variograma experimental de los residuos
plot(vresiduos_junio)
Código
# Ajustar el modelo de variograma para los residuos
auto_fit_residuos_junio <- autofitVariogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_junio + lst_junio, data_temp_filtered_junio)
#| results: hide print(auto_fit_residuos_junio)
plot(auto_fit_residuos_junio)
Código
# Definir manualmente los parámetros iniciales del modelo Ste (Stein)
nugget_ste_junio <- 1.25 # Estimación inicial del nugget
sill_ste_junio <- 4.5 # Estimación inicial del sill
psill_ste_junio <- sill_ste_junio - nugget_ste_junio # Cálculo del psill
range_ste_junio <- 70000 # Estimación del rango
# Crear el modelo de variograma Stein (Ste)
model_vgm_ste_junio <- vgm(psill = psill_ste_junio, model = "Ste", range = range_ste_junio, nugget = nugget_ste_junio, kappa = 0.3)
# Ajustar el modelo al variograma experimental
fit_v_ste_res_junio <- fit.variogram(vresiduos_junio, model_vgm_ste_junio)
print(fit_v_ste_res_junio) # Ver resultados del ajuste model psill range kappa
1 Nug 0.0000000 0.00 0.0
2 Ste 0.3623495 17517.03 0.3Código
# Graficar el variograma ajustado residuos
plot(vresiduos_junio, fit_v_ste_res_junio, main = "Variograma Ajustado Residuos- Mes: Junio",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_ste_res_junio, "data/procesada/taller3/fit_variograma_ste_res_junio.rds")
#5.3 Residuos del mes de septiembre
# Filtrar los datos para el mes de septiembre
data_temp_filtered_septiembre <- data_temp[data_temp$mes == 9, ] # Filtrar datos para SEPTIEMBRE
# Crear el variograma experimental para los residuos del modelo con predictores
vresiduos_septiembre <- variogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_septiembre + lst_septiembre, data_temp_filtered_septiembre, cutoff = 80000, width = 7000)
# Graficar el variograma experimental de los residuos
plot(vresiduos_septiembre)
Código
# Ajustar el modelo de variograma para los residuos
auto_fit_residuos_septiembre <- autofitVariogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_septiembre + lst_septiembre, data_temp_filtered_septiembre)
#| results: hide print(auto_fit_residuos_septiembre)
plot(auto_fit_residuos_septiembre)
Código
# Definir manualmente los parámetros iniciales del modelo Ste (Stein)
nugget_ste_septiembre <-0.5 # Estimación inicial del nugget
sill_ste_septiembre <- 0.7 # Estimación inicial del sill
psill_ste_septiembre <- sill_ste_septiembre - nugget_ste_septiembre # Cálculo del psill
range_ste_septiembre <- 60000 # Estimación del rango
# Crear el modelo de variograma Stein (Ste)
model_vgm_ste_septiembre <- vgm(psill = psill_ste_septiembre, model = "Ste", range = range_ste_septiembre, nugget = nugget_ste_septiembre, kappa = 0.3)
# Ajustar el modelo al variograma experimental
fit_v_ste_res_septiembre <- fit.variogram(vresiduos_septiembre, model_vgm_ste_septiembre)
print(fit_v_ste_res_septiembre) # Ver resultados del ajuste model psill range kappa
1 Nug 0.4935584 0.0 0.0
2 Ste 0.6222198 567742.5 0.3Código
# Graficar el variograma ajustado residuos
plot(vresiduos_septiembre, fit_v_ste_res_septiembre, main = "Variograma Ajustado Residuos- Mes: Septiembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_ste_res_septiembre, "data/procesada/taller3/fit_variograma_ste_res_septiembre.rds")
#5.4 Residuos del mes de diciembre
# Filtrar los datos para el mes de diciembre
data_temp_filtered_diciembre <- data_temp[data_temp$mes == 12, ] # Filtrar datos para DICIEMBRE
# Crear el variograma experimental para los residuos del modelo con predictores
vresiduos_diciembre <- variogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_diciembre + lst_diciembre, data_temp_filtered_diciembre, cutoff = 80000, width = 6000)
# Graficar el variograma experimental de los residuos
plot(vresiduos_diciembre)
Código
# Ajustar el modelo de variograma para los residuos
auto_fit_residuos_diciembre <- autofitVariogram(temp ~ dem + dem + aspect + slope + dist_costa + ndvi_diciembre + lst_diciembre, data_temp_filtered_diciembre)
#| results: hide print(auto_fit_residuos_diciembre)
plot(auto_fit_residuos_diciembre)
Código
# Definir manualmente los parámetros iniciales del modelo Ste (Stein)
nugget_ste_diciembre <- 0.5 # Estimación inicial del nugget
sill_ste_diciembre <- 2 # Estimación inicial del sill
psill_ste_diciembre <- sill_ste_diciembre - nugget_ste_diciembre # Cálculo del psill
range_ste_diciembre <- 70000 # Estimación del rango
# Crear el modelo de variograma Stein (Ste)
model_vgm_ste_diciembre <- vgm(psill = psill_ste_diciembre, model = "Ste", range = range_ste_diciembre, nugget = nugget_ste_diciembre, kappa = 0.3)
# Ajustar el modelo al variograma experimental
fit_v_ste_res_diciembre <- fit.variogram(vresiduos_diciembre, model_vgm_ste_diciembre)
print(fit_v_ste_res_diciembre) # Ver resultados del ajuste model psill range kappa
1 Nug 0.08663936 0.00 0.0
2 Ste 1.16105555 29159.78 0.3Código
# Graficar el variograma ajustado residuos
plot(vresiduos_diciembre, fit_v_ste_res_diciembre, main = "Variograma Ajustado Residuos- Mes: Diciembre",
xlab = "Distancia", ylab = "Semivarianza") # Graficar el variograma ajustado
Código
# Guardar el modelo ajustado en un archivo .rds
saveRDS(fit_v_ste_res_diciembre, "data/procesada/taller3/fit_variograma_ste_res_diciembre.rds")8 Conclusion
En el taller 3 se llevó a cabo un análisis variográfico de los residuos para los meses de marzo, junio, septiembre y diciembre, con el objetivo de evaluar la autocorrelación espacial en los errores de un modelo previamente ajustado. A través de los variogramas experimentales y ajustados, se observó la relación entre la semivarianza y la distancia entre las observaciones, lo que permite identificar la dependencia espacial no capturada por el modelo.
En los gráficos de los diferentes meses, se nota que el nugget, que refleja la variabilidad no explicada a distancias cortas, varía de un mes a otro. En algunos periodos, el nugget es pequeño, lo que sugiere que la variabilidad a distancias cortas es baja, posiblemente indicando que el modelo está capturando mejor la estructura espacial en esos meses. Sin embargo, en otros meses, como marzo, el nugget es más elevado, lo que indica una mayor variabilidad no explicada, probablemente debido a ruido en los datos o efectos locales no capturados por el modelo.
El sill, que marca el punto donde la autocorrelación espacial deja de ser significativa, también presenta variaciones a lo largo de los meses. En los meses con un sill más elevado, la autocorrelación espacial persiste a mayores distancias, lo que sugiere que el modelo no está capturando correctamente los patrones espaciales en esas áreas. Asimismo, el rango, que indica la distancia máxima hasta la cual las observaciones están correlacionadas espacialmente, difiere entre los meses. En algunos periodos, el rango es más corto, lo que implica que la autocorrelación espacial se limita a áreas más pequeñas, mientras que en otros meses, la correlación se extiende a distancias más amplias.
En conclusión, el análisis variográfico de los residuos revela que la estructura espacial de los errores del modelo no es constante durante todo el año. En meses como marzo y diciembre, se detecta una mayor autocorrelación espacial en los residuos, lo que sugiere que el modelo podría no estar capturando adecuadamente la variabilidad espacial en esos periodos. Por el contrario, en meses como junio, la menor autocorrelación espacial indica que el modelo es más preciso para representar la estructura espacial. Esto sugiere la necesidad de ajustar el modelo de acuerdo con la estacionalidad para mejorar su capacidad predictiva y reducir la autocorrelación en los residuos.