Datos de Temperatura en Antioquia

Autores/as

Daniel Villa

Raquel Latorre

Ivan Rojas

Juan Esteban Restrepo

Introducción

Antioquia es un departamento con diversidad geográfica y climática, que se caracteriza por tener variaciones de temperatura dependiendo su subregión. Durante los meses de febrero y marzo y agosto- septiembre el departamento se caracteriza por presentar altas temperaturas a lo largo del año, en subregiones como Urabá, Magdalena Medio, Bajo cauca, y nordeste en los límites con el departamento de Bolívar con temperaturas promedio que oscila entre los 32 °C y los 39°C. En esa misma temporada climática también se experimentan bajas temperaturas, que se han visto en los municipios de la subregión norte del departamento, en donde la temperatura ha llegado a bajar entre 1° y 3.5°. Estos fenómenos climáticos han generado una serie de impactos en la salud, el medio ambiente y la calidad de vida de las personas.

Marco teórico.

El departamento de Antioquia se encuentra ubicado en una región geográfica de Colombia conocida como la Cordillera de los Andes, y se caracteriza por presentar una gran diversidad climática y topográfica. En términos generales, la región presenta temperaturas cálidas durante la mayor parte del año, con una media anual de 24°C. Sin embargo, en los últimos años se ha observado un aumento en la frecuencia y la intensidad de las olas de calor en Antioquia, lo que ha generado preocupación entre los habitantes y las autoridades locales. Estas altas temperaturas pueden tener efectos negativos sobre la salud humana, la agricultura, la ganadería, la biodiversidad y otros aspectos clave del desarrollo económico y social de la región. Desde el punto de vista climático, este fenómeno puede estar relacionado con factores como el cambio climático, la variabilidad natural del clima y la interacción entre la atmósfera y la topografía local. En este sentido, diversos estudios científicos han analizado la evolución de las temperaturas en Antioquia y sus posibles causas. Estudios realizados por la Universidad de Antioquia y la gobernación de Antioquia se evidenció un aumento de la temperatura promedio en todo el departamento, con aumentos máximos cercanos a los 3 °C y mínimos de 0.5 °C. Las subregiones de Urabá, Bajo Cauca y Magdalena Medio es donde se espera aumente en mayor medida la temperatura media (3 °C), destacando que actualmente es en estas subregiones donde se presenta la mayor temperatura media. Para el Valle de Aburrá, se espera aumentos cercanos entre 2.1 °C y 2.2 °C, lo que sugiere un patrón de calentamiento a largo plazo. Otros estudios han señalado la En términos de las implicaciones para la sociedad, es importante destacar que las altas temperaturas pueden tener efectos negativos sobre la salud humana, especialmente en grupos vulnerables como los niños, los ancianos y las personas con enfermedades crónicas. Sin embargo, además de estos factores, también existe otro que es particularmente preocupante en Antioquia, la minería ilegal. Esta actividad se ha convertido en una práctica recurrente en la región, lo que ha generado un impacto negativo tanto en el medio ambiente como en la salud de las comunidades locales. La minería ilegal implica la extracción de minerales sin cumplir con las normas ambientales y sociales necesarias para garantizar la sostenibilidad de la actividad. Como resultado, se generan graves daños ambientales, como la deforestación, la contaminación del agua y la degradación del suelo. Estos impactos ambientales pueden tener consecuencias negativas para el clima de la región, lo que se traduce en un aumento de las temperaturas. Además, las olas de calor pueden afectar la productividad de los trabajadores y generar mayores costos en el sector de la salud. En el ámbito agrícola, las altas temperaturas pueden reducir el rendimiento de los cultivos y afectar la calidad de los alimentos. Por último, el aumento de las temperaturas puede tener impactos negativos sobre la biodiversidad, ya que muchas especies animales y vegetales dependen de condiciones climáticas específicas para su supervivencia.

Pregunta de investigación:

¿Qué regiones del departamento de Antioquia presentan las mayores temperaturas y cuáles son los factores que influyen en este fenómeno climático?

Fuente: Gobernación de Antioquia

Caso Valle de aburrá

Las altas temperaturas pueden ser explicadas a través de diversos factores, entre los que se encuentran la ubicación geográfica de la región, la presencia de grandes centros urbanos y la contaminación atmosférica. De acuerdo con estudios previos su ubicación en una cuenca rodeada por montañas, genera un efecto invernadero que aumenta la temperatura en la zona.

