Chapter 4 Part 2

Harold Nelson

10/10/2021

Setup

Get the libraries and re-create the data from Part 1.

library(sf)           # Objects and functions for geospatial data

## Linking to GEOS 3.8.1, GDAL 3.2.1, PROJ 7.2.1

library(rgdal)        # Functions for spatial data input/output

## Loading required package: sp

## rgdal: version: 1.5-23, (SVN revision 1121)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.2.1, released 2020/12/29
## Path to GDAL shared files: /Library/Frameworks/R.framework/Versions/4.1/Resources/library/rgdal/gdal
## GDAL binary built with GEOS: TRUE 
## Loaded PROJ runtime: Rel. 7.2.1, January 1st, 2021, [PJ_VERSION: 721]
## Path to PROJ shared files: /Library/Frameworks/R.framework/Versions/4.1/Resources/library/rgdal/proj
## PROJ CDN enabled: FALSE
## Linking to sp version:1.4-5
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading rgdal.
## Overwritten PROJ_LIB was /Library/Frameworks/R.framework/Versions/4.1/Resources/library/rgdal/proj

library(ggplot2)      # Graphing functions
library(dplyr)        # Functions for processing tabular data

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(tidyr)        # Functions for processing tabular data
library(scales)       # Additional graphics functions
library(RColorBrewer) # Color ramps for graphs and maps
library(gridExtra)    # Functions for arranging multiple plots on a page

## 
## Attaching package: 'gridExtra'

## The following object is masked from 'package:dplyr':
## 
##     combine

library(readr)

## 
## Attaching package: 'readr'

## The following object is masked from 'package:scales':
## 
##     col_factor

Get the data.

sdcounty <- st_read(dsn = "SD_counties.shp")

## Reading layer `SD_counties' from data source `/Users/haroldnelson/Dropbox/Wimberly/gdswr_chapter4_data/SD_counties.shp' using driver `ESRI Shapefile'
## Simple feature collection with 66 features and 18 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -167412.2 ymin: -264522.5 xmax: 456081.5 ymax: 118174.9
## Projected CRS: Custom

sdclimate <- read_csv(file = "sd_climate.csv")

## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   .default = col_double()
## )
## ℹ Use `spec()` for the full column specifications.

sdclimate <- mutate(sdclimate, fips = as.character(fips))
sdcounty <- left_join(sdcounty, sdclimate, by = c("FIPS" = "fips"))
sdcntrd = st_centroid(sdcounty)

## Warning in st_centroid.sf(sdcounty): st_centroid assumes attributes are constant
## over geometries of x

Practice Problem 1

Experiment with different versions of the West Nile virus incidence map. Try different color ramps from the RColorBrewer package and different numbers and breakpoints for classes in the choropleth maps.

First let’s examine the basic map code.

sdcounty <- mutate(sdcounty, 
                   inc2003 = 100000 * C2003/P2003, 
                   inc2003_c1 = cut(inc2003, breaks = quantile(inc2003,probs = seq(0, 1, 0.2)),
                                    include.lowest = T))
sdcntrd <- mutate(sdcntrd,
                  pop2003_c1 = cut(P2003, breaks = c(0, 10000, 50000, 200000),
                                   include.lowest = T, dig.lab=6))

ggplot() +
  geom_sf(data = sdcounty, aes(fill = inc2003_c1), color = "black", size = 0.25) + 
  scale_fill_brewer(name="Incidence", palette = "BuPu",
                    guide = guide_legend(override.aes = list(linetype = "blank", shape = NA))) +
  geom_sf(data = sdcntrd, aes(size = pop2003_c1), color = "black", show.legend = "point") +
  scale_size_discrete(name="Population") + 
  theme_bw() +
  labs(title="2003 West Nile Virus Incidence in South Dakota")

## Warning: Using size for a discrete variable is not advised.

First, let’s talk about how we could change the cuts in incidence.

https://www.youtube.com/watch?v=QmP2iiTw9fQ

Now, let’s talk about changing the palette.

https://www.youtube.com/watch?v=hThc_8X1KwI

Practice Problem 2

Generate mutli-panel map of mean monthly temperatures. The monthly temperature variables are t3001, t3002… t3012.

Basically we need to take the code for precipitation and repurpose it for temperature.

Here’s the precip code.

Here’s a video of me changing it to temperature.

https://www.youtube.com/watch?v=47TiKCVU1No

monthly_temp <- sdcounty %>%
  pivot_longer(sdcounty,
               cols = t3001:t3012,
               values_to = "temp",
               names_to = "monthvar") %>%
  select(FIPS, monthvar, temp, geometry) %>%
  mutate(monthvar = factor(monthvar, labels = month.name))

## Warning in gsub(paste0("^", names_prefix), "", names(cols)): argument 'pattern'
## has length > 1 and only the first element will be used

## Warning in val_cols[col_id] <- unname(as.list(data[cols])): number of items to
## replace is not a multiple of replacement length

glimpse(monthly_temp)

## Rows: 792
## Columns: 4
## $ FIPS     <chr> "46105", "46105", "46105", "46105", "46105", "46105", "46105"…
## $ monthvar <fct> January, February, March, April, May, June, July, August, Sep…
## $ temp     <dbl> -279.557, 3.004, 537.272, 1376.880, 2032.190, 2547.500, 3037.…
## $ geometry <POLYGON [m]> POLYGON ((40.80017 105012.2..., POLYGON ((40.80017 10…

Graph

We’ll just take the graph for precipation and use it for temperature.

You can watch at https://www.youtube.com/watch?v=XUvNMsn8jbg

monthly_temp <- st_as_sf(monthly_temp)

ggplot() +
  geom_sf(data = monthly_temp, aes(fill = temp), color = "black", size = 0.25) + 
  scale_fill_distiller(name="Temp", palette = "Blues", breaks = pretty_breaks())+
  facet_wrap(~ monthvar, ncol = 4) +
  labs(title="Monthly Climatology of Temperature in South Dakota")