SIG con R

Contacto: wa.me/+573217850490

En este nuevo vídeo tutorial de mi canal de YouTube (youtube.com/fabiocastrollanos) damos un paso más allá del mapa temático tradicional y respondemos una pregunta que difícilmente se puede contestar con una sola capa: ¿la infraestructura vial de un municipio guarda relación con la pobreza de los territorios que lo rodean?

Para resolverlo aplicamos un análisis de Moran bivariado entre la longitud de vías y el Índice de Pobreza Multidimensional (IPM) a nivel municipal para toda Colombia, usando R. La variable de longitud de vías la construimos en un tutorial anterior, donde cuantificamos la red vial de cada municipio; aquí la ponemos a «conversar» con la pobreza.

1. El trasfondo

Un mapa de coropletas nos dice dónde hay valores altos o bajos de una variable. El Moran univariado da un paso más y mide si esos valores se agrupan en el espacio, es decir, si municipios «parecidos» tienden a ser vecinos. Pero en el análisis territorial la pregunta interesante casi siempre involucra dos fenómenos a la vez.

Ahí entra el Moran bivariado. En lugar de comparar una variable consigo misma, compara el valor de una variable en un municipio con el promedio de otra variable en sus vecinos (lo que técnicamente llamamos rezago espacial). El resultado es una clasificación local (LISA) que ubica a cada municipio en uno de cuatro cuadrantes:

• Alto – Alto: municipios con alta longitud de vías rodeados de vecinos con alta pobreza.
• Bajo – Bajo: poca vía, rodeada de vecinos con baja pobreza.
• Alto – Bajo: buena dotación vial en un entorno de baja pobreza.
• Bajo – Alto: poca infraestructura vial inmersa en un entorno de alta pobreza.

Este último cuadrante, el «Bajo – Alto», suele ser el más relevante para la política pública: señala bolsones de pobreza con conectividad deficiente, exactamente el tipo de territorio donde una inversión en infraestructura de transporte podría tener mayor impacto social.

Un detalle importante: el mapa final solo colorea los municipios estadísticamente significativos. La significancia se estima mediante un proceso de permutaciones (aleatorización), lo que nos permite distinguir los patrones reales de los que podrían aparecer por simple azar. Todo lo que no supera esa prueba se descarta como «Sin significancia».

En resumen, el objetivo del ejercicio no es hacer un mapa bonito, sino detectar de forma reproducible dónde coinciden espacialmente la carencia vial y la pobreza, para pasar de la intuición a la evidencia.

2. La preparación

Antes de correr una sola línea, conviene tener el entorno y los insumos ordenados. Este flujo se apoya por completo en R y en un puñado de librerías del ecosistema geoespacial:

• sf para leer y manipular los datos vectoriales.
• terra como soporte del manejo espacial.
• tidyverse para toda la limpieza y transformación de las tablas.
• rgeoda, la protagonista de este análisis: es la interfaz de GeoDa en R y contiene toda la maquinaria de autocorrelación espacial y estadísticos LISA.
• ggplot2 junto con ggspatial y patchwork para la composición cartográfica.
• cptcity y RColorBrewer para las paletas de color.

En cuanto a los datos de entrada, necesitamos combinar dos mundos:

• Capas vectoriales: la capa municipal que ya contiene el atributo de longitud de vías (derivada de datos de INVIAS), los límites departamentales del IGAC y una capa mundial de contexto obtenida con geodata.
• Datos tabulares: el IPM municipal en formato CSV, tomado del Departamento Nacional de Planeación (DNP).

El puente entre ambos es el código DANE del municipio (MPIO_CCNCT). Aquí un pequeño pero clásico dolor de cabeza: algunos códigos pierden el cero inicial al leerse (típicamente Antioquia y Atlántico), así que hay que restituir ese 0 antes de unir las tablas, o el cruce fallará silenciosamente. Por último, todas las capas se llevan a un mismo sistema de referencia (EPSG:4326) para que se alineen correctamente.

3. Comandos clave

Sin entrar línea por línea, estos son los comandos que sostienen todo el análisis y vale la pena entender a nivel conceptual:

