Este cuaderno muestra cómo hacer un mapa interactivo de la producción agrícola por municipio en un departamento de Colombia, usando R.
Se toman los datos de las Evaluaciones Agropecuarias Municipales (EVA) de la UPRA, se filtran para el departamento de Cundinamarca, el año 2022 y el grupo de cultivo frutales, se unen a los polígonos de los municipios del Marco Geoestadístico Nacional (MGN) del DANE y se visualizan con la librería Leaflet.
Se basa en dos materiales del curso de Geomática Básica: el cuaderno How to read and subset .xlsx agricultural statistics data (Lizarazo, 2023) y el Tutorial de RLeaflet (Campos, 2019).
Se definen el departamento (con su código DANE), el grupo de cultivo y el año que se van a analizar:
mi_depto <- "Cundinamarca"
mi_cod_depto <- "25"
mi_grupo <- "FRUTALES"
mi_anio <- 2022Se cargan las librerías:
library(readr) # leer archivos csv
library(dplyr) # manipular datos (filtrar, agrupar, resumir, unir)
library(janitor) # limpiar nombres de columnas
library(stringi) # quitar tildes
library(sf) # datos geoespaciales (características simples)
library(leaflet) # mapas interactivos
library(htmlwidgets) # guardar el mapa como htmlComo los nombres pueden venir con o sin tildes y en mayúsculas o minúsculas, se crea una función que los deja todos iguales para poder compararlos:
normalizar <- function(x) toupper(trimws(stri_trans_general(x, "Latin-ASCII")))Los datos de la EVA 2019-2025 (Base Agrícola) están publicados por la UPRA en el Portal de Datos Abiertos de Colombia (datos.gov.co). La primera vez se descargan como archivo csv; las siguientes veces se lee el archivo ya guardado.
archivo <- "EVA_2019_2025.csv"
if (!file.exists(archivo)) {
options(timeout = 600)
download.file("https://www.datos.gov.co/resource/uejq-wxrr.csv?$limit=1000000",
archivo, mode = "wb")
}
eva <- read_csv(archivo, show_col_types = FALSE) %>% clean_names()
dim(eva)## [1] 166732 18Los nombres originales de las columnas tienen tildes y caracteres especiales. Se seleccionan solo las columnas necesarias y se les ponen nombres sencillos. El código del municipio se deja como texto de 5 dígitos, que es el formato del código DANE:
columna <- function(patron) grep(patron, names(eva), value = TRUE)[1]
eva <- eva %>%
transmute(
COD_MUN = sprintf("%05d", as.integer(.data[[columna("(cod|c_digo).*mun")]])),
DPTO = .data[[columna("^departamento$")]],
MUNICIPIO = .data[[columna("^municipio$")]],
GRUPO = .data[[columna("^grupo")]],
YEAR = as.integer(.data[[columna("^(ano|a_o|anio|year)$")]]),
PRODUCCION = as.numeric(.data[[columna("^produc")]])
)
head(eva)## # A tibble: 6 × 6
## COD_MUN DPTO MUNICIPIO GRUPO YEAR PRODUCCION
## <chr> <chr> <chr> <chr> <int> <dbl>
## 1 05001 Antioquia Medellín Hortalizas 2021 128
## 2 05001 Antioquia Medellín Hortalizas 2021 152
## 3 05001 Antioquia Medellín Hortalizas 2022 127.
## 4 05001 Antioquia Medellín Hortalizas 2022 152
## 5 05001 Antioquia Medellín Hortalizas 2023 133
## 6 05001 Antioquia Medellín Hortalizas 2023 152Primero se filtran los datos del departamento y el año:
depto_eva <- eva %>%
filter(normalizar(DPTO) == normalizar(mi_depto), YEAR == mi_anio)
depto_eva %>% count(GRUPO, sort = TRUE)## # A tibble: 8 × 2
## GRUPO n
## <chr> <int>
## 1 Frutales 652
## 2 Hortalizas 495
## 3 Cereales 350
## 4 Leguminosas 329
## 5 Raíces y tubérculos 317
## 6 Cultivos tropicales tradicionales 156
## 7 Cultivos para condimentos, bebidas medicinales y aromáticas 38
## 8 Oleaginosas 6Luego se selecciona el grupo de cultivo y se suma la producción de cada municipio:
produccion <- depto_eva %>%
filter(normalizar(GRUPO) == normalizar(mi_grupo)) %>%
group_by(COD_MUN, MUNICIPIO) %>%
summarize(prod_total = sum(PRODUCCION, na.rm = TRUE), .groups = "drop") %>%
arrange(desc(prod_total))
produccion## # A tibble: 114 × 3
## COD_MUN MUNICIPIO prod_total
## <chr> <chr> <dbl>
## 1 25649 San Bernardo 88442.
## 2 25386 La Mesa 56226.
## 3 25035 Anapoima 41141.
## 4 25123 Cachipay 29687.
## 5 25269 Facatativá 28609.
## 6 25524 Pandi 27271.
## 7 25740 Sibaté 24359
## 8 25599 Apulo 21331
## 9 25815 Tocaima 20511
## 10 25878 Viotá 18955.
## # ℹ 104 more rowsSe guarda la tabla en formato csv:
write.csv(produccion, "produccion_municipios.csv", row.names = FALSE)Los límites municipales se toman del Marco Geoestadístico
Nacional MGN2025 del DANE, nivel municipio, en formato
GeoPackage. Se busca el archivo en la carpeta de trabajo y se lee como
un objeto sf:
(mgn <- list.files(".", pattern = "MGN.*MPIO.*\\.gpkg$", recursive = TRUE, full.names = TRUE))## [1] "./Proyecto_ciudades/MGN2025_MPIO/MGN_ADM_MPIO_GRAFICO (Geopackage).gpkg"colombia <- st_read(mgn[1])## Reading layer `MGN_ADM_MPIO_GRAFICO (Geopackage)' from data source
## `C:\GEOMATICA2026\Restudio\Proyecto_ciudades\MGN2025_MPIO\MGN_ADM_MPIO_GRAFICO (Geopackage).gpkg'
## using driver `GPKG'
## Simple feature collection with 1122 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -81.73562 ymin: -4.229406 xmax: -66.84722 ymax: 13.39473
## Geodetic CRS: WGS 84names(colombia)## [1] "dpto_ccdgo" "mpio_ccdgo" "mpio_cdpmp" "dpto_cnmbr" "mpio_cnmbr"
## [6] "mpio_crslc" "mpio_tipo" "mpio_narea" "mpio_nano" "shape_Leng"
## [11] "shape_Area" "geom"Se seleccionan los municipios del departamento por su código DANE:
mpios <- colombia %>% filter(dpto_ccdgo == mi_cod_depto)
nrow(mpios)## [1] 116El MGN trae el código DANE de 5 dígitos de cada
municipio en la columna mpio_cdpmp, que es el mismo código
de la columna COD_MUN de la EVA. Por eso la unión se hace
por código y no por nombre, lo que evita problemas con tildes o nombres
escritos de forma distinta.
mapa_datos <- mpios %>% left_join(produccion, by = c("mpio_cdpmp" = "COD_MUN"))
# Municipios con producción que no encontraron su polígono (debe estar vacío)
setdiff(produccion$COD_MUN, mpios$mpio_cdpmp)## character(0)# Municipios sin producción de frutales en el año (quedan en gris)
sum(is.na(mapa_datos$prod_total))## [1] 2Leaflet trabaja con coordenadas geográficas WGS84 (EPSG 4326), así que se verifica y, si hace falta, se transforma:
mapa_datos <- st_transform(mapa_datos, 4326)Pocos municipios producen mucho más que los demás, así que una escala continua dejaría casi todo el mapa del mismo color. Por eso se usan quintiles: los municipios se dividen en 5 grupos con la misma cantidad cada uno. Los municipios en gris no reportaron producción de frutales en 2022.
pal <- colorQuantile("viridis", domain = mapa_datos$prod_total,
n = 5, na.color = "#cccccc")
popup <- paste0(
"<b>Municipio: </b>", mapa_datos$mpio_cnmbr, "<br>",
"<b>Código DANE: </b>", mapa_datos$mpio_cdpmp, "<br>",
"<b>Producción ", tolower(mi_grupo), " ", mi_anio, ": </b>",
ifelse(is.na(mapa_datos$prod_total), "sin datos",
paste(format(round(mapa_datos$prod_total), big.mark = ",", trim = TRUE), "t"))
)
mapa <- leaflet(mapa_datos) %>%
addProviderTiles(providers$Esri.WorldTerrain) %>%
addPolygons(
fillColor = ~pal(prod_total),
fillOpacity = 0.6,
color = "#444444",
weight = 1,
popup = popup,
highlightOptions = highlightOptions(weight = 3, color = "white", bringToFront = TRUE)
) %>%
addLegend(
position = "topright",
pal = pal,
values = ~prod_total,
title = paste0("Producción ", tolower(mi_grupo), "<br>", mi_anio, " (quintiles)")
)
mapa