Install mapping packages

#install.packages("tmap")
#install.packages("scales")
#install.packages("leaflet.extras")
#install.packages("rio")
#install.packages("htmlwidgets")
#install.packages("sf")
#install.packages("sp")
#install.packages("leaflet")
#install.packages("raster")
#install.packages("rgdal")
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ ggplot2 3.4.0      ✔ purrr   0.3.5 
## ✔ tibble  3.1.8      ✔ dplyr   1.0.10
## ✔ tidyr   1.2.1      ✔ stringr 1.4.1 
## ✔ readr   2.1.3      ✔ forcats 0.5.2 
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
library(raster)
## Loading required package: sp
## 
## Attaching package: 'raster'
## 
## The following object is masked from 'package:dplyr':
## 
##     select
library(rgdal) 
## Please note that rgdal will be retired by the end of 2023,
## plan transition to sf/stars/terra functions using GDAL and PROJ
## at your earliest convenience.
## 
## rgdal: version: 1.5-32, (SVN revision 1176)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.4.2, released 2022/03/08
## Path to GDAL shared files: /Library/Frameworks/R.framework/Versions/4.2/Resources/library/rgdal/gdal
## GDAL binary built with GEOS: FALSE 
## Loaded PROJ runtime: Rel. 8.2.1, January 1st, 2022, [PJ_VERSION: 821]
## Path to PROJ shared files: /Library/Frameworks/R.framework/Versions/4.2/Resources/library/rgdal/proj
## PROJ CDN enabled: FALSE
## Linking to sp version:1.5-0
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading sp or rgdal.
library(tmap)
library(tmaptools)
library(sf)
## Linking to GEOS 3.10.2, GDAL 3.4.2, PROJ 8.2.1; sf_use_s2() is TRUE
library(leaflet)
library(leaflet.extras)
library(dplyr)
library(rio)
library(sp)
library(scales)
## 
## Attaching package: 'scales'
## 
## The following object is masked from 'package:purrr':
## 
##     discard
## 
## The following object is masked from 'package:readr':
## 
##     col_factor
library(htmlwidgets)

Set various values needed, including names of files and FIPS codes for New Hampshire and South Carolina

nhdatafile <- "NHD2016.xlsx"
nhdatafilecsv <- "NHD2016.csv"
usshapefile <- "cb_2014_us_county_5m/cb_2014_us_county_5m.shp"
nhfipscode <- "33"
scdatafile <- "SCGOP2016.csv"
scfipscode <- "45"

Step 1: Get election results data

Start with the 2016 New Hampshire Democratic primary results

setwd("/Users/KathyOchoa/Documents/DATA 110/GIS")
nhdata <- import(nhdatafilecsv)

Select just the County, Clinton and Sanders columns

nhdata <- nhdata[,c("County", "Clinton", "Sanders")]

Step 2: Decide what data to map

Add columns for candidates’ margins of victory (or loss) and percent of the vote

nhdata$SandersMarginVotes <- nhdata$Sanders - nhdata$Clinton
nhdata$SandersPct <- (nhdata$Sanders) / (nhdata$Sanders + nhdata$Clinton)
# Will use formatting later to multiply by a hundred
nhdata$ClintonPct <- (nhdata$Clinton) / (nhdata$Sanders + nhdata$Clinton)
nhdata$SandersMarginPctgPoints <- nhdata$SandersPct - nhdata$ClintonPct

Step 3: Get geographic data files

setwd("/Users/KathyOchoa/Documents/DATA 110/GIS")
usgeo <- shapefile("cb_2014_us_county_5m/cb_2014_us_county_5m.shp")
## Warning: [vect] Z coordinates ignored
## Warning in .getSpatDF(x@ptr$df, ...): NAs introduced by coercion to integer
## range

## Warning in .getSpatDF(x@ptr$df, ...): NAs introduced by coercion to integer
## range

Do a quick qtm

qtm(usgeo)

#view(usgeo)

Subset just the NH data from the US shapefile

nhgeo <- usgeo[usgeo$STATEFP==nhfipscode,]

qtm of the New Hampshire data

qtm(nhgeo)

Step 3.5: Get the data into the right format

Check that the county names match up

str(nhgeo$NAME)
##  chr [1:10] "Grafton" "Hillsborough" "Coos" "Belknap" "Rockingham" ...
str(nhdata$County)
##  chr [1:10] "Belknap" "Carroll" "Cheshire" "Coos" "Grafton" "Hillsborough" ...
# Change the county names to plain characters in nhgeo:
nhgeo$NAME <- as.character(nhgeo$NAME) 
nhgeo$NAME
##  [1] "Grafton"      "Hillsborough" "Coos"         "Belknap"      "Rockingham"  
##  [6] "Cheshire"     "Strafford"    "Merrimack"    "Carroll"      "Sullivan"