queen_weights() → Es el punto de partida. Construye la matriz de vecindad bajo el criterio de la reina (dos municipios son vecinos si comparten un borde o incluso un solo vértice). Sin una definición de «quién es vecino de quién», no existe lo espacial en la estadística espacial: esta matriz es el cimiento de todo lo demás.

local_bimoran() → El corazón del ejercicio. Recibe la matriz de vecindad y las dos variables (vías e IPM) y calcula el Moran bivariado local para cada municipio, entregando su cuadrante y su nivel de significancia. Es la función que traduce la pregunta conceptual en un resultado cuantitativo.

lisa_clusters(), lisa_labels() y lisa_colors() → El trío que traduce el resultado a algo interpretable. La primera asigna la categoría a cada municipio, la segunda recupera las etiquetas estándar (Alto–Alto, Bajo–Alto, etc.) y la tercera aporta la paleta canónica de GeoDa, para que el mapa hable el mismo lenguaje visual que la literatura del tema.

inner_join() → El enlace entre lo tabular y lo geométrico. Es la operación que pega el IPM del DNP a la geometría municipal mediante el código DANE. Un cruce mal hecho aquí arruina todo el análisis, por eso la limpieza previa de los códigos es tan importante.

annotation_scale() y annotation_north_arrow() → Los elementos que convierten un gráfico en un mapa. La barra de escala y la flecha de norte no son adornos: son lo que le permite a cualquier lector interpretar correctamente la información. Es un paso que siempre insisto en no saltarse.

Conclusión

El Moran bivariado es una de esas técnicas que cambian la forma de leer un territorio: pasamos de describir variables por separado a detectar cómo se relacionan en el espacio. Aplicado a las vías y al IPM, nos entrega un mapa que no solo se ve bien, sino que señala con evidencia estadística dónde la falta de conectividad y la pobreza se dan la mano en Colombia.

Si quieres ver el paso a paso completo, con el código en pantalla y las decisiones metodológicas explicadas, te invito al vídeo en mi canal: youtube.com/fabiocastrollanos. Y si necesitas asesoría personalizada en temas de SIG y análisis espacial con R, escríbeme.

ontacto: Fabio Castro-Llanos — wa.me/+573217850490

🗺️ ¿Cómo medir la infraestructura vial a nivel municipal con R?

Acabo de subir un nuevo video donde pasamos directo a la práctica para resolver esto. Les muestro paso a paso cómo armar el flujo de trabajo para cuantificar la infraestructura de transporte usando datos geoespaciales.

En el tutorial vemos cómo:

Filtrar la red vial (tipos 1 a 4) para quedarnos solo con lo importante.

• Hacer una intersección espacial entre las vías y los municipios usando el paquete terra.

• Automatizar el cálculo para sacar la longitud total (en kilómetros) por cada municipio.

• Rematar con un mapa de calor para visualizar súper rápido qué zonas del país están mejor conectadas.

Si trabajan con SIG y quieren aprender a automatizar este tipo de geoprocesos en R, el tutorial te va a gustar.

Si me deseas contactar: wa.me/+573217850490

Los códigos los puedes descargar aquí:

Si deseas asesoría me puedes escribir: wa.me/+573217850490

Saludos

Contacto: wa.me/+573217850490

En este vídeo tutorial explicamos cómo realizar un mapa bivariado entre las variables de rezago escolar y trabajo infantil a nivel municipal para Colombia, los datos son tomados del Departamento Nacional de Planeación de Colombia.

En el análisis de datos geográficos, a menudo no basta con ver una sola variable en el mapa. A veces, la pregunta real es: ¿Cómo interactúan dos fenómenos en un mismo territorio? Para responder a esto, una de las herramientas más potentes y visualmente atractivas es el Mapa Bivariado. En este post, les explicaré la lógica detrás de un script desarrollado en R para analizar la relación entre el Rezago Escolar y el Trabajo Infantil.

¿Qué es un mapa bivariado?

