How to combine ggplot2 and autoplot.OpenStreetMap

Question

I am having trouble displaying an osm map under a series of points.

The shapefile is contained within a dataframe (df).

when I run this code,

df <- st_as_sf(df, coords = c("Longitude", "Latitude")) 
ggplot(df) + geom_sf() + theme_void()

I get this result:

and when I run this code,

dc <- openmap(c(19.353014580954277,-99.14118184937192), 
                    c(19.287683266903073,-99.03186392856783), 
                    zoom = NULL,
                    type = "osm", 
                    mergeTiles = TRUE)
autoplot.OpenStreetMap(dc) + theme_void()

I get this result:

I've tried multiple ways and read a lot of tutorials and cheat sheets and still don't understand how to combine these two to create one image,

My personal choice for this is maptiles + tidyterra for map tiles, see https://stackoverflow.com/questions/75214294/how-to-properly-combine-a-polygon-map-on-top-of-a-raster-map-in-r I think is easier — dieghernan, Apr 01 '23 at 15:02
thanks for the reply, but that post is even more confusing and besides, I'm not working in the US so those US specific functions don't work for this :/ — larrote, Apr 03 '23 at 13:08
I added an answer on this. See how to switch from OpenStreetMap to maptiles + tidyterra. Note also that OpenStreetMap package has not been updated in the last 4 years, so it is not actively maintained and IMHO outdated. — dieghernan, Apr 18 '23 at 06:49

score 1 · Answer 1 · answered Apr 18 '23 at 06:45

Expanding my comment, I do believe the approach with maptiles + tidyterra is easier for plotting basemaps. Note that you can retrieve also map tiles for OpenStreetMaps with maptiles, but you would get a SpatRaster (spatial raster) that is truly a geographic object. Let's start by creating your df with mock data from OpenStreetMap using nominatimlite (you didn't provide this data on the example so I had to get some points in your study area):

# Mock your df using restaurants on the bbox
thebbox <- c(
  -99.14118184937192, 19.353014580954277,
  -99.03186392856783, 19.287683266903073
)

library(nominatimlite)
df <- geo_amenity(thebbox,
  amenity = "restaurant", limit = 50,
  lat = "Latitude", lon = "Longitude"
)

Now we start to follow your code. One thing t note is that when you convert df to a sf object you should add the corresponding CRS to the object itself. As Longitude/Latitude are geographic coordinates, it is safe enough to consider crs = "+proj=lonlat":

# Start here
library(ggplot2)
library(sf)

df <- st_as_sf(df,
  coords = c("Longitude", "Latitude"),
  # Add crs
  crs = "+proj=lonlat"
)

ggplot(df) +
  geom_sf() +
  theme_void()

Good. Now, for getting base maps we can use maptiles with provider = "OpenStreetMap".

An important thing to bear in mind is that map tiles services usually provides the images in Mercator projection (i.e. Google Maps, etc.). So you may want to retrieve and plot the base tiles in Mercator.

Using a different CRS would cause the tile to be resampled and the results would not be as clean as the ones you would get using Mercator:


# Now use maptiles instead, with Mercator to get better results
df_merc <- st_transform(df, 3857)

library(maptiles)
dc <- get_tiles(df_merc, provider = "OpenStreetMap", zoom = 13)

#> Data and map tiles sources:
#> © OpenStreetMap contributors

Last step: Use ggplot2 and tidyterra to plot the RGB tile. tidyterra add additional geoms to ggplot2 for ploting SpatRaster objects:


# And finally plot with ggplot2 + tidyterra
library(tidyterra)

ggplot() +
  geom_spatraster_rgb(data = dc) +
  geom_sf(data = df) +
  coord_sf(crs = 3857)

^{Created on 2023-04-18 with reprex v2.0.2}

And that's all. A short version of the code (assuming you already have df) is:

library(ggplot2)
library(sf)
library(maptiles)
library(tidyterra)

df <- st_as_sf(df,
  coords = c("Longitude", "Latitude"),
  # Add crs
  crs = "+proj=lonlat"
)
df_merc <- st_transform(df, 3857)
dc <- get_tiles(df_merc, provider = "OpenStreetMap", zoom = 13)
ggplot() +
  geom_spatraster_rgb(data = dc) +
  geom_sf(data = df) +
  coord_sf(crs = 3857)

score 0 · Answer 2 · edited Apr 03 '23 at 18:29

instead of geom_sf() I used geom_point().
I also had to set the projection for the osm raster

dc <- openmap(c(19.353014580954277,-99.14118184937192), 
                         c(19.287683266903073,-99.03186392856783), 
                         zoom = NULL,
                         type = "osm", 
                         mergeTiles = TRUE)
     
#2.4: plotear mapas
dc <- openproj(dc, projection = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")

autoplot.OpenStreetMap(dc) + geom_point(data=df, aes(x=Longitude, y=Latitude)) + theme_void()

enter image description here

How to combine ggplot2 and autoplot.OpenStreetMap

2 Answers2