Objectives

There is no single package for mapping data in R. I’ve found that different packages and functions perform best for different types of maps. In this tutorial, I will walk you through several packages, data sources, and plotting tools. I realize this may be a bit more confusing than focusing on a single package, but after quite a bit of thought on the matter, I decided the most useful thing would be to introduce you to the best tools I’ve found to date and let you decide which package/toolkit is the most appropriate for your visualization needs. With this in mind, our objectives for this tutorial are as follows:


My assumptions


Set-up

Make sure you’ve installed and loaded the following:

library(tidyverse)
library(sf)
library(tmap) # lots to install here, sorry, but it's worth it.
# before installing tmap, do Tools --> Check for package updates and install all updates
library(spData)

In what follows, I’m going to walk you through the basics of map-making in R. We’ll start by learning how to access some great basemaps (not to be used for analysis, just for visualization), the slowly add shapefiles and rasters. We’ll learn how to make some common maps (i.e. choropleth) in R and walk through a couple of packages that simplify the whole map-making process. We won’t have time to get to interactive maps, which is one of the cooler things you can do with R, but if you’re interested in this, then you’ll need to learn more about Leaflet and Shiny. You can see a few examples of interactive maps made in R here.


Mapping your own data

You can download each of the data we’ll be working with here:

data("us_states") # from spData package
ga <- us_states %>% filter(NAME == "Georgia")
roads <- st_read("./data/tl_2016_us_primaryroads.shp") 
## Reading layer `tl_2016_us_primaryroads' from data source `C:\Users\eburchf\OneDrive - Emory University\Teaching\ENVS 270 Environmental Data Science\03_ESDA\data\tl_2016_us_primaryroads.shp' using driver `ESRI Shapefile'
## Simple feature collection with 12509 features and 4 fields
## geometry type:  LINESTRING
## dimension:      XY
## bbox:           xmin: -170.8313 ymin: -14.3129 xmax: -65.64866 ymax: 49.00209
## geographic CRS: NAD83
counties <- st_read("./data/Counties%20Georgia.shp", quiet=T) # add quiet argument to supress outcome

ANYTIME YOU ARE WORKING WITH MULTIPLE SPATIAL DATASETS, FIRST CHECK TO MAKE SURE THEY ARE ALL IN THE SAME PROJECTION!

In our case, each dataset is in a different projection. Let’s change them all to a good projection for Georgia that I found through some Googlin’.

proj <- "+proj=lcc +lat_1=31.41666666666667 +lat_2=34.28333333333333 +lat_0=0 +lon_0=-83.5 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +to_meter=0.3048006096012192 +no_defs"
roads <- st_transform(roads, proj)
ga <- st_transform(ga, proj)
counties <- st_transform(counties, proj)

roads <- st_crop(roads, ga)

Now that our data are in the same projection, we can quickly visualize all three data in a single plot. Since they are in the same projection, this also allows us to use the st_crop() function to crop the national road dataset to the boundaries of Georgia! If they weren’t in the same projection, you’d need to reproject them to a single projection before doing any visualization.

plot(st_geometry(ga))
plot(st_geometry(counties), add = T, col = "red")
plot(st_geometry(roads), add = T, col = "blue")

You can actually get fairly far with the base R plot() function. The code gets dense quite quickly and the visualizations are nothing like what tmap and ggplot can produce (see below). I recommend using plot() only to quickly visualize your data. Data visualizations you intend to publish or share should be made with one of the packages I describe below. If you’re interested in learning more about using plot(), I recommend Chapter 3 of Applied Spatial Data Analysis in R by Roger Bivand and colleagues. Below, I’ve inserted a table from p. 68 of this text that reviews useful arguments to be passed to the plot() function when working with spatial objects. If you want to learn even more about plot(), go here and here and/or here and/or here.


ggplot2

Though it wasn’t originally designed for mapping, you can do quite a bit of mapping with our friend ggplot2. Check this out:

ggplot(ga) +
  geom_sf()

From there, you can do all of the usual ggplot tuning to make things look a bit better:

ggplot(ga) +
  geom_sf(fill = "beige") +
  labs(title = "The fine state of Georgia") +
  theme_minimal()

What if we want to add multiple objects in one map?

ggplot() +
  geom_sf(data = counties, aes(fill = totpop10)) +
  geom_sf(data = roads, color = "orange") +
  labs(title = "All roads lead to the ATL",
       fill = "Population") +
  theme_minimal()

Don’t forget to change the legend title! More about manipulating legends with ggplot2 found here. If you look at this code all at once, it may seem overwhelming, but the big advantage of ggplot is that it’s modular. You can - and should - add pieces to your visualization one step at a time.

This post has some other useful tricks including using the ggspatial package to add a scale bar and compass:

library(ggspatial)
ggplot() +
  geom_sf(data = counties, aes(fill = totpop10)) +
  geom_sf(data = roads, color = "orange") +
  labs(title = "All roads lead to the ATL",
       fill = "Population") +
  theme_minimal() +
  annotation_scale(location = "bl", width_hint = 0.5) +
  annotation_north_arrow(location = "tr", which_north = "true",
                         style = north_arrow_fancy_orienteering())


tmap

