Plotting static base map underneath a sf object

Question

I'm trying to plot a static base map underneath my sf object (for print). When using ggmap I first get a lot of errors, then I can't seem to figure out how to chain the base map to my ggplot2 object with a geom_sf.

library(sf)
# devtools::install_github("tidyverse/ggplot2")
library(ggplot2)
library(ggmap) 

nc <- st_read(system.file("shape/nc.shp", package="sf"))
nc_map <- get_map(location = "North Carolina, NC", zoom = 7)

ggmap(nc_map)

nc_centers <- st_centroid(nc)

nc_centers %>%
  ggplot() +
  geom_sf(aes(color = SID79, size = BIR74),
          show.legend = "point") +
  coord_sf(datum = NA) +
  theme_minimal()

I also rather use the source = "osm" as style but those will always return '400 Bad Request'.

Is there maybe another good package for base maps?

Flintlock

leaflet is a good package imo - don't know though if you'll be able to complete your task with it. – Pyrophosphate

leaflet for interactive maps right? I'm looking for static print solutions. – Flintlock

You can use mapview::mapshot to save leaflet maps as a static file (jpeg, png etc.). Maybe that suits your needs – Gillian

Answer 1

You might consider reprojecting your data but the following code seems to work for me.
See here for an explanation about why you need inherit.aes = FALSE and see here for an alternative solution with base plots.

library(sf)
#> Linking to GEOS 3.5.1, GDAL 2.1.3, proj.4 4.9.2
# devtools::install_github("r-lib/rlang")
library(ggplot2)
library(ggmap) 

nc <- st_read(system.file("shape/nc.shp", package="sf"))
#> Reading layer `nc' from data source `/home/gilles/R/x86_64-pc-linux-gnu-library/3.4/sf/shape/nc.shp' using driver `ESRI Shapefile'
#> Simple feature collection with 100 features and 14 fields
#> geometry type:  MULTIPOLYGON
#> dimension:      XY
#> bbox:           xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
#> epsg (SRID):    4267
#> proj4string:    +proj=longlat +datum=NAD27 +no_defs
nc_map <- get_map(location = "North Carolina, NC", zoom = 7)
#> Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=North+Carolina,+NC&zoom=7&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
#> Information from URL : http://maps.googleapis.com/maps/api/geocode/json?address=North%20Carolina,%20NC&sensor=false
nc_centers <- st_centroid(nc)
#> Warning in st_centroid.sfc(st_geometry(x), of_largest_polygon =
#> of_largest_polygon): st_centroid does not give correct centroids for
#> longitude/latitude data

ggmap(nc_map) +
    geom_sf(data = nc_centers, 
            aes(color = SID79, size = BIR74),
            show.legend = "point", inherit.aes = FALSE) +
    coord_sf(datum = NA) +
    theme_minimal()
#> Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Created on 2018-04-03 by the reprex package (v0.2.0).

Sierrasiesser

@GillesThanks for the excellent solution. Could you please add text annotation too? I have tried geom_sf_label () but it seems to fail with this solution. – Botheration

Answer 2

I've been working on a package recently which someone may find useful. The ggmap answer now require an API key for use with google maps which adds an extra bit of faff.

basemapR also allows you a bit more flexibility in setting the extents for the basemap by using a bounding box.

#devtools::install_github('Chrisjb/basemapR')
library(basemapR)
library(sf)
library(ggplot2)


nc <- st_read(system.file("shape/nc.shp", package="sf"))

nc_centers <- st_centroid(nc)

# create bbox from our nc layer and expand it to include more area above/below
bbox <- expand_bbox(st_bbox(nc_centers), X = 0, Y = 150000)

ggplot() +
  base_map(bbox, increase_zoom = 2, basemap = 'google-terrain') +
  geom_sf(data = nc_centers, 
          aes(color = SID79, size = BIR74),
          show.legend = "point", inherit.aes = FALSE) +
  coord_sf(datum = NA,
           xlim = c(bbox['xmin'], bbox['xmax']),
           ylim = c(bbox['ymin'], bbox['ymax'])) +
  theme_minimal() +
  labs(caption = 'map data \uA9 2020 Google')

Remember to cite google maps in the caption or elsewhere in the map

Willmert

Answer 3

You can also use the package ggspatial which provides the "map tile" annotation layer.

ggplot(nc_centers) +
    annotation_map_tile(zoom = 7) +
    geom_sf(aes(color = SID79, size = BIR74),
            show.legend = "point", inherit.aes = FALSE) +
    coord_sf(datum = NA) +
    theme_minimal()

Mostly