Visualiza datos del Censo 2024 en mapas a nivel de manzana con R

Ahora que recientemente lanzaron los datos 2024 del Censo de Población y Vivienda de Chile, me di unos minutos para visualizarlos. En este tutorial veremos dos formas de visualizar datos censales a nivel de manzana: con mapas estáticos en {ggplot2}, y con mapas interactivos en {mapgl}.

Datos censales

Descarga los datos cartográficos del Censo 2024 desde la página del INE, entrando a Cartografía Censal y luego descargando el archivo Cartografía País Censo 2024 (geoparquet).

Cargar cartografías censales

Los datos censales vienen en formato geoparquet, que es un formato moderno y eficiente para datos espaciales.

Podemos cargar datos en geoparquet con la función read_parquet() del paquete {arrow}.

Si no tenemos {arrow} instalado, lo instalamos con install.packages("arrow").

library(arrow)
manzanas <- read_parquet("Cartografia_censo2024_Pais/Cartografia_censo2024_Pais_Manzanas.parquet")

Con la función glimpse() de {dplyr} podemos echar un vistazo a las columnas que vienen en la tabla, pero no las mostraremos aquí porque son más de 200 😂

library(dplyr)
glimpse(manzanas)

Ahora que tenemos los datos, vamos a ver cómo visualizarlos en mapas!

Mapas estáticos

Para generar mapas a partir de cartografía censal, lo que haremos será: filtrar la base de datos para obtener la unidad territorial que nos interesa, luego convertir la cartografía a formato {sf}, y finalmente visualizar con {ggplot2}.

Preparar datos

Primero filtremos una comuna del país:

# filtrar datos
manzanas_comuna <- manzanas |> 
  filter(COMUNA == "LA FLORIDA")

Luego, convertimos a formato {sf}, para poder usar en R los polígonos que vienen en la cartografía, especificando un sistema de referencia de coordenadas (CRS) adecuado:

library(sf)

# convertir a sf
manzanas_comuna_sf <- manzanas_comuna |> 
  st_as_sf(crs = 4326) 

Ahora podemos ver la tabla de datos, que vienen con sus características espaciales incluídas:

manzanas_comuna_sf

Simple feature collection with 3584 features and 217 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -70.61478 ymin: -33.57023 xmax: -70.51814 ymax: -33.50389
Geodetic CRS:  WGS 84
# A tibble: 3,584 × 218
   OBJECTID   CUT COD_REGION REGION        COD_PROVINCIA PROVINCIA COMUNA AREA_C
 *    <int> <int>      <int> <chr>                 <int> <chr>     <chr>  <chr> 
 1   148637 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 2   148638 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 3   148639 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 4   148640 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 5   148641 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 6   148642 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 7   148643 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 8   148644 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
 9   148645 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
10   148646 13110         13 METROPOLITAN…           131 SANTIAGO  LA FL… URBANO
# ℹ 3,574 more rows
# ℹ 210 more variables: MANZENT <dbl>, DISTRITO <chr>, COD_DISTRITO <int>,
#   COD_LOCALIDAD <int>, COD_ZONA <int>, LOCALIDAD <chr>, COD_ENTIDAD <int>,
#   COD_MANZANA <int>, ENTIDAD <chr>, TIPO_MZ <chr>, COD_CATEGORIA <int>,
#   CATEGORIA <chr>, MZ_BASE_CENSO <int>, ID_ENTIDAD <dbl>, ID_LOCALIDAD <dbl>,
#   ID_DISTRITO <dbl>, ID_ZONA <dbl>, n_per <dbl>, n_hombres <dbl>,
#   n_mujeres <dbl>, n_edad_0_5 <dbl>, n_edad_6_13 <dbl>, n_edad_14_17 <dbl>, …

Este dataframe espacial posee una columna especial (SHAPE) que contiene los polígonos territoriales.