I’ve saved the best for last. The tmap package makes it easy to visualize sf objects in R. It’s far more elegant than the other options I’ve shown you and makes creating interactive maps very simple. It also makes all of the DOGTAILS work straightforward. I am going to walk through some of the tools in the tmap package as applied to our Utah dataset. The full package vignette is fantastic and found here.

Let’s start with the qtm() or quick thematic map tool. This quickly creates maps of any shape object (similar to the base R plot() function):

library(tmap)
tmap_mode("plot") # set to "view" for interactive visualizations, see below
qtm(roads)

We can use qtm() to build fast choropleth maps:

qtm(counties, fill = "totpop10", fill.title = "Total population")

Notice how without specification, qtm() creates a legend! As you bang your head against your desk expressing frustration over the last four packages we’ve gone through when, ultimately, you’ll probably just use this package, remember that knowledge is POWER! I wanted to introduce you to the MANY tools out there, you never know when you might get stuck and need one of them. qtm() has a number of cool arguments:

  • text allows you to specify text to go over each polygon/attribute in your shapefile (i.e. state names, abbreviations, whatever)
  • text.size controls text size.
  • format - the package comes with a few pre-defined map formats that are worth checking out. There are also a few predefined style options included.
  • fill.title titles the legend
  • text.root determines how text size is increased
  • fill.textNA allows you to specify the name of attributes with value NA

In addition to this awesome tool, tmap has a number of other tools for specific spatial data types. I’ve reproduced the full list from the tmap vignette here and recommend you go read this vignette if you plan on plotting with tmap.

Think of these as similar to the geo elements you add to a ggplot. Like ggplot, you can use tmap to add complexity and layers to a map. Hooray! Here’s an example with our Georgia data sets. First first element is tm_shape(), which specifies the shape object. Then we can add additional information. Here’s a very simple example:

tm_shape(counties) +
  tm_fill()

Very informative map. Basically here I tell tmap that I am working with the tracts shapefile (tm_shape(tracts)) and I just want to fill in the whole thing in the default color, a bleak shade of gray. We can do better than this:

tm_shape(counties) +
  tm_fill(col = "totpop10", title = "Population") +
  tm_compass() +
  tm_scale_bar()

Oooohhh! We can add more information to tm_fill() to specify exactly how we want to fill polygons. Here’s one example:

tm_shape(counties) +
  tm_fill(col = "totpop10", convert2density = T, style = "jenks", title = "Population") +
  tm_compass() +
  tm_scale_bar()

You know what to do to get more info about fill options (?tm_fill). Notice how the style we chose significantly changes our visualization (i.e. maps can lie). Here, I’ve used jenks, which minimizes each class’s average deviation from the class mean. I’ve also asked tmap to plot population density, i.e. totpop10 divided by the area of each polygon (which tmap calculates in the background… sniff sniff, so beautiful).

And to add the tract polygons we can add tm_borders() and specify the darkness of borders using alpha:

tm_shape(counties) +
  tm_fill(col = "totpop10", convert2density = T, style = "jenks", title = "Population (x100)") +
  tm_borders(alpha = 0.3) +
  tm_compass() +
  tm_scale_bar()

Did you notice that tmap makes adding our DOGTAILS painless? It’s a lot more painful in ggplot and spplot, which is why I didn’t show it in this tutorial (trust me, painful). Each tm_shape() and it’s specifications are referred to as a group. We can stack multiple groups (i.e. specific visualization of tm_shape(roads), tm_shape(counties), and tm_shape(ga)). A crazy thing about tmap is that each shapefile you intend to plot can have different projections and extents. tmap will automatically reproject stacked groups of data to the projection of the first group (counties in our case).

tm <- tm_shape(counties) +
  tm_fill(col = "totpop10", convert2density = T, style = "jenks", title = "Population (x100)") +
  tm_borders(alpha = 0.3) +
  tm_compass() +
  tm_scale_bar() +
tm_shape(roads) +
  tm_lines(col = "black", alpha = 0.4)

tm

You can alter the format, shape, and legend, among other things, but digging a bit more into the tmap elements and code. I strongly recommend you review the [tmap in a Nutshell Vignette](https://cran.r-project.org/web/packages/tmap/vignettes/tmap-nutshell.html) and the full package documentation that explains all of the cool things you can do. I’ll quickly show you how to fix our map to make the legend a bit more legible:

tm <- tm +
  tm_legend(text.size=1,
    title.size=1.2,
    position = c("left","bottom"), 
    bg.color = "white", 
    bg.alpha=.5, 
    frame="gray50", 
    height=.6)
tm

Want to save this beautiful map we made?

save_tmap(tm, "./fig/my_awesome_map.png", width = 800, height = 1000)

The last amazing, mindbogglingly cool thing we can do with tmap is quickly make interactive visualizations. I won’t have time to cover Leaflet and Shiny, the main tools for interactive mapping in R. If, however, you want to QUICKLY build a cool interactive map, tmap can help!

tmap_mode("view")
tm
# interactive map!

There are millions of (fairly easy) ways you can create beautiful maps with tmap. A few more examples can be found here:

  1. Interactive multivariate world map
  2. Fancy world map
  3. World facets
  4. US choropleth
  5. London crimes
  6. Adding insets for Hawaii and Alaska

Advanced mapping tools

…that we won’t cover, but that are very cool and that I encourage you to explore:


Additional resources