MAR 536 Lab 8: plotting spatial data

2023-03-29

Today: Advanced plotting

• Plot review
• Layouts & multi-panel plots
• Mapping
• Spatial analyses

3 Things required for ggplots:

1. Call to ggplot()
2. Aesthetics using aes()
3. Geometries e.g. geom_point(), geom_bar(), geom_line(), ect.

Multi-panel plot options

ggplot

• facet_grid()
• grid and gridExtra packages
• patchwork package
• par()
• lattice package other examples
• ggridges
• cowplot

facet_grid() example: Laengelmavesi data (lab 2)

facet_grid()

ggplot(data, aes(y=length, col=species)) +
  geom_boxplot() +
  facet_grid(cols = vars(species))

facet_grid()

ggplot(data, aes(length, col=species)) +
  geom_histogram(aes(fill=species)) +
  facet_grid(rows = vars(species))

grid and gridExtra examples

The gridExtra package arranges “grobs”

grob = graphical object

gridExtra functions can also arrange gtables and ggplot objects

library(gridExtra)
library(grid)

Some gapminder example plots to work with:

library(gapminder)
library(gghighlight)

plot1 <- filter(gapminder, year == 2007) |>
  ggplot() + 
  aes(x = lifeExp) + 
  geom_histogram(bins = 10, col = "white")
plot1

continents with low life expectancy

plot2 <- ggplot(gapminder) +
  aes(x=year, 
      y=lifeExp, 
      group=country, 
      color=country) +
  geom_line() +
  gghighlight(min(lifeExp) < 30) +
  facet_wrap(~continent)+
  theme_minimal() +
  theme(legend.position = "none")
plot2

distribution of life expectancy over years for two countries

plot3 <- gapminder |>
  filter(country %in% c("Peru","Iceland")) |>
ggplot() +
  aes(x=lifeExp, color = country) +
  geom_density()
plot3

grid and gridExtra examples continued…

Passing plots to grid.arrange() and specifying either the number of rows or columns gives a simple layout.

grid.arrange(plot1, plot2, plot3, nrow = 1)

grid and gridExtra examples continued…

Grobs may also be placed in a list and arranged using customized formats using the argument layout_matrix.

# Create a new grob
newgrob <- textGrob("This is a \n text grob")
groblist <- list(plot1, plot2, plot3, newgrob)
# Create a custom layout
grid.arrange(grobs = groblist,
             widths = c(2,2,4,1),
             layout_matrix = rbind(c(2,2,1,1),
                                   c(2,2,3,4)))

Patchwork

The patchwork package provides a shorthand method to plot multiple ggplot objects together

library(patchwork)
plot1 + plot2 + plot3

Patchwork example continued…

To structure plot layouts further use the plot_layout() function

plot1 + {plot2 + plot3 + plot_layout(ncol=2, widths = c(3,1))} + plot_layout(ncol = 1)

Note that the {} indicate a nested plot

Patchwork example continued…

This package also uses | to indicate plots adjacent to one another and / to indicate vertical stacking

(plot1 | plot3) / plot2

For more examples see: https://gotellilab.github.io/GotelliLabMeetingHacks/NickGotelli/ggplotPatchwork.html

Exercise 1

1. Use patchwork or grid+gridExtra to create a 4-panel plot grid with the following characteristics:

• space for 4 ggplots
• 2 columns & 3 rows
• the first column should be twice as wide as the second column
• the first plot should appear in the entire first column
• plots 2-4 should fill the second column

Mapping

There are n (where n is large) ways to produce maps using R.

Many of these are beyond this course and delve into the dark abstractedly projected rabbithole that is cartography.

However, we can produce very nice looking maps in R quite easily.
(because learning GIS can be hard)

Using {ggplot2} and the {sf} package.

We can also use {sf} to do lots of spatial data operations & analyses, all within R.

Country coastlines

The package {rnaturalearth} provides a map of countries of the entire world.

library("rnaturalearth")
library("rnaturalearthdata")
#library("rnaturalearthhires") # for scale = large
world <- ne_countries(scale = "medium", returnclass = "sf")
class(world)
## [1] "sf"         "data.frame"

simple map

ggplot(data = world) +
    geom_sf()

map color

ggplot(data = world) + 
    geom_sf(color = "black", fill = "lightgreen")

chloropleth maps

plotting a data variable with our geometries

ggplot(data = world) +
    geom_sf(aes(fill = pop_est)) +
    scale_fill_viridis_c(option = "plasma", trans = "sqrt")

Projection and extent

coord_sf() deals with the coordinate system Used to change the map projections, and the extent. e.g.

ggplot(data = world) +
    geom_sf() +
    coord_sf(xlim = c(-102.15, -74.12), ylim = c(7.65, 33.97), expand = FALSE)

census data example

library(tidycensus)
#options(tigris_use_cache = TRUE)

