Disclaimer: The contents of this document come from 8 Making maps with R of Geocomputation with R (Lovelace, Nowosad, & Muenchow, 2019). This document is prepared for CP6521 Advanced GIS, a graduate-level city planning elective course at Georgia Tech in Spring 2019. For any question, contact the instructor, Yongsung Lee, Ph.D. via yongsung.lee(at)gatech.edu.
This document is also published on RPubs.
library(tidyverse)
library(sf)
library(raster)
library(spData)  
library(spDataLarge)

library(tmap)    # for static and interactive maps
library(leaflet) # for interactive maps
library(mapview) # for interactive maps
library(ggplot2) # tidyverse vis package
library(shiny)   # for web applications

Amateur-looking maps can undermine your audience’s ability to understand important information and weaken the presentation of a professional data investigation.

8.2 Static maps

8.2.1 tmap basics

‘grammar of graphics’

  • a separation between the input data and the aesthetics (how data are visualised)
  • each input dataset can be ‘mapped’ in a range of different ways,
  • which include location on the map (defined by data’s geometry), color, and other visual variables.
# Add fill layer to nz shape
tm_shape(nz) +
  tm_fill() 
# Add border layer to nz shape
tm_shape(nz) +
  tm_borders() 
# Add fill and border layers to nz shape
tm_shape(nz) +
  tm_fill() +
  tm_borders()
  • the common task of adding new layers is undertaken by the addition operator +, followed by tm_*()
  • the asterisk (*) refers to a wide range of layer types which have self-explanatory names including fill, borders (demonstrated above), bubbles, text and raster (see tmap-element for a full list)
8.2.2 Map objects
  • With tmap, we are able to store objects representing maps
map_nz = tm_shape(nz) + tm_polygons()
class(map_nz)
#> [1] "tmap"
  • map_nz can be plotted later:
    • by adding additional layers (as shown below) or
    • simply running map_nz in the console, which is equivalent to print(map_nz).
map_nz1 = map_nz +
  tm_shape(nz_elev) + tm_raster(alpha = 0.7) # alpha makes the layer semi-transparent
  • New shapes can be added with + tm_shape(new_obj).
  • In this case new_obj represents a new spatial object to be plotted on top of preceding layers.
  • When a new shape is added in this way, all subsequent aesthetic functions refer to it, until another new shape is added.
  • There is no limit to the number of layers or shapes that can be added to tmap objects.
  • Also, the same shape can even be used multiple times.
nz_water = st_union(nz) %>% st_buffer(22200) %>% 
  st_cast(to = "LINESTRING")
map_nz2 = map_nz1 +
  tm_shape(nz_water) + tm_lines()
map_nz3 = map_nz2 +
  tm_shape(nz_height) + tm_dots()
  • multiple map objects can be arranged in a single ‘metaplot’ with tmap_arrange().
tmap_arrange(map_nz1, map_nz2, map_nz3)
8.2.3 Aesthetics
  • There are two main types of map aesthetics:
    • those that change with the data
    • those that are constant
  • Unlike ggplot2, tmap accepts aesthetic arguments that are either variable fields (based on column names) or constant values
  • The most commonly used aesthetics for fill and border layers include color, transparency, line width and line type, set with col, alpha, lwd, and lty arguments, respectively
ma1 = tm_shape(nz) + tm_fill(col = "red")
ma2 = tm_shape(nz) + tm_fill(col = "red", alpha = 0.3)
ma3 = tm_shape(nz) + tm_borders(col = "blue")
ma4 = tm_shape(nz) + tm_borders(lwd = 3)
ma5 = tm_shape(nz) + tm_borders(lty = 2)
ma6 = tm_shape(nz) + tm_fill(col = "red", alpha = 0.3) +
  tm_borders(col = "blue", lwd = 3, lty = 2)
tmap_arrange(ma1, ma2, ma3, ma4, ma5, ma6)
  • Unlike the base R code below, tmap aesthetic arguments will not accept a numeric vector
plot(st_geometry(nz), col = nz$Land_area)  # works
tm_shape(nz) + tm_fill(col = nz$Land_area) # fails
#> Error: Fill argument neither colors nor valid variable name(s)
  • Instead, col (and other aesthetics that can vary such as lwd for line layers and size for point layers) requires a character string naming an attribute associated with the geometry to be plotted
tm_shape(nz) + tm_fill(col = "Land_area")

  • title sets the title of the associated legend
legend_title = expression("Area (km"^2*")") # expression() helps create superscript text
map_nza = tm_shape(nz) +
  tm_fill(col = "Land_area", title = legend_title) + tm_borders()
8.2.4 Color settings
  • The default setting uses ‘pretty’ breaks, described below.
  • breaks allows you to manually set the breaks.
  • n sets the number of bins into which numeric variables are categorized.
  • palette defines the color scheme, for example BuGn.