Order each data set by county name

nhgeo <- nhgeo[order(nhgeo$NAME),]
nhdata <- nhdata[order(nhdata$County),]
# Are the two county columns identical now? They should be:
identical(nhgeo$NAME,nhdata$County)
## [1] TRUE

Step 4: Merge spatial and results data

nhmap <- merge(nhgeo, nhdata, by.x = "NAME", by.y = "County")
# See the new data structure with
#str(nhmap)

Step 5: Create a static map

qtm(nhmap, "SandersMarginVotes") 
## Some legend labels were too wide. These labels have been resized to 0.63, 0.63, 0.63, 0.58, 0.54. Increase legend.width (argument of tm_layout) to make the legend wider and therefore the labels larger.

qtm(nhmap, "SandersMarginPctgPoints")

For more control over the map’s colors, borders and such, use the tm_shape() function

tm_shape(nhmap) +
  tm_fill("SandersMarginVotes", title="Sanders Margin, Total Votes", palette = "PRGn") +
  tm_borders(alpha=.5) +
  tm_text("NAME", size=0.8)
## Some legend labels were too wide. These labels have been resized to 0.63, 0.63, 0.63, 0.58, 0.54. Increase legend.width (argument of tm_layout) to make the legend wider and therefore the labels larger.

Same code as above, but store the static map in a variable, and change the theme to “classic” style:

nhstaticmap <- tm_shape(nhmap) +
  tm_fill("SandersMarginVotes", title="Sanders Margin, Total Votes", palette = "viridis") +
  tm_borders(alpha=.5) +
  tm_text("NAME", size=0.8) +
  tm_style("classic") 
nhstaticmap
## Some legend labels were too wide. These labels have been resized to 0.63, 0.63, 0.63, 0.58, 0.54. Increase legend.width (argument of tm_layout) to make the legend wider and therefore the labels larger.

Save the map to a jpg file

tmap_save(nhstaticmap, filename="nhdemprimary.jpg") 
## Map saved to /Users/KathyOchoa/Documents/DATA 110/nhdemprimary.jpg
## Resolution: 1501.336 by 2937.385 pixels
## Size: 5.004452 by 9.791282 inches (300 dpi)

Part 6: Create palette and pop-ups for interactive map

Create a palette

# "Blues" is a range of blues from ColorBrewer and domain is the data range of the color scale
clintonPalette <- colorNumeric(palette = "Blues", domain=nhmap$ClintonPct)

And a pop-up window

nhpopup <- paste0("County: ", nhmap$NAME,
"Sanders ", percent(nhmap$SandersPct), " - Clinton ", percent(nhmap$ClintonPct))

Step 7: Generate an interactive map

nhmap_projected <- sp::spTransform(nhmap, "+proj=longlat +datum=WGS84")
leaflet(nhmap_projected) %>%
  addProviderTiles("CartoDB.Positron") %>%
  addPolygons(stroke=FALSE, smoothFactor = 0.2, fillOpacity = .8, popup=nhpopup, color=~clintonPalette(nhmap$ClintonPct))

Step 8: Add palettes for a multi-layer map

South Carolina data

setwd("/Users/KathyOchoa/Documents/DATA 110/GIS")
scdata <- rio::import(scdatafile)

SC shapefile and quickplots

scgeo <- usgeo[usgeo@data$STATEFP=="45",]
qtm(scgeo)

Add a column with percent of votes for each candidate

# Candidates are in columns 2-7:
candidates <- colnames(scdata[2:7])
for(i in 2:7){
  j = i + 7
  temp <- scdata[[i]] / scdata$Total
  scdata[[j]] <- temp
  colnames(scdata)[j] <- paste0(colnames(scdata)[i], "Pct")
}
winner <- colnames(scdata[2:7])

Get winner in each precinct

for(i in 1:nrow(scdata)) {
  scdata$winner[i] <- names(which.max(scdata[i,2:7])) 
}

Import spreadsheet with percent of adult population holding at least a 4-yr college degree

setwd("/Users/KathyOchoa/Documents/DATA 110/GIS")
sced <- rio::import("SCdegree.xlsx")

Check if county names are in the same format in both files

str(scgeo$NAME)
##  chr [1:46] "Edgefield" "Lee" "Horry" "Allendale" "Marion" "Dorchester" ...
str(scdata$County)
##  chr [1:46] "Abbeville" "Aiken" "Allendale" "Anderson" "Bamberg" "Barnwell" ...
# Change the county names to plain characters in scgeo:
scgeo$NAME <- as.character(scgeo$NAME)
# Order each data set by county name
scgeo <- scgeo[order(scgeo$NAME),]
scdata <- scdata[order(scdata$County),]
# Are the two county columns identical now? They should be:
identical(scgeo$NAME,scdata$County ) 
## [1] TRUE

Add the election results and rename county column

scmap <- merge(scgeo, scdata, by.x = "NAME", by.y = "County")

Color by strength of win with multiple layers

# Use same intensity for all - get minimum and maximum for the top 3 combined
minpct <- min(c(scdata$`Donald J TrumpPct`, scdata$`Marco RubioPct`, scdata$`Ted CruzPct`
))
maxpct <- max(c(scdata$`Donald J TrumpPct`, scdata$`Marco RubioPct`, scdata$`Ted CruzPct`
))

Create leaflet palettes for each layer of the map:

trumpPalette <- colorNumeric(palette = "Purples", domain=c(minpct, maxpct))
rubioPalette <- colorNumeric(palette = "Reds", domain = c(minpct, maxpct))
cruzPalette <- colorNumeric(palette = "Oranges", domain = c(minpct, maxpct))
winnerPalette <- colorFactor(palette=c("#984ea3", "#e41a1c"), domain = scmap$winner)
edPalette <- colorNumeric(palette = "Blues", domain=scmap$PctCollegeDegree)

Create a pop-up:

scpopup <- paste0("<b>County: ", scmap$NAME, "<br />Winner: ", scmap$winner, "</b><br />< br />Trump: ", percent(scmap$`Donald J TrumpPct`), "<br />Rubio: ", percent(scmap$`Marco RubioPct`), "<br />Cruz: ", percent(scmap$`Ted CruzPct`), "<br /><br />Pct w college ed: ", sced$PctCollegeDegree, "% vs state-wide avg of 25%")

Add the projection we know

scmap <- sp::spTransform(scmap, "+proj=longlat +datum=WGS84")

Basic interactive map showing winner in each county:

leaflet(scmap) %>%
  addProviderTiles("CartoDB.Positron") %>%
  addPolygons(stroke=TRUE, weight=1, smoothFactor = 0.2, fillOpacity = .75, popup=scpopup, color= ~winnerPalette(scmap$winner), group="Winners" ) %>%
  addLegend(position="bottomleft", colors=c("#984ea3", "#e41a1c"), labels=c("Trump", "Rubio"))

Step 9: Add map layers and controls

Put top 3 candidates in their own layers and add education layer, store in scGOPmap variable

scGOPmap <- leaflet(scmap) %>%
  addProviderTiles("CartoDB.Positron") %>%
  addPolygons(stroke=TRUE, weight=1, smoothFactor = 0.2, fillOpacity = .75, popup=scpopup, color= ~winnerPalette(scmap$winner), group="Winners" ) %>%
  addLegend(position="bottomleft", colors=c("#984ea3", "#e41a1c"), labels=c("Trump", "Rubio")) %>%
  
  addPolygons(stroke=TRUE, weight=1, smoothFactor = 0.2, fillOpacity = .75, popup=scpopup, color= ~trumpPalette(scmap$`Donald J TrumpPct`), group="Trump") %>%
  addPolygons(stroke=TRUE, weight=1, smoothFactor = 0.2, fillOpacity = .75, popup=scpopup, color= ~rubioPalette(scmap$`Marco RubioPct`), group="Rubio") %>% 
  addPolygons(stroke=TRUE, weight=1, smoothFactor = 0.2, fillOpacity = .75, popup=scpopup, color= ~cruzPalette(scmap$`Ted CruzPct`), group="Cruz") %>%
  
  addPolygons(stroke=TRUE, weight=1, smoothFactor = 0.2, fillOpacity = .75, popup=scpopup, color= ~edPalette(sced$PctCollegeDegree), group="College degs") %>%
  
  addLayersControl(baseGroups=c("Winners", "Trump", "Rubio", "Cruz", "College degs"),
 position = "bottomleft",
 options = layersControlOptions(collapsed = FALSE))

# Now display the map
scGOPmap

Step 10: Save your interactive map

#htmlwidgets::saveWidget(scGOPmap2, file="scGOPwidget2.html")
# save as an HTML file with dependencies in another directory:
#htmlwidgets::saveWidget(widget=scGOPmap2, file="scGOPprimary_withdependencies.html", selfcontained=FALSE, libdir = "js")