Difficulty plotting boroughs in R

Question

I'm trying to plot the boroughs of Montreal in R and show the name of each borough on the map. However, the following code doesn't seem to plot the map properly as it doesn't show all boroughs:

library(osmdata)
library(ggplot2)

# Get boundary data for Montreal
montreal_boundaries <- getbb("Montreal") %>%
  opq(timeout = 50) %>%
  add_osm_feature(key = "boundary", value = "administrative") %>%
  osmdata_sf()

# Plot the map using ggplot2
ggplot() +
  geom_sf(data = montreal_boundaries$osm_multipolygons, fill = "grey", color = "black") +
  coord_sf() +
  theme_void() +
  labs(title = "Montreal Boroughs")

Could someone help me correct the code so that all Montreal boroughs are properly displayed on the map, and each borough is labeled?

Thank you!

Answer 1

The issue is that your OSM data contains data on multiple administrative levels which overlap and depending on the order of the data are plotted on top of each other. Hence, to show all boroughs you have to filter (i.e. get rid of admin_level == "") and arrange your data properly:

library(osmdata)
library(ggplot2)

# Get boundary data for Montreal
montreal_boundaries <- getbb("Montreal") %>%
  opq(timeout = 50) %>%
  add_osm_feature(
    key = "boundary",
    value = "administrative"
  ) %>%
  osmdata_sf()

boundaries <- montreal_boundaries$osm_multipolygons |>
  subset(
    select = c(osm_id, name, admin_level, geometry)
  )
# Order by admin level
boundaries <- boundaries[
  order(as.numeric(boundaries$admin_level), decreasing = FALSE),
]

ggplot(data = boundaries) +
  geom_sf(
    # Drop admin_level == ""
    data = ~ subset(., !admin_level %in% c("")),
    aes(
      fill = factor(admin_level),
    ),
  ) +
  theme_void() +
  labs(title = "Montreal Boroughs")