How to combine a raster layer with sf objects in ggplot?

[kbzsl]

Sorry, for my novice question, but I cannot figure out what’s gone wrong when combining a raster background and other sf objects using the basemaps and ggplot2 packages. Thank you.

Here is my code, it should draw a square around Margaret Island, Budapest

library(tidyverse)
library(sf)
library(basemaps)

set_defaults(map_service = "esri", map_type = "world_street_map")

p1 = st_point(c(19.040767, 47.540033))
p2 = st_point(c(19.065183, 47.541033))
p3 = st_point(c(19.047933, 47.509817))
p4 = st_point(c(19.029983, 47.512850))

l1 = st_sfc(st_linestring(c(p1, p2)))
l2 = st_sfc(st_linestring(c(p2, p3)))
l3 = st_sfc(st_linestring(c(p3, p4)))
l4 = st_sfc(st_linestring(c(p4, p1)))

bp_margit <- st_as_sf(c(l1, l2, l3, l4), crs = 4326)

ggplot() + 
  geom_sf(data = bp_margit, color = "red") 

ggplot() + 
  basemap_gglayer(bp_margit) + 
  coord_sf() +
  scale_fill_identity() 

ggplot() + 
  basemap_gglayer(bp_margit) + 
  geom_sf(data = bp_margit, color = "red") +
  coord_sf() +
  scale_fill_identity()

[16EAGLE]

This is a projection issue. basemap returns maps in the Web Mercator proejction (EPSG: 3857), the standard used by many tiling services (see section Value under ?basemap). That means you are trying to plot two objects with different CRSs on top of each other.

Just transform your polygon to Web Mercator first, e.g. by extracting the CRS form the pulled map like this

bp_margit <- st_transform(
  bp_margit,
  crs = st_crs(basemap_stars(bp_margit))
)

or manually setting EPSG 3857

bp_margit <- st_transform(bp_margit,  crs = st_crs(3857)
)

before plotting

ggplot() + 
  basemap_gglayer(bp_margit) + 
  geom_sf(data = bp_margit, color = "red") +
  coord_sf() +
  scale_fill_identity()