Cada fila del dataframe corresponde a una unidad territorial, y en la columna SHAPE se indica la forma o el polígono que ocupa esa unidad en el territorio.

Visualización

Vamos a generar mapas estáticos con {ggplot2}, un paquete muy poderoso y personalizable para visualización de datos.

Para visualizar estos datos en {ggplot2}, usamos la función de geometría geom_sf(), que está diseñada para trabajar con datos espaciales en formato {sf}. Esta función detecta la columna que contiene los datos espaciales y otras características espaciales, y permite visualizarlas como mapas.

library(ggplot2)

manzanas_comuna_sf |> 
  ggplot() +
  geom_sf() +
  theme_minimal(base_size = 10)

Con ésto anterior obtenemos el mapa base, generado a partir de los polígonos de la columna SHAPE.

Ahora podemos especificar una variable censal para visualizar, definiéndola en aes(). Luego podemos personalizar el gráfico agregando paletas de colores, ajustes visuales, y textos:

manzanas_comuna_sf |> 
  ggplot() +
  aes(fill = n_edad_60_mas) +
  geom_sf(color = "white", linewidth = 0.02) +
  scale_fill_fermenter(palette = 13) +
  theme_minimal(base_size = 10) +
  theme(axis.text.x = element_text(angle = 90, vjust = .5)) +
  guides(fill = guide_legend(title = "Población",
                             position = "top")) +
  labs(title = "Población de 60 años o más por manzana",
       subtitle = "Comuna de La Florida",
       caption = "Fuente: Censo 2024, INE")

Reutilizar mapa

Podemos repetir todo el proceso, cambiando sólo un par de líneas, para obtener el mapa de otra variable en una comuna distinta:

manzanas |> 
  # filtrar
  filter(COMUNA == "PUENTE ALTO") |> # comuna
  # convertir
  st_as_sf(crs = 4326) |> 
  # graficar
  ggplot() +
  aes(fill = n_discapacidad) + # variable
  geom_sf(color = "white", linewidth = 0.01) +
  scale_fill_fermenter(palette = 3) +
  theme_minimal(base_size = 10) +
  theme(axis.text.x = element_text(angle = 90, vjust = .5)) +
  guides(fill = guide_legend(title = "Población",
                             position = "top")) +
  labs(title = "Población con discapacidad por manzana",
       subtitle = "Comuna de Puente Alto",
       caption = "Fuente: Censo 2024, INE")

Lo bueno de disponer del código es que podemos automatizar la creación de mapas para distintas comunas o regiones, simplemente cambiando los filtros y las variables que queremos graficar, o podemos cambiar el nivel territorial cargando uno de los otros archivos, o bien, realizando una operación de resumen sobre los datos.

Mapa por comunas

Solamente como prueba, filtremos una región del país, luego agrupemos por comunas, y sumemos la población de una de las variables en cada comuna con summarize() para pasar desde un mapa por manzanas a un mapa por comunas:

sf_use_s2(FALSE)

comunas <- manzanas |> 
  # filtrar
  filter(REGION == "METROPOLITANA DE SANTIAGO") |> 
  # convertir
  st_as_sf(crs = 4326) |> 
  # resumir por comuna
  group_by(COMUNA) |> 
  summarize(n_ocupado = sum(n_ocupado), # sumar datos por grupo
            SHAPE = st_union(SHAPE)) # unir polígonos por grupo

Obtenemos como resultado una tabla de datos por comuna:

comunas