A diferencia de un mapa de coropletas tradicional (que usa una escala de un solo color), el mapa bivariado utiliza una matriz de colores. Esto nos permite identificar no solo dónde los valores son altos o bajos para una variable, sino dónde se combinan ambas situaciones.

La lógica del proceso

Para lograr este resultado, el flujo de trabajo se divide en tres grandes pilares:

1. Preparación de los datos espaciales

El primer paso consiste en cargar la información base (capas departamentales y bases de datos georeferenciadas). En este punto, es vital definir los límites geográficos (bbox) para que el mapa se centre exactamente en el área de interés, eliminando ruidos visuales innecesarios.

2. Clasificación bivariada (El corazón del mapa)

Aquí es donde ocurre la «magia». Utilizando la librería biscale, los datos se clasifican en grupos (generalmente terciles o cuantiles). El script toma dos variables independientes:

• Eje X: Rezago Escolar.
• Eje Y: Trabajo Infantil.

El resultado es una nueva columna en nuestros datos que combina ambas dimensiones en categorías únicas (por ejemplo: «Bajo-Bajo», «Alto-Bajo», «Alto-Alto»).

3. Composición y Diseño Cartográfico

Un buen mapa debe ser tan informativo como estético. La lógica del diseño incluye:

• Capas de contexto: Se añade un fondo global y los límites departamentales para dar ubicación espacial.
• Simbolización: Se aplica una paleta de colores específica (como DkCyan) que permite diferenciar claramente los cuadrantes de la matriz.
• Elementos técnicos: Se incorporan la flecha de norte, la escala gráfica y etiquetas precisas para que cualquier lector pueda interpretar la escala.

El resultado final: La Leyenda

Lo más importante en un mapa bivariado es su leyenda en forma de diamante. El script no solo genera el mapa, sino que construye una leyenda independiente que se superpone al lienzo principal. Esta matriz nos dice, por ejemplo, que las zonas en los tonos más oscuros son aquellas donde ambos problemas sociales convergen con mayor intensidad.

Conclusión

Visualizar datos de esta manera nos permite identificar patrones críticos de forma inmediata, facilitando la toma de decisiones basada en datos geográficos. R, a través de potentes librerías como terra, sf y ggplot2, sigue siendo la herramienta líder para este tipo de análisis avanzado.

El mapa que creamos en el vídeo tutorial

Sin embargo, es posible replicarlo para cualquier ciudad del mundo, entre tanto hacemos uso de datos libres de Open Street Maps

Contacto: wa.me/+573217850490

En este vídeo tutorial explico como hacer uso de la librería geoBounds dentro de R para acceder a información de shapefiles político administrativo a distintos niveles. Y podemos lograr mapas como este:

Blog original: https://dieghernan.github.io/202602_geobounds/

Contacto: Fabio Castro – Llanos Whatsapp: wa.me/+573217850490

Contacto

Fabio Alexander Castro Llanos: fabioalexandercastro@gmail.com wa.me/+573217850490

En este pequeño proyecto desarrollé una aplicación interactiva en R Shiny para visualizar el Palmer Drought Severity Index (PDSI) a nivel municipal en Colombia. La idea es bastante simple: poder explorar cómo cambian las condiciones de sequía o humedad mes a mes usando un mapa interactivo.

El resultado final es una app donde el usuario puede seleccionar año y mes, y automáticamente el mapa se actualiza mostrando el valor del índice en cada municipio.

La aplicación está construida usando sf, leaflet, tidyverse y Shiny, junto con algunos ajustes visuales usando bslib para darle una apariencia moderna

📦 Carga de librerías y datos

El primer paso consiste en cargar las librerías necesarias y leer los datos que vamos a usar.

En este caso se utilizan dos datasets principales:

• Un shapefile simplificado de los municipios de Colombia
• Una tabla con los valores mensuales del PDSI

require(sf)
require(tidyverse)
require(shiny)
require(leaflet)
require(bslib)
require(RColorBrewer)
mpio <- st_read("./www/mpios_simplify.gpkg", quiet = TRUE)
znal <- read_csv(
  "./www/zonal_pdsi_year-monthly_gather.csv"
) |>
  select(MpCodigo, MpNombre, Depto, year, month, value)