tm_shape(nz) + tm_polygons(col = "Median_income")
breaks = c(0, 3, 4, 5) * 10000
tm_shape(nz) + tm_polygons(col = "Median_income", breaks = breaks)
tm_shape(nz) + tm_polygons(col = "Median_income", n = 10)
tm_shape(nz) + tm_polygons(col = "Median_income", palette = "BuGn")
  • tmap also allows users to specify algorithms to automatically create breaks with the style argument
    • style = pretty, the default setting, rounds breaks into whole numbers where possible and spaces them evenly.
    • style = equal divides input values into bins of equal range, and is appropriate for variables with a uniform distribution (not recommended for variables with a skewed distribution as the resulting map may end-up having little color diversity).
    • style = quantile ensures the same number of observations fall into each category (with the potential down side that bin ranges can vary widely).
    • style = jenks identifies groups of similar values in the data and maximizes the differences between categories.
    • style = cont (and order) present a large number of colors over continuous color field, and are particularly suited for continuous rasters (order can help visualize skewed distributions).
    • style = cat was designed to represent categorical values and assures that each category receives a unique color.
  • palette defines the color ranges associated with the bins and determined by the breaks, n, and style arguments

  • add a - as prefix to reverse the palette order

  • Categorical palettes consist of easily distinguishable colors and are most appropriate for categorical data without any particular order such as state names or land cover classes
  • Sequential palettes follow a gradient, for example from light to dark colors (light colors tend to represent lower values), and are appropriate for continuous (numeric) variables.

  • Sequential palettes can be single (Blues go from light to dark blue, for example) or multi-color/hue (YlOrBr is gradient from light yellow to brown via orange, for example)

tm_shape(nz) + tm_polygons("Population", palette = "Blues")
tm_shape(nz) + tm_polygons("Population", palette = "YlOrBr")
  • Diverging palettes, typically range between three distinct colors (purple-white-green below) and are usually created by joining two single-color sequential palettes with the darker colors at each end.
  • These are used to visualize the difference from an important reference point, e.g., a certain temperature, the median household income or the mean probability for a drought event. The reference point’s value can be adjusted in tmap using the midpoint argument

Two important principles

  • Perceptibility
    • colors on maps should match our perception (e.g., green - vegetation or lowlands, blue - water or cool)
    • easy to understand to effectively convey information
    • e.g., it should be clear which values are lower and which are higher, and colors should change gradually
  • Accessibility (e.g, don’t use rainbow colors for this reason)
    • use colorblind friendly palettes as often as possible
8.2.5 Layouts
  • The map layout refers to the combination of all map elements into a cohesive map.
map_nz + 
  tm_compass(type = "8star", position = c("left", "top")) + # north arrows 
  tm_scale_bar(breaks = c(0, 100, 200), size = 1)           # scale bars
# see args(tm_layout) or ?tm_layout
map_nz + tm_layout(title = "New Zealand")
map_nz + tm_layout(scale = 5)
map_nz + tm_layout(bg.color = "lightblue")
map_nz + tm_layout(frame = FALSE)
  • Frame width (frame.lwd) and an option to allow double lines (frame.double.line).
  • Margin settings including outer.margin and inner.margin.
  • Font settings controlled by fontface and fontfamily.
  • Legend settings including binary options such as legend.show (whether or not to show the legend) legend.only (omit the map) and legend.outside (should the legend go outside the map?), as well as multiple choice settings such as legend.position.
  • Default colors of aesthetic layers (aes.color), map attributes such as the frame (attr.color).

  • Color settings controlling sepia.intensity (how yellowy the map looks) and saturation (a color-grayscale).

  • Beyond the low-level control over layouts and colors, tmap also offers high-level styles, using the tm_style() function

map_nza + tm_style("bw")
map_nza + tm_style("classic")
map_nza + tm_style("cobalt")
map_nza + tm_style("col_blind")
8.2.6 Faceted maps
  • Facets enable the visualization of how spatial relationships change with respect to another variable, such as time.
  • All individual facets in a faceted map contain the same geometry data repeated multiple times, once for each column in the attribute data.
  • However, facets can also represent shifting geometries such as the evolution of a point pattern over time.
urb_1970_2030 = urban_agglomerations %>% 
  filter(year %in% c(1970, 1990, 2010, 2030))
tm_shape(world) + tm_polygons() + 
  tm_shape(urb_1970_2030) + tm_symbols(col = "black", border.col = "white",
                                       size = "population_millions") +
  tm_facets(by = "year", nrow = 2, free.coords = FALSE)
  • Shapes that do not have a facet variable are repeated (the countries in world in this case).
  • The by argument which varies depending on a variable (year in this case).
  • nrow/ncol setting specifying the number of rows and columns that facets should be arranged into.
  • The free.coords-parameter specifying if each map has its own bounding box.
8.2.7 Inset maps
  • An inset map is a smaller map rendered within or next to the main map.

  • It could serve many different purposes, including providing a context or bringing some non-contiguous regions closer to ease their comparison.

  • In the example below, the first step is to define the area of interest, which can be done by creating a new spatial object, nz_region.

nz_region = 
  st_bbox(c(xmin = 1340000, xmax = 1450000,
            ymin = 5130000, ymax = 5210000),
          crs = st_crs(nz_height)) %>%
  st_as_sfc()
  • In the second step, we create a base map showing the New Zealand’s Southern Alps area.
