Tutorial: Mapa en {ggplot2} con calles desde Open Street Map en R

En este tutorial crearemos un mapa de una región de Chile, y sobre el polígono geográfico aplicaremos otros elementos geográficos como calles, avenidas y carreteras, obtenidos desde Open Street Map (proveedor de mapas online abierto y comunitario).

Esto puede servir para crear visualizaciones de datos espaciales minimalistas que de todos modos entreguen elementos urbanos de referencia para que les usuaries puedan ubicarse mejor espacialmente.

Usaremos {dplyr} para manipular los datos, el paquete {ggplot2} para visualización de datos, {sf} para tratamiento de elementos espaciales, {rnaturalearth} para obtener mapas de cualquier país o región del mundo, y {osmdata} para obtener datos estaciales desde Open Street Map (OSM) por medio de su API pública.

Cargar paquetes

library(dplyr) #manipulación de datos
library(ggplot2) #gráficos
library(sf) #elementos espaciales
library(rnaturalearth) #mapas del mundo
library(osmdata) #datos open street map

Obtener mapa de país y de región

El primer paso es obtener un mapa de Chile, a nivel de regiones, usando la función ne_states() del paquete {rnaturalearth}. Luego filtramos el dataframe resultante para enfocarnos en una sola región del país.

# obtener mapa pero a nivel de regiones
mapa_1 <- rnaturalearth::ne_states(country = "Chile") |> 
    select(name, name_es, iso_3166_2, code_hasc, geometry)

mapa_2 <- mapa_1 |> 
    filter(name == "Los Lagos")

mapa_2

Simple feature collection with 1 feature and 4 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -74.40722 ymin: -44.045 xmax: -71.65926 ymax: -40.27304
Geodetic CRS:  WGS 84
       name   name_es iso_3166_2 code_hasc                       geometry
1 Los Lagos Los Lagos      CL-LL     CL.LG MULTIPOLYGON (((-71.76926 -...

Previsualizamos la región obtenida usando {ggplot2}:

mapa_2 |> 
    ggplot() +
    geom_sf(linewidth = .4,
            color = "grey80", fill = "grey95") +
    theme_classic()

Luego realizamos el primer paso para obtener los datos desde Open Street Map, que es definir las coordenadas de la caja delimitadora o bounding box del lugar del que vamos a solicitar datos. Por ejemplo, una ciudad específica. Si bien esto se puede hacer con la función osmdata::getbb(), los resultados no siempre coinciden con lo que esperamos (la región o ciudad puede salir más recortada de lo esperado), por lo que preferimos usar directamente la caja del mapa que queremos visualizar, usando sf::st_bbox().

# obtener perímetro o área abarcada por el mapa
caja <- mapa_2 |> st_bbox()

# opcionalmente, se podría hacer con el nombre del lugar
# area <- getbb("Región de Los Lagos, Chile")

caja

     xmin      ymin      xmax      ymax 
-74.40722 -44.04500 -71.65926 -40.27304 

Obtener datos geográficos desde Open Street Map

Teniendo estas coordenadas, se solicitan los datos haciendo un query a la API Overpass de OSM. En esta solicitud se especifica la característica del mapa que se va a obtener, y sus valores; en este caso, queremos obtener las calles de una región, por lo que la feature corresponde a highway. Dividimos la query en varias partes, para tener en objetos separados las calles de distintos tamaños (a grandes rasgos: autopistas y avenidas, calles secundarias, y calles interiores).

Considerar que cada query puede tomar un tiempo, dependiendo de la cantidad de datos, y que OSM es un proyecto soportado por voluntarios, por lo que se recomienda no abusar de la API (hacer la mínima cantidad de solicitudes, no hacer muchas solicitudes seguidas, etc.)

# obtener calles dentro de la caja
calles_grandes <- caja |> 
    opq(timeout = 500) |> 
    add_osm_feature(key = "highway", 
                    value = c("motorway", "primary", "motorway_link", "primary_link")) |> 
    osmdata_sf()

calles_medianas <- caja |> 
    opq(timeout = 500) |> 
    add_osm_feature(key = "highway", 
                    value = c("secondary", "tertiary", "secondary_link", "tertiary_link")) |> 
    osmdata_sf()

calles_chicas <- caja |>
    opq(timeout = 500) |>
    add_osm_feature(key = "highway",
                    value = c("residential", "living_street")) |>
    osmdata_sf()

calles_grandes

Object of class 'osmdata' with:
                 $bbox : -44.0449981164699,-74.4072159499999,-40.2730445289999,-71.6592646899999
        $overpass_call : The call submitted to the overpass API
                 $meta : metadata including timestamp and version numbers
           $osm_points : 'sf' Simple Features Collection with 29125 points
            $osm_lines : 'sf' Simple Features Collection with 3270 linestrings
         $osm_polygons : 'sf' Simple Features Collection with 0 polygons
       $osm_multilines : NULL
    $osm_multipolygons : NULL

Habiendo obtenido los datos geográficos desde OSM, recibimos objetos que contienen los elementos espaciales. Extraemos las líneas de cada uno de los objetos, para quedarnos sólo con las calles. Opcionalmente, podemos recortar los elementos espaciales obtenidos con la figura del mapa, usando el mapa como figura de recorte para que las líneas de las calles no se salgan de los márgenes del mapa.

# extraer toda la geometría
calles_grandes_a <- calles_grandes$osm_lines
calles_medianas_a <- calles_medianas$osm_lines
calles_chicas_a <- calles_chicas$osm_lines

# recortar geometría para que calce dentro del polígono del mapa
calles_grandes_b <- st_intersection(calles_grandes$osm_lines, mapa_2)
calles_medianas_b <- st_intersection(calles_medianas$osm_lines, mapa_2)
calles_chicas_b <- st_intersection(calles_chicas$osm_lines, mapa_2)

Visualizar mapa con calles

Finalmente, podemos aplicar las capas obtenidas a nuestro gráfico de {ggplot2}. Podemos ajustar la apariencia de cada tipo de calle por separado dado que las ubicamos en capas distintas.

    ggplot() +
    # capa del país entero (como fondo)
    geom_sf(data = mapa_1, 
            linewidth = .4,
            color = "grey95", fill = "grey98") +
    # capa del mapa de la región
    geom_sf(data = mapa_2,
            linewidth = .4, 
            color = "grey80", fill = "grey95") +
    # capas osm
    geom_sf(data = calles_medianas_b, 
            linewidth = .3, alpha = .3) +
    geom_sf(data = calles_grandes_b, 
            linewidth = .4, alpha = .3) +
    # recortar vista
    coord_sf(xlim = caja[c(1, 3)],
             ylim = c(-42.4, -40.3)) +
    theme_classic()

Otro ejemplo del mapa visto desde más cerca, incluyendo la capa de calles interiores:

ggplot() +
    # capa del mapa
    geom_sf(data = mapa_2,
            fill = "grey95") +
    # capas osm
    geom_sf(data = calles_grandes_b, 
            linewidth = .3, alpha = 1) +
    geom_sf(data = calles_medianas_b, 
            linewidth = .3, alpha = .4) +
    geom_sf(data = calles_chicas_b, 
            linewidth = .2, alpha = .2) +
    # recortar vista
    coord_sf(xlim = c(-73.18, -73.09),
             ylim = c(-40.6, -40.55)) +
    theme_classic()