Algo importante aquí es asegurarse de que el código del municipio tenga el mismo tipo de dato en ambas tablas para poder hacer correctamente el join.

mpio$MPIO_CCNCT <- as.character(mpio$MPIO_CCNCT)
znal$MpCodigo <- as.character(znal$MpCodigo)

🎛️ Interfaz de usuario (UI)

La interfaz se construyó con Bootstrap 5 usando el paquete bslib. Esto permite aplicar temas modernos y tipografías personalizadas.

theme = bs_theme(
  version = 5,
  bootswatch = "darkly",
  primary = "#38bdf8"
)

En el panel izquierdo se agregaron dos controles simples:

• selección de año
• selección de mes

selectInput("yr", "Año", choices = yrs_avl)
selectInput("mt", "Mes", choices = mth_avl)

Estos inputs serán los que controlan la actualización del mapa.

En el panel principal se incluye el mapa interactivo:

leafletOutput("map_pdsi")

Además agregué algunos detalles visuales como tarjetas (cards), tipografía moderna y enlaces a redes sociales.

⚡ Lógica reactiva en el servidor

La magia de Shiny ocurre en el server, donde se construye el flujo reactivo de la aplicación.

Primero se filtran los datos del índice según el año y el mes seleccionados:

tbl_pdsi <- reactive({
  znal |>
    filter(
      year == input$yr,
      month == input$mt
    )
})

Luego se realiza la unión con el shapefile municipal:

mpio_pdsi <- reactive({
  mpio |>
    left_join(
      tbl_pdsi(),
      by = c("MPIO_CCNCT" = "MpCodigo")
    )
})

De esta forma cada municipio queda asociado con su valor del índice PDSI.

🗺️ Construcción del mapa interactivo

Primero se inicializa el mapa base usando Leaflet:

output$map_pdsi <- renderLeaflet({
  leaflet() %>%
    addProviderTiles(providers$CartoDB.Positron) %>%
    setView(lng = -74, lat = 4.5, zoom = 6)
})

Después se actualizan los polígonos dinámicamente cada vez que cambia el año o el mes.

pal <- colorNumeric(
  palette = brewer.pal(11, "BrBG"),
  domain  = dat$value
)

Esta paleta permite representar:

• valores negativos → sequía
• valores positivos → humedad

Finalmente se dibujan los municipios en el mapa.

addPolygons(
  fillColor = ~pal(value),
  fillOpacity = 0.8,
  label = ~paste0(MpNombre, "<br>PDSI: ", round(value,2))
)

Además se agrega una leyenda dinámica que ayuda a interpretar los colores.

🌎 Resultado

El resultado es una pequeña aplicación interactiva que permite explorar rápidamente las condiciones de sequía o humedad en Colombia a escala municipal.

Este tipo de herramientas resulta especialmente útil para:

• monitoreo climático
• análisis agroclimático
• exploración de datos espaciales en tiempo real

Y demuestra lo poderoso que puede ser R para construir aplicaciones geoespaciales interactivas con relativamente pocas líneas de código.

Si requieres asesoría en temas de SIG contactame a mi correo fabioalexandercastro@gmail.com o whatsapp https://wa.me/+573217850490

En este nuevo tutorial muestro cómo calcular la estadística zonal del PDSI para todos los municipios de Colombia y generar un archivo espacial final listo para análisis multitemporal.

Trabajamos con datos previamente descargados de TerraClimate, específicamente el Palmer Drought Severity Index (PDSI), y realizamos el procesamiento usando un flujo reproducible en R.

🔎 ¿Qué hace el script?

El código realiza tres procesos fundamentales:

1️⃣ Lectura y organización de rasters

Se listan todos los .tif del PDSI recortados a Colombia y se ordenan correctamente por año.

fles <- dir_ls('./tif/pdsi', regexp = '.tif$')
fles <- grep('COL', fles, value = T)
fles <- mixedsort(fles)

2️⃣ Cálculo de estadística zonal

Utilizamos exactextractr::exact_extract() para calcular el promedio mensual del PDSI por municipio:

znal <- exactextractr::exact_extract(
  x = rstr, 
  y = st_as_sf(mpio), 
  fun = 'mean'
)

El resultado genera:

• 12 columnas (mean.PDSI_1 a mean.PDSI_12)
• Una fila por municipio
• Un identificador temporal por año

Todo queda organizado en una tabla estructurada.

3️⃣ Creación del archivo espacial final

Luego hacemos el join con el shapefile municipal y exportamos:

st_write(znal, './gpkg/zonal_pdsi.gpkg')

El resultado es un GeoPackage listo para:

• Mapas temáticos
• Dashboards en Shiny
• Modelos econométricos espaciales
• Análisis multitemporal por municipio

🎥 Mira el paso a paso completo

En el video del canal Un Geógrafo en YouTube explico:

• Por qué usar exactextractr en lugar de terra::extract
• Cómo estructurar datos multitemporales
• Buenas prácticas para zonal statistics
• Errores comunes y cómo evitarlos

👉 Te invito a ver el tutorial completo en el canal.

Si requieres asesoría en temas de SIG contactame a mi correo fabioalexandercastro@gmail.com o whatsapp https://wa.me/+573217850490

En este tutorial te explico cómo descargar y procesar datos climáticos del dataset TerraClimate directamente desde R, automatizando la descarga de múltiples años, recortando por país y exportando los resultados en formato GeoTIFF listo para análisis espacial.

El flujo que implementamos es completamente reproducible y escalable, ideal para proyectos de modelación climática, análisis de sequías o estudios de cambio climático.

🌎 ¿Qué es TerraClimate?

TerraClimate es un conjunto de datos climáticos mensuales de alta resolución (~4 km) que combina información de estaciones, reanálisis y climatologías históricas. Es ampliamente utilizado para estudios de:

• Balance hídrico
• Índices de sequía (como PDSI)
• Temperatura y precipitación
• Modelación ecológica
• Agricultura y seguridad alimentaria

En este ejemplo trabajamos con el Palmer Drought Severity Index (PDSI).

🧠 Objetivo del flujo en R

El script realiza cuatro procesos clave:

1. Genera automáticamente las URLs para múltiples años.
2. Descarga archivos NetCDF (.nc) desde el servidor THREDDS.
3. Recorta y enmascara los datos usando el límite oficial de Colombia.
4. Exporta cada resultado como archivo GeoTIFF.

Todo esto se ejecuta con una función aplicada iterativamente.

📦 Librerías utilizadas

Trabajamos con un stack moderno y eficiente

require(pacman)
p_load(terra, fs, sf, tidyverse, gtools, glue)

🔗 Construcción automática de URLs

Generamos las rutas directamente hacia el servidor:

urls <- glue('http://thredds.northwestknowledge.net:8080/thredds/fileServer/TERRACLIMATE_ALL/data/TerraClimateV1.1_PDSI_{1995:2020}.nc')

Esto evita descargar manualmente año por año y permite escalar fácilmente el rango temporal.

🗺 Descarga del límite oficial de Colombia

Utilizamos geodata::gadm() para obtener el shapefile del país:

col0 <- geodata::gadm(country = 'COL', level = 0, path = './tmpr')

Este objeto se usa posteriormente como máscara espacial.

⚙️ Función principal de descarga y procesamiento

La función ejecuta tres operaciones:

• Descarga del archivo .nc
• Recorte espacial (crop)
• Enmascaramiento (mask)
• Exportación a .tif

make.down <- function(urle){
  wrtn <- glue('{dout}/{basename(urle)}')
  download.file(
    url = urle, 
    destfile = wrtn, 
    mode = 'wb'
  )
  rstr <- rast(wrtn)
  rstr <- terra::crop(rstr, col0)
  rstr <- terra::mask(rstr, col0)
  wrtn <- glue('{dirname(wrtn)}/COL_{basename(wrtn)}')
  wrtn <- gsub('.nc', '.tif', wrtn)
  terra::writeRaster(x = rstr, filename = wrtn, overwrite = TRUE)
}

Luego simplemente:

map(urls, make.down)

📁 Resultado final