Además, la presencia de grandes ciudades como Medellín y el incremento de la actividad industrial, el transporte y la construcción, contribuyen a la emisión de gases de efecto invernadero y a la contaminación del aire, lo que agrava aún más la situación. Esto ha llevado a la implementación de políticas públicas y estrategias de adaptación y mitigación para reducir la vulnerabilidad ante el cambio climático y sus impactos en la región. En el valle de aburrá, la temperatura más alta generalmente ocurre durante los meses de febrero-marzo-abril y agosto-septiembre, que es la época más seca del año. Durante estos meses, la temperatura diurna promedio puede alcanzar los 27°C, aunque en días extremadamente calurosos puede superar los 30°C. Sin embargo, es importante tener en cuenta que el clima puede ser impredecible y la temperatura puede variar de un día a otro, incluso durante la misma temporada. El centro de Medellín, ubicado en la Comuna 10 La Candelaria, registra las temperaturas más altas de la ciudad, con una temperatura promedio de 32 grados frente a los 28 grados en zonas como El Poblado o Belén. La temperatura de la superficie llega a ser hasta cuatro grados centígrados más caliente que la del resto de la ciudad. La causa es que allí se concentran más construcciones, tanto de viviendas como comerciales, y también un mayor tráfico vehicular y peatonal. El lugar más frío del área metropolitana es la cuenca suroccidental, donde están San Antonio de Prado y la Estrella, y la suroriental, que alberga a El Poblado. Las altas temperaturas en Medellín pueden tener varios efectos en la salud y en el medio ambiente. Algunos de estos efectos pueden incluir:

  1. Deshidratación: Las altas temperaturas pueden causar una pérdida rápida de líquidos en el cuerpo, lo que puede provocar deshidratación y otros problemas de salud.

  2. Golpe de calor: El golpe de calor es una condición grave que puede ocurrir cuando el cuerpo no puede regular su temperatura correctamente debido a una exposición prolongada a altas temperaturas. Esto puede llevar a síntomas como confusión, calambres, estrés , mareo y pérdida del conocimiento.

  3. Problemas respiratorios: Las altas temperaturas pueden empeorar los síntomas de enfermedades respiratorias como el asma y la enfermedad pulmonar obstructiva crónica (EPOC).

  4. Contaminación del aire: Aumenta la contaminación del aire al aumentar la formación de ozono en la atmósfera, lo que puede provocar problemas respiratorios y otros problemas de salud.

  5. Sequía: Las altas temperaturas pueden aumentar la tasa de evaporación de la humedad del suelo y los cuerpos de agua, afectado el caudal de importantes quebradas como la Iguaná y Santa Elena y lo que puede provocar sequía y escasez de agua.

  6. Incendios forestales: Aumenta el riesgo de incendios forestales debido a que los materiales vegetales se vuelven más inflamables y secos en el valle de Aburrá en las zonas de interface, es decir en las zonas entre lo urbano y lo rural. En la ciudad de Medellín los sectores como Manrique (comuna 3), Villa Hermosa (comuna 8) , Robledo (Comuna 7) y Belén (comuna 16) y los Corregimientos de San Cristóbal y Altavista, que es donde se encuentra el Cerro de las Tres Cruces, en temporadas con temperatura alta tienen afectaciones de hectáreas al igual que partes altas de los barrios, donde ya se empieza a volver más rural y al igual que otros los municipios de valle de aburrá como Bello, Caldas, Envigado, Itagüí, La Estrella y Sabaneta.

Fuente: Siata

El pavimento, el concreto, las canchas sintéticas, la ausencia de árboles y la pérdida de vegetación por incendios contribuyen al efecto “isla de calor”, zonas donde la sensación térmica es superior. Por otro lado hay barrios donde su ubicación y su la altura en la que se encuentra, incrustado en una de las montañas laterales del valle, son zonas relativamente frías.

Algunos paquetes de interés

Código 
library(geoR); library(gstat) ;library(sf);library(raster)
library(rgdal)
library(vegan)
library(ade4)
library(tidyverse)
require(magrittr)
library(tmaptools) # Geocodificación y visualización
library(tidygeocoder) # Geocodificación
library(leaflet) # Visualización 
library(sm)
require(gstat)
require(MASS)
require(gt)
set.seed(1231)