Simple feature collection with 52 features and 2 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -71.51187 ymin: -34.05721 xmax: -70.22923 ymax: -32.96357
Geodetic CRS:  WGS 84
# A tibble: 52 × 3
   COMUNA          n_ocupado                                               SHAPE
   <chr>               <dbl>                                  <MULTIPOLYGON [°]>
 1 ALHUÉ                2542 (((-71.24133 -34.05452, -71.2413 -34.05457, -71.24…
 2 BUIN                47409 (((-70.85738 -33.81385, -70.85736 -33.81381, -70.8…
 3 CALERA DE TANGO      5893 (((-70.83452 -33.66046, -70.83477 -33.66032, -70.8…
 4 CERRILLOS           36188 (((-70.73703 -33.51005, -70.73702 -33.51005, -70.7…
 5 CERRO NAVIA         52273 (((-70.76262 -33.42757, -70.76262 -33.42757, -70.7…
 6 COLINA              62631 (((-70.71874 -33.32434, -70.71884 -33.32437, -70.7…
 7 CONCHALÍ            54640 (((-70.68942 -33.38985, -70.68916 -33.3904, -70.68…
 8 CURACAVÍ            10081 (((-71.15734 -33.41664, -71.15734 -33.41663, -71.1…
 9 EL BOSQUE           64497 (((-70.68728 -33.5783, -70.68727 -33.57831, -70.68…
10 EL MONTE            13596 (((-71.05384 -33.69958, -71.05354 -33.69968, -71.0…
# ℹ 42 more rows

Ahora podemos visualizar un mapa de la población ocupada por comuna en la Región Metropolitana:

comunas |> 
  ggplot() +
  aes(fill = n_ocupado) + # variable
  geom_sf(color = "white", linewidth = 0.01) +
  colorspace::scale_fill_continuous_sequential(palette = "Sunset", na.value = "white",
                                               labels = scales::label_number(big.mark = "."),
                                               breaks = c(0, 50000, 100000, 150000, 200000, 250000)) +
  theme_minimal(base_size = 10) +
  theme(axis.text.x = element_text(angle = 90, vjust = .5),
        legend.key.width = unit(3, "mm")) +
  guides(fill = guide_colorbar(title = "Ocupados/as", position = "right")) +
  coord_sf(xlim = c(-70.83, -70.47), 
           ylim = c(-33.3, -33.68)) +
  labs(title = "Población ocupada",
       subtitle = "Región Metropolitana",
       caption = "Fuente: Censo 2024, INE")

Warning in prettyNum(.Internal(format(x, trim, digits, nsmall, width, 3L, :
'big.mark' and 'decimal.mark' are both '.', which could be confusing

Mapas interactivos

Para visualizar mapas interactivos en R existen varias opciones, como el paquete {leaflet}, pero yo quería usar el paquete {mapgl}, desarrollado por Kyle Walker que permite crear mapas interactivos usando MapLibre, una biblioteca de código abierto para mapas web.

Este paquete se caracteriza por el acceso que ofrece a usuarixs de R a mapas online interactivos de alta calidad y alto rendimiento.

Primero instalamos {mapgl}:

install.packages("mapgl")

Luego lo cargamos:

library(mapgl)

Y probamos su funcionamiento con un mapa vacío (prueba haciendo zoom o navegando el mapa con el cursor)

maplibre(center = c(-71.5, -33.0), 
         zoom = 3)

Ahora usemos {mapgl} con los datos censales. Podemos repetir el proceso de filtrado de la comuna o región que nos interese:

# filtrar y seleccionar
manzanas_comuna <- manzanas |> 
  filter(COMUNA == "VALPARAÍSO") |> 
  select(REGION, n_edad_60_mas, SHAPE)

Luego hacemos la misma conversión a {sf} para trabajar mejor con datos espaciales:

# convertir a sf
manzanas_comuna_sf <- manzanas_comuna |> 
  st_as_sf(crs = 4326)

Finalmente podemos visualizar los datos en un mapa interactivo con la función maplibre_view(), especificando la variable que queremos graficar en el mapa:

manzanas_comuna_sf |> 
  maplibre_view(column = "n_edad_60_mas")

¡Así de simple! 💜

Lo interesante es que {mapgl} también permite hacer cosas mucho más increíbles usando Mapbox, una plataforma de mapas de primer nivel que tiene muchas más opciones de personalización y estilos de mapas, incluyendo figuras tridimensionales, edificios, efectos visuales y más.