Obtendrás una carpeta con archivos como:

🚀 Ventajas del enfoque

✔ Automatización completa
✔ Flujo reproducible
✔ Uso de terra (más eficiente que raster)
✔ Escalable a cualquier variable de TerraClimate
✔ Compatible con flujos avanzados (SPEI, modelación de especies, análisis agrícola)

🎥 Mira el video completo

En el video del canal Un Geógrafo en YouTube explico cada paso con detalle

En este video presento un flujo completo y reproducible en R para identificar zonas climáticamente idóneas para el cultivo de café bajo sombra, inspirado en el artículo científico:

“A coffee corridor for biodiversity and livelihoods: climatic feasibility of shade coffee cultivation in western Rwanda”
Trees, Forests and People 🌱📄

Aunque el estudio original se desarrolla en África oriental, aquí adapto su enfoque metodológico al contexto colombiano, utilizando datos climáticos globales y herramientas SIG en R.

🧭 Área de estudio y datos utilizados

El análisis se centra en los departamentos de Quindío, Risaralda y Caldas, una de las principales regiones cafeteras de Colombia ☕🇨🇴.

Se emplean datos oficiales y abiertos:

• 🌡️🌧️ Variables bioclimáticas (BIO) de WorldClim v2.1
• 🗺️ Límites administrativos de GADM
• 📦 Paquetes especializados de R para análisis espacial (terra, sf, tmap, geodata)

📊 Variables climáticas consideradas

Siguiendo el artículo de referencia, se seleccionan cinco variables bioclimáticas clave para el cultivo de Coffea arabica:

• BIO1: Temperatura media anual 🌡️
• BIO2: Rango diurno de temperatura
• BIO5: Temperatura máxima del mes más cálido 🔥
• BIO6: Temperatura mínima del mes más frío ❄️
• BIO12: Precipitación anual 🌧️

Cada variable se recorta y enmascara a la zona de estudio para asegurar coherencia espacial.

🧮 Reclasificación climática (idoneidad)

Las variables climáticas se reclasifican en cinco niveles de idoneidad, desde:

• ❌ No apto
• ⚠️ Marginal
• ✅ Apto
• 🌿 Muy apto

Esto se realiza mediante matrices de reclasificación, construidas a partir de rangos climáticos reportados en la literatura agroecológica. El resultado es un conjunto de capas raster comparables entre sí.

⚖️ Análisis multicriterio (Weighted Linear Combination)

Para integrar las variables se aplica un Análisis Multicriterio por Combinación Lineal Ponderada (WLC):

• Todas las variables reciben el mismo peso (0.25)
• Se asume igual importancia relativa de cada factor climático
• El resultado es un mapa continuo de idoneidad climática

Finalmente, el índice se normaliza a un rango de 0 a 1, donde:

• 0 = climáticamente no apto
• 1 = condiciones óptimas para café bajo sombra

📈 Este enfoque permite sintetizar información climática compleja en un solo indicador espacial interpretable.

🗺️ Visualización cartográfica

El resultado se representa mediante un mapa temático en tmap, incorporando:

• 🎨 Paleta de colores gradual por clases de idoneidad
• 🧭 Rosa de los vientos
• 📏 Barra de escala
• 🏷️ Etiquetas departamentales
• 🖼️ Exportación en alta resolución (300 dpi)

El mapa final resume visualmente las zonas con mayor potencial climático para el café, bajo un enfoque reproducible y transparente.

🔁 Reproducibilidad y aplicación

Todo el proceso se desarrolla 100% en R, lo que permite:

• 🔁 Repetir el análisis en otras regiones
• 🌱 Adaptarlo a otros cultivos
• 🌍 Incorporar escenarios climáticos futuros
• 📊 Integrarlo en estudios de planificación territorial

Este tipo de análisis es ampliamente utilizado en agroecología, SIG ambiental y estudios de cambio climático, y aquí se presenta de forma accesible para la comunidad técnica.

¿Me deseas contactar?

Escribeme a mi whatsapp https://wa.me/+573217850490

El vídeo tutorial de YouTube aquí lo puedes encontrar: