Making maps with ggplot2 and sf

Recently, the newest version of the popular ggplot2 graphics package was announced, and it has some nifty mapping features that I was keen to try out (read more here). Mainly, I was interested in the support for sf, or “simple features”, objects. This class of objects were created as part of a wider R package designed to make mapping and spatial analyses far easier.

The latest update of ggplot2 not only makes plotting from sf objects trivial, but also means that some quite nice map figures can be made with relatively little effort, as you’ll hopefully see below.

# first you'll need to update your version of ggplot to the latest version
# install.packages("ggplot2")
library(ggplot2)

# now lets load some other packages that we'll need to load and manipulate
# spatial data in R
library(mapdata)
library(sf)
library(lwgeom)

Assuming you managed to get those packages installed and loaded ok (if you didn’t, you almost certainly have encountered some dependency issues), we can now start playing with some data. To begin, let’s load a world map to play with.

mapBase <- map("worldHires", fill = T, plot = F)

# now we need to coerce it to an "sf" object, and fix any
mapBase <- st_as_sf(mapBase)

# now let's try cropping it to a region of Europe
cropMap <- st_crop(mapBase, xmin = -15, xmax = 30, ymin = 30, ymax = 60)
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## Error in CPL_geos_op2(op, st_geometry(x), st_geometry(y)): Evaluation error: TopologyException: Input geom 0 is invalid: Self-intersection at or near point -6.7722639937312801 41.090800036327813 at -6.7722639937312801 41.090800036327813.
# note we get an error message about Self-intersection...
# we can fix this using the lwgeom library...
mapBase <- st_make_valid(mapBase)
cropMap <- st_crop(mapBase, xmin = -15, xmax = 30, ymin = 30, ymax = 60)
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## Warning: attribute variables are assumed to be spatially constant
## throughout all geometries
# not it works fine...

Now that we have a valid sf object to plot, we can start using the geom_sf function to start making some nice maps.

# basic map to start with
ggplot(cropMap) + geom_sf()

plot of chunk unnamed-chunk-3

Cool right? Notice a few neat things? Firstly, ggplot draws a nice graticule for you, complete with the correct “degrees” symbol and N, S, E, or W to denote the correct hemisphere. Secondly, go ahead and resize the map… What do you notice? Hopefully, you’ll see that the aspect ratio of the map has been fixed so that your map always projects correctly. This is a great feature and saves you some time faffing around trying to manually correct the aspect ratio.

Having explored some basic features, we can now start to tailor our map and add features to it using other packages, or data.

# let's say we want to highlight the location of Spain, clarify the axes,
# and remove those annoying grid lines
spainMap <- ggplot(cropMap,
    aes(fill = factor(ifelse(cropMap$ID == "Spain", 1, 2)))) +
  geom_sf() +
  labs(x = "Longitude", y = "Latitude", fill = "") +
  scale_fill_manual(values = c("darkgrey", "lightgrey"),
    labels = c("Spain", "Not Spain")) +
  theme_bw() +
  theme(panel.grid = element_line(colour = "transparent"),
    axis.title = element_text(size = 18),
    axis.text = element_text(size = 16),
    legend.text = element_text(size = 14),
    legend.title = element_text(size = 14))

spainMap

plot of chunk unnamed-chunk-4

Pretty easy so far right? But what if we want to add some points to the map, for example to highlight sampling locations? Well, this is really easy too!

# simulate some random coordinates within the limits of our map
locations <- data.frame(lat = runif(25, 30, 60), long = runif(25, -15, 30))

# now convert this dataframe into an sf dataframe
sfPoints <- st_as_sf(locations, coords = c("long", "lat"), crs = 4326)