age_data <- get_acs(
  geography = "tract",
  variables = "B19013_001", #"B01002_001",
  state = "MA",
  geometry = TRUE
)
## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |=                                                                     |   2%
  |                                                                            
  |==                                                                    |   2%
  |                                                                            
  |==                                                                    |   3%
  |                                                                            
  |===                                                                   |   4%
  |                                                                            
  |===                                                                   |   5%
  |                                                                            
  |====                                                                  |   6%
  |                                                                            
  |=====                                                                 |   7%
  |                                                                            
  |======                                                                |   9%
  |                                                                            
  |=======                                                               |  10%
  |                                                                            
  |========                                                              |  12%
  |                                                                            
  |=========                                                             |  12%
  |                                                                            
  |=========                                                             |  13%
  |                                                                            
  |===========                                                           |  16%
  |                                                                            
  |============                                                          |  17%
  |                                                                            
  |=============                                                         |  19%
  |                                                                            
  |===============                                                       |  22%
  |                                                                            
  |================                                                      |  23%
  |                                                                            
  |=================                                                     |  24%
  |                                                                            
  |=================                                                     |  25%
  |                                                                            
  |===================                                                   |  28%
  |                                                                            
  |=====================                                                 |  29%
  |                                                                            
  |=======================                                               |  33%
  |                                                                            
  |=========================                                             |  35%
  |                                                                            
  |==========================                                            |  37%
  |                                                                            
  |============================                                          |  40%
  |                                                                            
  |=============================                                         |  41%
  |                                                                            
  |================================                                      |  46%
  |                                                                            
  |=================================                                     |  47%
  |                                                                            
  |==========================================                            |  60%
  |                                                                            
  |===========================================                           |  62%
  |                                                                            
  |==============================================                        |  66%
  |                                                                            
  |==================================================                    |  72%
  |                                                                            
  |===================================================                   |  72%
  |                                                                            
  |====================================================                  |  75%
  |                                                                            
  |======================================================                |  77%
  |                                                                            
  |========================================================              |  81%
  |                                                                            
  |===========================================================           |  84%
  |                                                                            
  |============================================================          |  86%
  |                                                                            
  |================================================================      |  91%
  |                                                                            
  |======================================================================| 100%
ggplot(age_data, aes(fill = estimate)) + 
  geom_sf(color = NA) + 
  theme_void() + 
  scale_fill_viridis_c(option = "magma") + 
  labs(title = "Massachusetts median income",
       subtitle = "2017-2021 ACS",
       fill = "Estimate")

looking at the sf object

age_data
## Simple feature collection with 1620 features and 5 fields (with 4 geometries empty)
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -73.50814 ymin: 41.23796 xmax: -69.92839 ymax: 42.88659
## Geodetic CRS:  NAD83
## First 10 features:
##          GEOID                                                  NAME   variable
## 1  25025050200       Census Tract 502, Suffolk County, Massachusetts B19013_001
## 2  25023540103  Census Tract 5401.03, Plymouth County, Massachusetts B19013_001
## 3  25025100500      Census Tract 1005, Suffolk County, Massachusetts B19013_001
## 4  25005652200      Census Tract 6522, Bristol County, Massachusetts B19013_001
## 5  25023505104  Census Tract 5051.04, Plymouth County, Massachusetts B19013_001
## 6  25017330201 Census Tract 3302.01, Middlesex County, Massachusetts B19013_001
## 7  25021404100      Census Tract 4041, Norfolk County, Massachusetts B19013_001
## 8  25023511000     Census Tract 5110, Plymouth County, Massachusetts B19013_001
## 9  25025091200       Census Tract 912, Suffolk County, Massachusetts B19013_001
## 10 25025082000       Census Tract 820, Suffolk County, Massachusetts B19013_001
##    estimate   moe                       geometry
## 1     71402 10472 MULTIPOLYGON (((-71.04003 4...
## 2    167604 31302 MULTIPOLYGON (((-71.03657 4...
## 3     49644 23048 MULTIPOLYGON (((-71.0782 42...
## 4     69574 12716 MULTIPOLYGON (((-70.95382 4...
## 5    125833 37004 MULTIPOLYGON (((-70.75747 4...
## 6    136444 44601 MULTIPOLYGON (((-71.12136 4...
## 7    250001    NA MULTIPOLYGON (((-71.27219 4...
## 8     61875 10620 MULTIPOLYGON (((-71.01678 4...
## 9     71125 33529 MULTIPOLYGON (((-71.06743 4...
## 10    69041 27243 MULTIPOLYGON (((-71.08833 4...

just for Bristol County

age_data <- get_acs(
  geography = "tract",
  variables = "B19013_001",  #for median age = "B01002_001",
  state = "MA",
  county = "Bristol",
  geometry = TRUE
)

ggplot(age_data, aes(fill = estimate)) + 
  geom_sf(color = NA) + 
  theme_void() + 
  scale_fill_viridis_c(option = "magma") + 
  labs(title = "Bristol County median income",
       subtitle = "2017-2021 ACS",
       fill = "Estimate")

Exercise 2

Use ggplot to create a map of Cape Cod with the following features:

• Longitudes should range from \(-71^{\circ}\text{W}\) to \(-69^{\circ}\text{W}\)
• Latitudes should range from \(41.25^{\circ}\text{N}\) to \(43^{\circ}\text{N}\)
• Label axes
• Color land
• Add points indicating the locations of Woods Hole, Chatham, and Provincetown
• BONUS Change the map projection

For coastline use:

library("rnaturalearthhires") # for scale = large
world <- ne_countries(scale = "large", returnclass = "sf")