nz_height_map = 
  tm_shape(nz_elev, bbox = nz_region) +
  tm_raster(style = "cont", palette = "YlGn", legend.show = TRUE) +
  tm_shape(nz_height) + 
  tm_symbols(shape = 2, col = "red", size = 1) +
  tm_scale_bar(position = c("left", "bottom"))
  • The third step consists of the inset map creation.
nz_map = 
  tm_shape(nz) + 
  tm_polygons() +
  tm_shape(nz_height) + 
  tm_symbols(shape = 2, col = "red", size = 0.1) + 
  tm_shape(nz_region) + 
  tm_borders(lwd = 3)
  • Finally, we combine the two maps using the function viewport() from the grid package
    • the first arguments of which specify the center location (x and y) and a size (width and height) of the inset map.
library(grid)
nz_height_map
print(nz_map, vp = viewport(0.8, 0.27, width = 0.5, height = 0.5))
  • Inset maps are also used to create one map of non-contiguous areas.
  • Probably, the most often used example is a map of the United States, which consists of the contiguous United States, Hawaii and Alaska
# first create a map for the contiguous US with an appropriate projection 
us_states_map = tm_shape(us_states, projection = 2163) + tm_polygons() + 
  tm_layout(frame = FALSE)

# create two states on separate maps. 
hawaii_map = tm_shape(hawaii) + tm_polygons() + 
  tm_layout(title = "Hawaii", frame = FALSE, bg.color = NA, 
            title.position = c("LEFT", "BOTTOM"))
alaska_map = tm_shape(alaska) + tm_polygons() + 
  tm_layout(title = "Alaska", frame = FALSE, bg.color = NA)

# now, combine  and arrange these three maps 
us_states_map
print(hawaii_map, vp = grid::viewport(0.35, 0.1, width = 0.2, height = 0.1))
print(alaska_map, vp = grid::viewport(0.15, 0.15, width = 0.3, height = 0.3))
  • Alternative mapping scripts can be found here.

8.3 Animated maps

Skip

8.4 Interaction maps

tmap_mode("view")
map_nz 
map_nz + tm_basemap(server = "OpenTopoMap")
world_coffee = left_join(world, coffee_data, by = "name_long")
facets = c("coffee_production_2016", "coffee_production_2017")
tm_shape(world_coffee) + 
  tm_polygons(facets) + 
  tm_facets(nrow = 1, sync = TRUE)
# If you are not proficient with `tmap`, the quickest way to create interactive maps may be with `mapview`. 
mapview::mapview(nz)
trails %>%
  st_transform(st_crs(franconia)) %>%
  st_intersection(franconia[franconia$district == "Oberfranken", ]) %>%
  st_collection_extract("LINE") %>%
  mapview(color = "red", lwd = 3, layer.name = "trails") +
  mapview(franconia, zcol = "district", burst = TRUE) + # 'burst' datasets into multiple layers
  breweries # addition of multiple layers with + followed by the name of a geographic object

Other options:

leaflet
- the most mature and widely used interactive mapping package in R
- provides a relatively low-level interface to the Leaflet JavaScript library
- many of its arguments can be understood by reading the documentation of the original JavaScript library (check Leaflet API reference)

pal = colorNumeric("RdYlBu", domain = cycle_hire$nbikes)
leaflet(data = cycle_hire) %>% 
  addProviderTiles(providers$OpenStreetMap.BlackAndWhite) %>% 
  addCircles(col = ~pal(nbikes), opacity = 0.9) %>% 
  addPolygons(data = lnd, fill = FALSE) %>% 
  addLegend(pal = pal, values = ~nbikes) %>% 
  setView(lng = -0.1, 51.5, zoom = 12) %>% 
  addMiniMap()

8.5 Mapping applications

A few limitaions of interactive mapping approaches above

# this code chunk is not run here. Try in your R session. 
library(shiny)    # for shiny apps
library(leaflet)  # renderLeaflet function
library(spData)   # loads the world dataset 
ui = fluidPage(
  sliderInput(inputId = "life", "Life expectancy", 49, 84, value = 80),
      leafletOutput(outputId = "map")
  )
server = function(input, output) {
  output$map = renderLeaflet({
    leaflet() %>% addProviderTiles("OpenStreetMap.BlackAndWhite") %>%
      addPolygons(data = world[world$lifeExp < input$life, ])})
}
shinyApp(ui, server)
How the above scripts work

The user interface (ui) of lifeApp is created by fluidPage(). This contains input and output ‘widgets’ — in this case, a sliderInput() (many other *Input() functions are available) and a leafletOutput(). These are arranged row-wise by default, explaining why the slider interface is placed directly above the map (see ?column for adding content column-wise).

The server side (server) is a function with input and output arguments. output is a list of objects containing elements generated by render*() function — renderLeaflet() which in this example generates output$map. Input elements such as input$life referred to in the server must relate to elements that exist in the ui — defined by inputId = "life" in the code above.

The function shinyApp() combines both the ui and server elements and serves the results interactively via a new R process. When you move the slider in the map, you are actually causing R code to re-run, although this is hidden from view in the user interface.

Resources:

8.6 Other mapping packages

Skip