# and simply add another geom_sf layer to the plot to include the points
pointsMap <- ggplot() +
  geom_sf(data = cropMap,
    aes(fill = factor(ifelse(cropMap$ID == "Spain", 1, 2)))) +
  geom_sf(data = sfPoints, col = "red", size = 3) +
  labs(x = "Longitude", y = "Latitude", fill = "") +
  scale_fill_manual(values = c("darkgrey", "lightgrey"),
    labels = c("Spain", "Not Spain")) +
  theme_bw() +
  theme(panel.grid = element_line(colour = "transparent"),
    axis.title = element_text(size = 18),
    axis.text = element_text(size = 16),
    legend.text = element_text(size = 14),
    legend.title = element_text(size = 14))

pointsMap

plot of chunk unnamed-chunk-5

Again, super simple. Now, I’d like to create a smaller map, that just shows Spain and any points inside it, so I can make a nice panel figure with the two maps side by side. Also, don’t forget a North arrow and scale bar, which we can add using the ggsn package.

library(ggsn)

# work out which points are in the polygon of interest
spainPoints <- st_join(sfPoints, cropMap[cropMap$ID == "Spain", ],
  join = st_intersects)
## although coordinates are longitude/latitude, st_intersects assumes that they are planar
spainZoom <- ggplot() +
  geom_sf(data = cropMap[cropMap$ID == "Spain", ]) +
  scalebar(cropMap[cropMap$ID == "Spain", ], dist = 200, dd2km = TRUE,
    st.dist = 0.05) +
  north(cropMap[cropMap$ID == "Spain", ]) +
  geom_sf(data = spainPoints[spainPoints$ID == "Spain", ], colour = "red",
    size = 3) +
  theme_void() +
  theme(panel.grid = element_line(colour = "transparent"))

spainZoom

plot of chunk unnamed-chunk-6

So, now we have our two maps sorted, we can arrange them side by side using the cowplot package.

library(cowplot)

plot_grid(spainZoom, pointsMap, labels = "AUTO", rel_widths = c(0.6, 1))

plot of chunk unnamed-chunk-7

There you have it, combining the sf and newest ggplot2 packages allows you quickly and easily make some neat looking map figures!

Advertisements

4 thoughts on “Making maps with ggplot2 and sf

  1. Nice post, I got an error trying to reproduce the zomm on Spain map before cowplot:

    Error in gcDestination(lon = x, lat = y, bearing = 90 * direction, dist = dist, : argument “model” is missing, with no default

    Regards,
    Walter

    Like

    • Hi Walter,
      Many thanks for your comment, and sorry you are getting an error. Could you try upgrading your version of the “ggsn” package to the development version (see here: https://github.com/oswaldosantos/ggsn)

      The newest version will auto-detect the coordinate reference system used by the other sf objects being plotted, and will use this to add the scale bar. The older version requires you to specify this, and I think this is what may be causing your error.

      Hope this helps, and cheers for checking out my blog!
      Dave

      Like

  2. Very interesting plot and easy to build. Really good job! Here’s the question. What about if we want to add regions for any country? For instance, keep with Spain. We need to add the 17 regions and moreover, inside the regions all the provinces. Would it be possible? Another question, is it possible to add points for specific coordinates? I mean, I want to point the city of Barcelona. Do I need to know the longitude and latitude to locate it? Thanks in advance and once again, very useful map!!!

    Like

    • Hi Jose,
      Thanks for the kind comments, glad it is useful. You can certainly add boundary lines for regions within a country – you just have to find this data somewhere (e.g. as a shapefile, geoJSON etc). You can then use the sf package to import it into R, and add it as a geom_sf layer to your plot.

      To answer your 2nd question, you can add points for specific coordinates similarly to how I added them. If you know the coordinates, it’s a simple matter of converting them to an sf object (see here https://r-spatial.github.io/sf/reference/st_transform.html), or if you don’t know them, try a geo-coding approach (read more here: https://www.jessesadler.com/post/geocoding-with-r/), before then converting coordinates into sf format.

      Hope this is useful, and thanks for checking out the blog!
      Dave

      Like

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google+ photo

You are commenting using your Google+ account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s