Lectura de datos contorno de colombia (Archivo sp)

Se descarga de la base de datos del DANE la información necesaria para tener información del contorno de Colombia. En particular, se seleccionará el municipio de Antioquia. (https://geoportal.dane.gov.co/servicios/descarga-y-metadatos/descarga-mgn-marco-geoestadistico-nacional/)

En esta pagina, se encuentras diversos archivos en diferentes formatos espaciales. Para este caso, se utiliza la información de la pagina para crear el mapa del borde de “Colombia”.

Código 
col <- st_read("databases/MGN_DPTO_POLITICO.shp") # Archivo de interés
Reading layer `MGN_DPTO_POLITICO' from data source 
  `/home/daniel/OneDrive/Documentos (OneDrive-UNAL)/Semestre 2023-01/Spatial/proyecto_espacial/databases/MGN_DPTO_POLITICO.shp' 
  using driver `ESRI Shapefile'
Simple feature collection with 33 features and 0 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -81.73562 ymin: -4.229406 xmax: -66.84722 ymax: 13.39473
CRS:           NA
Código 
plot(col$geometry)

Código 
antioquia <- col[14,"geometry"]
border_ant <- st_coordinates(antioquia) #matriz del sistema coordenadas
border_ant <- as.data.frame(border_ant[,1:2])

¿Que se hizo en este paso?

Bueno, la base de datos general de Colombia (que por cierto, es muy pesada) está dividida por el numero de departamentos.

La fila 14 es el departamento de Antioquia y simplemente se le asigna a el vector homónimo los datos geométricos de dicho departamento.

Luego, se hacen dos cosas. Primero, se gráfica los bordes del departamento y segundo, se guarda en un dataframe la información de las columnas 1 y 2. Estas columnas son el sistema de coordenadas.

Cargue base de datos de interes

Sacado de : https://swat.tamu.edu/data/

La pagina suministra 2 archivos con 2 bases de datos diferentes. Una con datos de la precipitación y otro de la temperatura.

Ahora, con ambas bases de datos se pueden realizar preguntas de diversos tipos.

Código 
dt0 <- read_csv("databases/temperature.csv")
head(dt0)

Se escoje una fecha para trabajar

Debido a que no se trabajará con problemas espacio-temporales se seleccionó de manera arbitraria una fecha. Para este ejercicio se toman los datos de la manera ya mencionada para el 2019-01-10

Código 
dt2 <- dt0 %>% filter(date == "2019-01-10") %>% na.omit()
dt2 <- dt2[-c(12, 10, 11, 17, 28, 31, 29, 4, 1),]

Mapa de los puntos de extracción de la información

Graficamos los puntos donde se toma la información del día seleccionado

Código 
leaflet(antioquia) %>%
  addPolygons(color = "red") %>% 
  addTiles() %>% 
  addProviderTiles(providers$OpenTopoMap) %>% 
  addMarkers(lng = dt2$long, lat = dt2$lat,
             popup = paste0("[x = ", dt2$long,", y = ",
                            dt2$lat,"]"))

Sobre la autocorrelacion

Para el correcto análisis espacial es necesario probar autocorrelación espacial. Se procede a analizar las pruebas pertinentes.

Código 
datos2 <- as.matrix(data.frame(x=as.numeric(dt2$long),
                              y=as.numeric(dt2$lat)))

Creacion Dataframe Espacial

Código 
x1.2 <- as.numeric(dt2$tmin)
x2.2 <- as.numeric(dt2$tmax)
x3.2 <- as.numeric(dt2$elev)
df2 <- cbind(datos2, x1.2,x2.2,x3.2)
df2 <- data.frame(df2)
colnames(df2)<-c("X","Y", "Temp min", "Temp max", "elev")

geo_temp.max2 <- as.geodata(df2, coords = 1:2, data.col = 4, borders=TRUE)
geo_temp.max2$borders <- border_ant
datos2.sp <- df2
coordinates(datos2.sp) <- ~X+Y

Gráfico descriptivo para observar estructura espacial:

Código 
points(geo_temp.max2, col = "blue3", pt.divide = "equal",
       xlab = "Longitud", ylab = "Latitud")

Código 
plot(geo_temp.max2, lowess=TRUE)

Para ver si la normalidad gráfica cambia cuando se aplica una tendencia de segundo orden:

Código 
plot(geo_temp.max2, lowess=TRUE, trend = "2nd")

Prueba de normalidad

Código 
qqnorm(df2$`Temp max`, pch = 19)
qqline(df2$`Temp max`)

\[ H_0: La~temperatura~maxima~distribuye~de~forma~normal \] \

\[ H_1: La~temperatura~maxima~NO~distribuye~de~forma~normal \]

Código 
shapiro.test(df2$`Temp max`)

    Shapiro-Wilk normality test

data:  df2$`Temp max`
W = 0.90931, p-value = 0.04583

Dado que el valor p es de 0.04583, es tiene que la temperatura máxima en Antioquia registrada el 10 de Enero del 2019 no distribuye normal.

Ahora probamos con otra fecha para observar si los datos realmente en todas sus fechas no se acomodan a una distribución normal.

Otros datos

Suponga los datos de otra fecha (2019-10-21) para la temperatura máxima.

Gráfico descriptivo para observar estructura espacial:

Código 
dt <- dt0 %>% filter(date == "2019-10-21") %>% na.omit()
dt <- dt[-c(12, 10, 11, 17, 28, 31, 29, 4, 1),]

datos <- as.matrix(data.frame(x=as.numeric(dt$long),
                              y=as.numeric(dt$lat)))

x1 <- as.numeric(dt$tmin)
x2 <- as.numeric(dt$tmax)
x3 <- as.numeric(dt$elev)
df <- cbind(datos, x1,x2,x3)
df <- data.frame(df)
colnames(df)<-c("X","Y", "Temp min", "Temp max", "elev")

geo_temp.max <- as.geodata(df, coords = 1:2, data.col = 4, borders=TRUE)
geo_temp.max$borders <- border_ant

plot(geo_temp.max, lowess=TRUE, trend = "2nd")

Código 
# Conversión de coordenadas geográficas (LongLat)
# a coordenadas planas (UTM)
datos.sp <- df
coordinates(datos.sp) <- ~X+Y
proj4string(datos.sp) <- CRS("+init=epsg:4326")
datos.sp.utm <- spTransform(datos.sp, CRS("+init=epsg:9377"))

coordinates(border_ant) <- ~X+Y
proj4string(border_ant) <- CRS("+init=epsg:4326")
border_ant.utm <- spTransform(border_ant, CRS("+init=epsg:9377"))

geo_temp.max.utm <- datos.sp.utm %>% as.data.frame() %>%  as.geodata(., coords = 4:5, data.col = 2, borders=TRUE)
geo_temp.max.utm$borders <- border_ant.utm %>% as.data.frame()

plot(geo_temp.max.utm, lowess=TRUE, trend = "2nd")

Prueba Visual de Normalidad

Código 
qqnorm(df$`Temp max`, pch = 19)
qqline(df$`Temp max`)

Código 
shapiro.test(df$`Temp max`)

    Shapiro-Wilk normality test

data:  df$`Temp max`
W = 0.92992, p-value = 0.1222

Al realizar la misma prueba de hipótesis sobre la variable de temperatura máxima pero de otra fecha se tienen datos normales.

Lo que se desea enfatizar es que, a pesar de que el resultado de la primer muestra seleccionada ha rechazado el supuesto de normalidad de los datos, un fenómeno natural como la temperatura cabría esperar, por la naturaleza de su proceso, que siguiese un proceso normal. Dicha información puede ser constatada en otros artículos pero esta no es la preocupación principal del presente.

Es por esto que, para este caso,se asume que a pesar de la no normalidad presentada por el test de Shapiro-Wilk, se procede con los análisis estadísticos-espaciales.

A partir de este momento, se realizan análisis para dos fechas, dado que se busca demostrar las técnicas vistas en clase, se plantea entonces analizar los datos propuestos inicialmente y los datos para la fecha que si cumple los supuestos, esto con la intención de ver como la desviación de la normalidad podría afectar los resultados.

Código 
df.sp.utm <- datos.sp.utm %>% as.data.frame()

distancia1 <- dist(df.sp.utm[,4:5], diag=TRUE, upper=TRUE)

dif.temp1 <- dist(df.sp.utm[,2], diag=TRUE, upper=TRUE)^2

distancia2 <- dist(df2[,1:2], diag=TRUE, upper=TRUE)

dif.temp2 <- dist(df2[,4], diag=TRUE, upper=TRUE)^2

Grafico distancia vs disimilitud

Código 
datosmatriz1 <- data.frame(dif.prof=dif.temp1[lower.tri(dif.temp1)],
                           distancia1=distancia1[lower.tri(distancia1)])

datosmatriz2 <- data.frame(dif.prof=dif.temp2[lower.tri(dif.temp2)],
                           distancia2=distancia2[lower.tri(distancia2)])


par(mfrow = c(1,2))

# (2019-10-21) [Normal]
plot(x = datosmatriz1$distancia1, y = datosmatriz1$dif.prof,
     main = "2019-10-21 (Normal)", ylab = "Diferencia en temperatura",
     xlab = "Distancia (m)")
abline(lm(dif.prof~distancia1, data = datosmatriz1), lwd = 2,
       col = "red")

# (2019-01-10) [NO Normal]
plot(x = datosmatriz2$distancia2, y = datosmatriz2$dif.prof,
     main = "2019-01-10 (No Normal)", ylab = "Diferencia en temperatura",
     xlab = "Distancia (m/100k)")
abline(lm(dif.prof~distancia2, data = datosmatriz2), lwd = 2,
       col = "red")

Aunque no es muy marcado, se observa que, a distancias más pequeñas, la disimilitud es más pequeña. Se podría pensar existe una autocorrelación espacial para el gráfico de la izquierda en contra parte con el gráfico de la derecha que muestra una tendencia constante sugiriendo la no correlación espacial.

Test de Mantel

Se desea probar que:

\[ H_0: Hay~aleatoridad~espacial~~vs~~H_1:Hay~autocorrelación~espacial \]

Código 
mantel(distancia1, dif.temp1, na.rm = T)

Mantel statistic based on Pearson's product-moment correlation 

Call:
mantel(xdis = distancia1, ydis = dif.temp1, na.rm = T) 

Mantel statistic r: 0.1496 
      Significance: 0.042 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.104 0.140 0.178 0.219 
Permutation: free
Number of permutations: 999

Según este método a una significancia de 0.05 se rechaza la hipótesis nula acerca de aleatoriedad espacial. Entonces se concluye hay autocorrelación espacial.

Código 
mantel.rtest(distancia1, dif.temp1, nrepet=999)
Monte-Carlo test
Call: mantelnoneuclid(m1 = m1, m2 = m2, nrepet = nrepet)

Observation: 0.1496266 

Based on 999 replicates
Simulated p-value: 0.042 
Alternative hypothesis: greater 

     Std.Obs  Expectation     Variance 
1.8241208422 0.0002735359 0.0067037993

Con el test de Mantel se prueba que existe autocorrelación espacial; es decir, la temperatura es similar entre unidades geográficas próximas.

2019-01-10 (NO-Normal)

Código 
mantel(distancia2, dif.temp2)

Mantel statistic based on Pearson's product-moment correlation 

Call:
mantel(xdis = distancia2, ydis = dif.temp2) 

Mantel statistic r: 0.09854 
      Significance: 0.127 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.110 0.151 0.185 0.232 
Permutation: free
Number of permutations: 999
Código 
mantel.rtest(distancia2, dif.temp2)
Monte-Carlo test
Call: mantelnoneuclid(m1 = m1, m2 = m2, nrepet = nrepet)

Observation: 0.09854303 

Based on 99 replicates
Simulated p-value: 0.09 
Alternative hypothesis: greater 

     Std.Obs  Expectation     Variance 
 1.396905205 -0.001539747  0.005133160

Lamentablemente, en las dos pruebas no se cumple el test de mantel, por lo cual se dice que para el día 2019-01-10 existe evidencia para no autocorrelación espacial.

La hipótesis de autocorrelación espacial se vio afectada por la desviación ligera de normalidad, entonces se procede a trabajar con los datos del 2019-10-21 y se analizará la temperatura máxima (si es posible también la mínima) en antioquia para la fecha seleccionada.

Ahora se procede a extraer empíricamente la distancia máxima para los semivariogramas.

Código 
resumen <- summary(geo_temp.max.utm)
distancia <- resumen$distances.summary
dm <- distancia[2]/2; dm
     max 
177035.3

Calculamos los variogramas

Código 
vario1 <- variog(geo_temp.max.utm, option = "cloud")
variog: computing omnidirectional variogram
Código 
vario2 <- variog(geo_temp.max.utm, option = "cloud", max.dist = dm)
variog: computing omnidirectional variogram
Código 
par(mfrow=c(1,2))
plot(vario1, pch=21, bg="blue", lwd="black",
     main="Variograma sin ajuste\ncon regla empírica")
plot(vario2, pch=21, bg="red", lwd="black",
     main="Variograma ajuste\ncon regla empírica")

trend.l Ambos semivariogramas refuerzan la hipótesis inicial respecto a la correlación espacial, lo anterior es fácil de observar en el gráfico de la derecha pues para distancias pequeñas la semivarianza tiene valores más pequeños que para distancias intermedias y grandes, se observa un comportamiento en alguna medida “creciente” respecto a las distancias.

Variogramas tipo bin

Código 
bin1 <- variog(geo_temp.max.utm,option="bin")
variog: computing omnidirectional variogram
Código 
bin2 <- variog(geo_temp.max.utm,option="bin",max.dist=dm)
variog: computing omnidirectional variogram
Código 
par(mfrow=c(1,2))
plot(bin1, pch=21, bg="blue", lwd="black",
     main="Semivariograma sin ajuste\ncon regla empírica")
plot(bin2, pch=21, bg="red", lwd="black",
     main="Semivariograma ajuste\ncon regla empírica")

Se puede notar a simple vista que los semivariogramas no siguen una estructura de los modelos teóricos vistos. La idea será finalmente, tratar de ajustar uno con una medida aceptable para realizar predicciones.

Autocorrelacion con semivariograma

Código 
par(mfrow = c(1,1))
sm.variogram(datos.sp.utm@coords,
             datos.sp.utm@data$`Temp max`,model="independent",
             ylim = c(0,70))
Test of spatial independence: p =  0.195

Este semivariograma apunta que no existe una fuerte autocorrelación espacial puesto que la mayoría de puntos del semivariograma están dentro de las bandas de confianza, aun así, la prueba formal ya comprobó que si existe dicha correlación.

Código 
ini.vals1 <- c(10.5,168500)
fit_matern <- variofit(vario = bin2,ini.cov.pars = ini.vals1, fix.nugget =T,
                   weights = "npairs", min="optim", nugget = 0)
variofit: covariance model used is matern 
variofit: weights used: npairs 
variofit: minimisation function used: optim
Código 
plot(bin2, pch=21, bg="blue", lwd="black",
     main="Semivariograma con ajuste Matern")
lines(fit_matern, lwd = 2, col = "red")

En este caso, el vector ini.vals1 contiene dos valores iniciales para el modelo de covarianza de Matern.

  • El segundo valor (\(168500\)) se utiliza para el parámetro de rango, que determina la distancia a la cual la correlación entre los puntos de datos cae a 0.5.

  • El primer valor (\(10.5\)) se utiliza para el parámetro de sill, que controla la tasa a la que disminuye la correlación a medida que aumenta la distancia entre los puntos de datos.

Predicción Espacial

Código 
# Nueva regla empirica
resumen2 <- summary(geo_temp.max)
distancia2 <- resumen2$distances.summary
dm2 <- distancia2[2]/2

# Creación de nuevas bases de datos con la libreria sp
geo_df <- rename(df, Temperature = `Temp max`)

geo_df2 <- data.frame(geo_df)
geo_df3 <- data.frame(geo_df)
coordinates(geo_df2) <- ~X+Y

geo_df <- as.geodata(geo_df, coords = 1:2, var = 4)

# Creación de las geometrías
poly <- Polygon(border_ant)
polys <- Polygons(list(poly), 1)
sps <- SpatialPolygons(list(polys))

# Creación de los datos espaciales
sp_df <- SpatialPolygonsDataFrame(sps, data.frame(Temperature = 35.7))

# Creación de la maya con 10.000 puntos
grd <- makegrid(sp_df, n = 10000)
grd.ko <- as.data.frame(grd)
colnames(grd) <- c('X','Y')
grd_pts <- SpatialPoints(coords = grd, 
                         proj4string=CRS(proj4string(sp_df)))

grd_pts_in <- grd_pts[sp_df, ]
gdf <- as.data.frame(coordinates(grd_pts_in)) 

ggplot(gdf) +
  geom_point(aes(x=X,y=Y))+
  labs(title = "Malla de 10k puntos",
       subtitle = "Para hacer interpolación espacial")+
  theme_bw()

Código 
## Ajuste del varigrama en GeoR con un modelo matern
variog1 <- variog(geo_df, option="bin", max.dist=dm2)
variog: computing omnidirectional variogram
Código 
mod_ols1 <- variofit(variog1, ini.cov.pars = c(11.3, 1.6),
                     fix.nugget =T,weights = "npairs", 
                     min="optim", nugget = 0)
variofit: covariance model used is matern 
variofit: weights used: npairs 
variofit: minimisation function used: optim
Código 
# Kriging Ordinario en GeoR
pred_krig <- krige.conv(geo_df, locations = grd.ko,
                        borders = border_ant@coords,
                        krige = krige.control(nugget = 0.03, trend.d = "2nd",
                                              trend.l = "2nd",
                                              cov.pars = c(sigmasq = 3.9,
                                                           phi = 1.6)))
krige.conv: results will be returned only for prediction locations inside the borders
krige.conv: model with mean given by a 2nd order polynomial on the coordinates
krige.conv: Kriging performed using global neighbourhood
Código 
plot(border_ant@coords, main = "Kriging ordinario",
     xlab = "Longitud", ylab = "Latitud", pch = '-')
image(pred_krig, add=TRUE)
contour(pred_krig, add= TRUE, drawlabels = TRUE)

La influencia del mar y ríos en la temperatura es un factor importante a tener en cuenta en las regiones, debido a que el agua tiene una mayor capacidad térmica que la tierra, tiende a retener el calor durante más tiempo y liberarlo de manera más gradual. Esto puede resultar en temperaturas más suaves y moderadas en comparación con las áreas continentales cercanas. Además, la altitud del nivel del mar también puede tener un efecto en la temperatura. En general, a mayor altitud sobre el nivel del mar, las temperaturas tienden a disminuir debido a la disminución de la presión atmosférica y la menor retención de calor en la atmósfera. Por lo tanto, en regiones situadas a menor altitud sobre el nivel del mar, como las zonas costeras, es posible que se experimenten temperaturas más altas en comparación con áreas más elevadas. Sin embargo, es importante tener en cuenta que otros factores también pueden influir en las temperaturas, como la latitud, la dirección y fuerza de los vientos, la presencia de corrientes marinas, la cobertura nubosa, la vegetación y la urbanización, entre otros. Por lo tanto, es necesario realizar un análisis detallado y considerar múltiples variables para comprender completamente la influencia del mar y el nivel del mar en las temperaturas en una región específica.

Varianza de la predicción

Código 
plot(border_ant@coords, main = "Varianza (kriging ordinario)",
     xlab = "Longitud", ylab = "Latitud", pch = '-')
image(pred_krig, values = sqrt(pred_krig$krige.var),add=TRUE)
contour(pred_krig, values = sqrt(pred_krig$krige.var), add= TRUE,
        drawlabels = TRUE)

Validación Cruzada

Código 
valid.krigo <- xvalid(geo_temp.max, model = mod_ols1)
xvalid: number of data locations       = 22
xvalid: number of validation locations = 22
xvalid: performing cross-validation at location ... 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 
xvalid: end of cross-validation
Código 
summary(valid.krigo)
                Min.    1st Qu.      Median      Mean   3rd Qu.     Max.
errors     -3.675870 -1.0298969 -0.09482950 0.1341619 1.7769295 4.014368
std.errors -1.796986 -0.5284779 -0.05831335 0.0314131 0.8178424 1.679985
                 sd
errors     2.024287
std.errors 0.951569
Código 
plot(valid.krigo)

Según los resultados de la validación cruzada, se pueden sacar las siguientes conclusiones:

  • Sesgo (Mean): El sesgo promedio es de aproximadamente 0.134, lo que indica que, en promedio, el modelo de kriging ordinario tiende a sobreestimar ligeramente los valores de temperatura.

  • Error cuadrático medio (RMSE, sd): El error cuadrático medio es de aproximadamente 2.024, lo que indica la discrepancia promedio entre los valores observados y los valores predichos por el modelo. Un RMSE más bajo indica una mejor precisión del modelo.

  • Distribución de errores: Los errores en la validación cruzada varían desde -3.676 hasta 4.014. Los errores negativos indican que el modelo subestima la temperatura, mientras que los errores positivos indican sobreestimación. El rango y la desviación estándar (sd) de los errores muestran la dispersión de los errores en el conjunto de datos de validación.

En general, los resultados de la validación cruzada indican que el modelo de kriging ordinario tiene cierto sesgo hacia la sobreestimación y muestra una discrepancia promedio de alrededor de 2.024 entre los valores observados y los valores predichos.

Código 
# Grafico de los datos reales en antioquia
ggplot(data.frame(cbind(geo_temp.max$coords,
                        Temp.Max = geo_temp.max$data)))+
  geom_point(aes(X,Y, color = Temp.Max), size = 2)+
  coord_fixed(ratio =1)+
  scale_color_gradient(low = 'blue', high = 'orange')+
  geom_path(data = data.frame(geo_temp.max$border),
            aes(X,Y))+
  labs(title = "Datos de la temperatura maxima en Antioquia",
       subtitle = "Con su escala de color en grados celsius (C°)")+
  theme_bw()

Código 
variog2 <- variogram(Temperature ~ X+Y, data = geo_df2, cutoff = dm2)
mod_ols2 <- fit.variogram(variog2,
                          model = vgm(model = "Exp", psill = 11.3,
                                      range = 1.6),fit.method = 6)


# Crear un nuevo objeto SpatialPixelsDataFrame con las mismas coordenadas que grd.ko

pred_krigu <- krige(Temperature ~X+Y, ~X+Y, geo_df3, gdf,
                    model = mod_ols2)
[using universal kriging]
Código 
# Puntos de medición para el grafico de predicción
puntos_medicion <- SpatialPoints(coords = geo_df3[,1:2],
                                 proj4string = CRS(proj4string(sp_df)))

# Grafico de la interpolación espacial con kriging universal
p <- ggplot(aes(X,Y), data = pred_krigu)+
  geom_point(aes(color = var1.pred))+
  geom_point(aes(X,Y), data = geo_df3)+
  geom_polygon(aes(X,Y), fill = NA, color = "black",
               data = as.data.frame(border_ant@coords))+
  scale_color_gradient(low = "blue", high = "orange")+
  labs(title = "Predicciones con Kriging Universal",
     y = "latitud", x = "longitud", color = "C°")+
  theme_bw()

p

Varianza de la predicción

Código 
p <- ggplot(aes(X,Y), data = pred_krigu)+
  geom_point(aes(color = var1.var))+
  geom_point(aes(X,Y), data = geo_df3)+
  geom_polygon(aes(X,Y), fill = NA, color = "black",
               data = as.data.frame(border_ant@coords))+
  scale_color_gradient(low = "blue", high = "orange")+
  labs(title = "Predicciones con Kriging Universal",
     y = "latitud", x = "longitud", color = "C°")+
  theme_bw()

p

Conclusión

En este trabajo se ha analizado el comportamiento de las temperaturas en el departamento de Antioquia, Colombia, y los factores que influyen en este fenómeno climático.

Las subregiones del departamento de Antioquia que presentan las mayores temperaturas son el Bajo Cauca, el Magdalena Medio y el Urabá, los factores que influyen en este fenómeno climático son la altitud, la cercanía al mar, la vegetación, la humedad y las corrientes de aire. Además, estas subregiones tienen una altitud baja, lo que implica una menor presión atmosférica y una mayor radiación solar. Estas regiones también están cerca del mar Caribe o del río Magdalena, lo que les aporta humedad y calor.

  • Estas regiones tienen una vegetación escasa o degradada, lo que reduce la capacidad de evaporación y transpiración de las plantas. Por el contrario, las regiones que presentan las menores temperaturas son el Suroeste, el Oriente y el Occidente, con valores promedio inferiores. Estas regiones se caracterizan por tener una alta altitud, una baja radiación solar y una amplia cobertura vegetal. Estos resultados muestran que la temperatura en Antioquia está influenciada por factores geográficos, climáticos y ambientales que determinan la distribución espacial y temporal del calor.

Referencias