GIS Tutorial Practice

Preparation

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"

#```{r}

install.packages(“tmap”) install.packages(“tmaptools”) install.packages(“leaflet”) install.packages(“scales”) install.packages(“leaflet.extras”) install.packages(“rio”) install.packages(“htmlwidgets”) install.packages(“sf”) install.packages(“sp”) #```

library(tmap)
library(leaflet)
library(scales)
library(leaflet.extras)
library(rio)
library(htmlwidgets)
library(sf)
library(tidyverse)
library(sp)

Step 1: Read in the NH election results file:

setwd("C:/Users/ykim2/Downloads/MC/DATA110/GIS")

nhdata <- import(nhdatafilecsv)

Eliminate columns for minor candidates and just use County, Clinton and Sanders columns:

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

Step 2: Decide what data to map

Add columns for percents and margins

nhdata$SandersMarginVotes <- nhdata$Sanders - nhdata$Clinton
nhdata$SandersPct <- (nhdata$Sanders) / (nhdata$Sanders + nhdata$Clinton) 
nhdata$ClintonPct <- (nhdata$Clinton) / (nhdata$Sanders + nhdata$Clinton)
nhdata$SandersMarginPctgPoints <- nhdata$SandersPct - nhdata$ClintonPct
#view(nhdata)

Step 3: Get Geographic data files

Read in the shapefile for US states and countries:

#install.packages("raster")
#install.packages("rgdal")
library(sp)
library(raster)
## 
## 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-29, (SVN revision 1165M)
## 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: C:/Users/ykim2/OneDrive/Documents/R/win-library/4.1/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: C:/Users/ykim2/OneDrive/Documents/R/win-library/4.1/rgdal/proj
## PROJ CDN enabled: FALSE
## Linking to sp version:1.4-6
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading sp or rgdal.
## Overwritten PROJ_LIB was C:/Users/ykim2/OneDrive/Documents/R/win-library/4.1/rgdal/proj
usgeo <- shapefile("cb_2014_us_county_5m/cb_2014_us_county_5m.shp")
## Warning in rgdal::readOGR(dirname(x), fn, stringsAsFactors = stringsAsFactors, :
## Z-dimension discarded

Do a quick plot(qtm stands for quick thematic map) of the shapefile and check its structure:

qtm(usgeo)

#view(usgeo)

Subset just the NH data from the US shapefile

nhgeo <- usgeo[usgeo$STATEFP==nhfipscode,]
#view(nhgeo)

tmap test plot of the New Hampshire data

qtm(nhgeo)

### Structure of the object

str(nhgeo)
## Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots
##   ..@ data       :'data.frame':  10 obs. of  9 variables:
##   .. ..$ STATEFP : chr [1:10] "33" "33" "33" "33" ...
##   .. ..$ COUNTYFP: chr [1:10] "009" "011" "007" "001" ...
##   .. ..$ COUNTYNS: chr [1:10] "00873178" "00873179" "00873177" "00873174" ...
##   .. ..$ AFFGEOID: chr [1:10] "0500000US33009" "0500000US33011" "0500000US33007" "0500000US33001" ...
##   .. ..$ GEOID   : chr [1:10] "33009" "33011" "33007" "33001" ...
##   .. ..$ NAME    : chr [1:10] "Grafton" "Hillsborough" "Coos" "Belknap" ...
##   .. ..$ LSAD    : chr [1:10] "06" "06" "06" "06" ...
##   .. ..$ ALAND   : chr [1:10] "4425927252" "2269220216" "4648216798" "1036582289" ...
##   .. ..$ AWATER  : chr [1:10] "105375486" "41604851" "90773891" "177039345" ...
##   ..@ polygons   :List of 10
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71.8 43.9
##   .. .. .. .. .. .. ..@ area   : num 0.508
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:327, 1:2] -72.3 -72.3 -72.3 -72.3 -72.3 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71.8 43.9
##   .. .. .. ..@ ID       : chr "685"
##   .. .. .. ..@ area     : num 0.508
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71.7 42.9
##   .. .. .. .. .. .. ..@ area   : num 0.255
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:45, 1:2] -72.1 -72 -72 -72 -72 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71.7 42.9
##   .. .. .. ..@ ID       : chr "866"
##   .. .. .. ..@ area     : num 0.255
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71.3 44.7
##   .. .. .. .. .. .. ..@ area   : num 0.539
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:511, 1:2] -71.8 -71.8 -71.8 -71.7 -71.7 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71.3 44.7
##   .. .. .. ..@ ID       : chr "922"
##   .. .. .. ..@ area     : num 0.539
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71.4 43.5
##   .. .. .. .. .. .. ..@ area   : num 0.136
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:33, 1:2] -71.7 -71.7 -71.7 -71.7 -71.7 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71.4 43.5
##   .. .. .. ..@ ID       : chr "1100"
##   .. .. .. ..@ area     : num 0.136
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71.1 43
##   .. .. .. .. .. .. ..@ area   : num 0.208
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:123, 1:2] -71.5 -71.4 -71.4 -71.4 -71.4 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71.1 43
##   .. .. .. ..@ ID       : chr "1278"
##   .. .. .. ..@ area     : num 0.208
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -72.3 42.9
##   .. .. .. .. .. .. ..@ area   : num 0.208
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:161, 1:2] -72.6 -72.6 -72.6 -72.6 -72.6 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -72.3 42.9
##   .. .. .. ..@ ID       : chr "1877"
##   .. .. .. ..@ area     : num 0.208
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71 43.3
##   .. .. .. .. .. .. ..@ area   : num 0.11
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:177, 1:2] -71.2 -71.2 -71.2 -71.2 -71.1 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71 43.3
##   .. .. .. ..@ ID       : chr "2676"
##   .. .. .. ..@ area     : num 0.11
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71.7 43.3
##   .. .. .. .. .. .. ..@ area   : num 0.274
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:49, 1:2] -72.1 -72.1 -72 -72.1 -72.1 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71.7 43.3
##   .. .. .. ..@ ID       : chr "2773"
##   .. .. .. ..@ area     : num 0.274
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -71.2 43.9
##   .. .. .. .. .. .. ..@ area   : num 0.288
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:73, 1:2] -71.6 -71.4 -71.4 -71.3 -71.4 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -71.2 43.9
##   .. .. .. ..@ ID       : chr "3077"
##   .. .. .. ..@ area     : num 0.288
##   .. .. .. ..$ comment: chr "0"
##   .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
##   .. .. .. ..@ Polygons :List of 1
##   .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
##   .. .. .. .. .. .. ..@ labpt  : num [1:2] -72.2 43.4
##   .. .. .. .. .. .. ..@ area   : num 0.159
##   .. .. .. .. .. .. ..@ hole   : logi FALSE
##   .. .. .. .. .. .. ..@ ringDir: int 1
##   .. .. .. .. .. .. ..@ coords : num [1:118, 1:2] -72.5 -72.4 -72.4 -72.4 -72.4 ...
##   .. .. .. ..@ plotOrder: int 1
##   .. .. .. ..@ labpt    : num [1:2] -72.2 43.4
##   .. .. .. ..@ ID       : chr "3167"
##   .. .. .. ..@ area     : num 0.159
##   .. .. .. ..$ comment: chr "0"
##   ..@ plotOrder  : int [1:10] 3 1 9 8 2 6 5 10 4 7
##   ..@ bbox       : num [1:2, 1:2] -72.6 42.7 -70.7 45.3
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "x" "y"
##   .. .. ..$ : chr [1:2] "min" "max"
##   ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
##   .. .. ..@ projargs: chr "+proj=longlat +datum=NAD83 +no_defs"
##   .. .. ..$ comment: chr "GEOGCRS[\"NAD83\",\n    DATUM[\"North American Datum 1983\",\n        ELLIPSOID[\"GRS 1980\",6378137,298.257222"| __truncated__
##   ..$ comment: chr "TRUE"
# Check if county names are in the same format in both files
str(nhgeo$NAME)
##  chr [1:10] "Grafton" "Hillsborough" "Coos" "Belknap" "Rockingham" ...
str(nhdata$County)
##  chr [1:10] "Belknap" "Carroll" "Cheshire" "Coos" "Grafton" "Hillsborough" ...

Order each data set by county name

nhgeo <-nhgeo[order(nhgeo$NAME),]
nhdata <- nhdata[order(nhdata$County),]
if (identical(nhgeo$NAME, nhdata$County)) {nhmap <- merge(nhgeo, nhdata, by.x = "NAME", by.y = "County")} else {stop} # merge geo and vote datasets.

Step 4: Merge geo data with results data using the merge function

nhmap <- merge(nhgeo, nhdata, by.x = "NAME", by.y = "County")

Step 5: Create a static map with tmap’s qtm() funtion

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

qtm(nhmap, "SandersMarginPctgPoints")

For more control over look and feel, 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.62, 0.62, 0.62, 0.57, 0.53. 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.62, 0.62, 0.62, 0.57, 0.53. Increase legend.width (argument of tm_layout) to make the legend wider and therefore the labels larger.

# Save the map to a jpg file with tmap’s tmap_save():

tmap_save(nhstaticmap, filename="nhdemprimary.jpg")
## Map saved to C:\Users\ykim2\Downloads\MC\DATA110\GIS\nhdemprimary.jpg
## Resolution: 1501.336 by 2937.385 pixels
## Size: 5.004452 by 9.791282 inches (300 dpi)

Step 6

Next up : Code for a basic interactive map, this time for Clinton percentages in NH

Create a palette

clintonPalette <- colorNumeric(palette="Blues", domain=nhmap$ClintonPct)

and a pop-up window

library(scales)
nhpopup <- paste0("County: ", nhmap$NAME, "<br /><br />Sanders, ",percent(nhmap$SandersPct), " - clinton ", percent(nhmap$ClintonPct))

Step 7: Now generate the interactive map:

# re-project
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

scdata <- rio::import(scdatafile)
scdata
##          County Jeb Bush Ben Carson Ted Cruz John R Kasich Marco Rubio
## 1     Abbeville      236        305      876           157         740
## 2         Aiken     1238       2642     7138          1957        6886
## 3     Allendale       21         32       44            12          92
## 4      Anderson     2411       2959     8750          2038        7175
## 5       Bamberg       58         78      193            57         333
## 6      Barnwell      129        171      516            91         575
## 7      Beaufort     2185       1095     4407          6245        9126
## 8      Berkeley     2662       1641     6244          1555        6176
## 9       Calhoun      250        139      564           110         421
## 10   Charleston     6140       2467     9050          7301       15247
## 11     Cherokee      382        879     1978           308        1449
## 12      Chester      126        284      844           119         553
## 13 Chesterfield      180        246     1434           139         695
## 14    Clarendon      354        267     1036           243         793
## 15     Colleton      503        281     1103           223         936
## 16   Darlington      693        525     2594           328        1387
## 17       Dillon      169        112      912            65         445
## 18   Dorchester     1951       1310     5199          1418        5244
## 19    Edgefield      116        326     1148           151         900
## 20    Fairfield      204        145      631           164         496
## 21     Florence     1462       1267     5208           823        3357
## 22   Georgetown      959        609     1591          1023        2693
## 23   Greenville     6548       8036    22899          7683       22875
## 24    Greenwood      866        989     2408           697        2543
## 25      Hampton       83         92      365            88         314
## 26        Horry     3161       2638     8496          3430        9605
## 27       Jasper      119        157      566           275         695
## 28      Kershaw     1150        616     2772           579        1896
## 29    Lancaster      561       1136     2689          1047        2754
## 30      Laurens      765        791     2745           436        1580
## 31          Lee      129         60      328            53         173
## 32    Lexington     5298       4847    13167          3449       11419
## 33       Marion      215        140      740           100         436
## 34     Marlboro      164        128      409            57         254
## 35    McCormick       76        106      371           180         468
## 36     Newberry      619        522     1267           343        1156
## 37       Oconee     1391       1123     3164          1477        4002
## 38   Orangeburg      779        404     1303           353        1305
## 39      Pickens     2033       1714     6393          1448        4961
## 40     Richland     4484       2633     7964          4051       10955
## 41       Saluda      382        295      797           156         553
## 42  Spartanburg     2744       4036    11478          2457       10691
## 43       Sumter     1352        817     2474           532        2236
## 44        Union      269        189      884           114         752
## 45 Williamsburg      230        152      575            76         520
## 46         York     2036       3953     9169          2632        8089
##    Donald J Trump Total
## 1            1353  3667
## 2            9148 29009
## 3             160   361
## 4           10953 34286
## 5             480  1199
## 6            1106  2588
## 7            9912 32970
## 8            9406 27684
## 9             805  2289
## 10          14420 54625
## 11           3616  8612
## 12           1403  3329
## 13           1806  4500
## 14           1601  4294
## 15           2326  5372
## 16           3352  8879
## 17            972  2675
## 18           6949 22071
## 19           1366  4007
## 20            919  2559
## 21           6004 18121
## 22           4620 11495
## 23          24913 92954
## 24           3029 10532
## 25            721  1663
## 26          26445 53775
## 27           1169  2981
## 28           3673 10686
## 29           4190 12377
## 30           3321  9638
## 31            650  1393
## 32          16361 54541
## 33           1206  2837
## 34            737  1749
## 35            615  1816
## 36           2043  5950
## 37           4424 15581
## 38           2536  6680
## 39           7136 23685
## 40           9172 39259
## 41           1248  3431
## 42          15249 46655
## 43           3762 11173
## 44           1403  3611
## 45           1186  2739
## 46          12031 37910

Couth Carolina shapefile and Quick plot of scgeo SC geospatial object:

scgeo <- usgeo[usgeo@data$STATEFP=="45", ]
scgeo
## class       : SpatialPolygonsDataFrame 
## features    : 46 
## extent      : -83.35324, -78.54109, 32.0346, 35.2154  (xmin, xmax, ymin, ymax)
## crs         : +proj=longlat +datum=NAD83 +no_defs 
## variables   : 9
## names       : STATEFP, COUNTYFP, COUNTYNS,       AFFGEOID, GEOID,      NAME, LSAD,      ALAND,    AWATER 
## min values  :      45,      001, 01244251, 0500000US45001, 45001, Abbeville,   06, 1016985343, 108743547 
## max values  :      45,      091, 01252740, 0500000US45091, 45091,      York,   06,  987179097,   9746224
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

sced <- rio::import("SCdegree.xlsx") 
sced
##          County PctCollegeDegree
## 1     Abbeville             12.2
## 2         Aiken             24.0
## 3     Allendale             15.0
## 4      Anderson             19.3
## 5       Bamberg             17.5
## 6      Barnwell             11.2
## 7      Beaufort             37.6
## 8      Berkeley             21.3
## 9       Calhoun             16.9
## 10   Charleston             39.4
## 11     Cherokee             14.3
## 12      Chester             12.4
## 13 Chesterfield             12.6
## 14    Clarendon             13.0
## 15     Colleton             14.7
## 16   Darlington             17.0
## 17       Dillon              8.2
## 18   Dorchester             24.5
## 19    Edgefield             17.2
## 20    Fairfield             15.5
## 21     Florence             21.3
## 22   Georgetown             22.7
## 23   Greenville             31.2
## 24    Greenwood             22.1
## 25      Hampton             11.4
## 26        Horry             22.7
## 27       Jasper             12.5
## 28      Kershaw             19.3
## 29    Lancaster             18.7
## 30      Laurens             14.5
## 31          Lee              9.1
## 32    Lexington             28.7
## 33       Marion             13.6
## 34     Marlboro              8.7
## 35    McCormick             18.2
## 36     Newberry             18.9
## 37       Oconee             21.7
## 38   Orangeburg             18.6
## 39      Pickens             23.0
## 40     Richland             36.1
## 41       Saluda             12.4
## 42  Spartanburg             21.3
## 43       Sumter             19.2
## 44        Union             12.9
## 45 Williamsburg             12.2
## 46         York             28.5

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" ...
# Order each data set by county name

scgeo <- scgeo[order(scgeo$NAME), ]
scdata <- scdata[order(scdata$County), ]

#view(scdata)
# 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") 

Instead of just coloring the winner, let’s color by strength of win with multiple layers

# Use the 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`))
minpct
## [1] 0.1218837

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 from the NH map we’ll need for this data on a Leaflet map

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", "R
ubio"))

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))
scGOPmap 

MCPS High School FARMS Map

I will try to apply those mapping codes to MCPS data which was used for Tableau.

setwd("C:/Users/ykim2/Downloads/MC/DATA110/GIS")
mchs <- read.csv("MoCoHSdata18-19.csv")
#head(mchs)
library(tmap)
library(tmaptools)
library(raster)   # to open shapefile
mcgeo <- shapefile("High School Service Areas/geo_export_acd65c0b-d0d1-403d-9608-62b1d821e428.shp")
## Warning in OGRSpatialRef(dsn, layer, morphFromESRI = morphFromESRI, dumpSRS =
## dumpSRS, : Discarded datum WGS84 in Proj4 definition: +proj=longlat +ellps=WGS84
## +no_defs
qtm(mcgeo)

str(mcgeo$schoolname)
##  chr [1:25] "Seneca Valley HS" "Quince Orchard HS" "Rockville HS" ...
str(mchs$Schoolname)
##  chr [1:25] "Bethesda-Chevy Chase HS" "Blair HS" "Blake HS" "Churchill HS" ...
mcgeo <- mcgeo[order(mcgeo$schoolname),]
df_hs <- mchs[order(mchs$Schoolname),]
str(mcgeo$schoolname)
##  chr [1:25] "Bethesda-Chevy Chase HS" "Blair HS" "Blake HS" "Churchill HS" ...
str(mchs$Schoolname)
##  chr [1:25] "Bethesda-Chevy Chase HS" "Blair HS" "Blake HS" "Churchill HS" ...
mchsmap<- merge(mcgeo, mchs, by.x = "schoolname", by.y = "Schoolname")
#str(mchsmap)
qtm(mchsmap, "schoolname")

Because of too many schools, it is hard to identify each school.

tm_shape(mchsmap) +
  tm_fill("schoolname", title= "MCPS High Shcool", palette = "PRGn") +
  tm_borders(alpha = 0.7) +
  tm_text("schoolname", size = 0.5)

Interactive map using leaflet function will help giving information about each school.

hspopup <- paste0("<b>School Name: ", mchsmap$schoolname, 
                  "<br/><br/>FARMS: ", mchsmap$X..FARMS,"%",
                  "</b><br/>White: ", mchsmap$X..white, "%",
                  "<br/>Black: ", mchsmap$X..black, "%",
                  "<br/>Hispanic: ",mchsmap$X..hispanic.latinx,"%",
                  "<br/>Asian: ", mchsmap$X..asian, "%",
                  "<br/>Graduation Rate : ", mchsmap$X..graduation.rate,"%" )

wh_pop <-paste0("<b>School Name: ", mchsmap$schoolname, "</b><br/>White: ", mchsmap$X..white ,"%")
bl_pop <-paste0("<b>School Name: ", mchsmap$schoolname, "</b><br/>Black: ", mchsmap$X..black,"%")
his_pop <-paste0("<b>School Name: ", mchsmap$schoolname, "</b><br/>Hispanic: ", mchsmap$X..hispanic.latinx,"%")
as_pop <-paste0("<b>School Name: ", mchsmap$schoolname, "</b><br/>Asian: ", mchsmap$X..asian,"%")
library(leaflet)
library("ggsci") # add a new palette
hspal0 <- colorFactor(palette = "RdBu", domain = mchsmap$schoolname)   
# Don't use RColorBrewer palette in colorFactor() because there are too many schoolnames(n) so don't generate the popup boxes.
hspal <- colorFactor(palette = "viridis", domain = mchsmap$schoolname)
hspal2 <- colorFactor(viridis_pal(option = "C")(2), domain = mchsmap$schoolname)
hspal3 <- colorFactor(palette = pal_jco()(10), domain = mchsmap$schoolname)

# Create a continuous palette function
wh_pal <- colorNumeric(palette = "Reds" , domain = mchsmap$X..white)
bl_pal <- colorNumeric(palette = "Purples" , domain = mchsmap$X..black)
his_pal <- colorNumeric(palette = "Greens" , domain = mchsmap$X..hispanic.latinx)
as_pal <- colorNumeric(palette = "Blues" , domain = mchsmap$X..asian)
hspal_farms <- colorNumeric(palette = "Dark2", domain = mchsmap$X..FARMS)
leaflet(mchsmap) %>%
 addProviderTiles("CartoDB.Positron") %>% 
  addPolygons(stroke = TRUE,
              #color ="black" ,    # I can't apply line color.
              weight = 1,    # border line width
              opacity = 0.8,     # border line opacity
              smoothFactor = 0.1,    # border line smooth
              fillOpacity = .8, 
              popup= hspopup, 
              color =  ~hspal3(mchsmap$schoolname),
              highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE))  
leaflet(mchsmap) %>%
  addProviderTiles("CartoDB.Positron") %>%
  
  addPolygons(stroke= TRUE, #color = 'black', opacity = 1,
              weight=1, 
              smoothFactor = 0.2, 
              fillOpacity = 1, 
              popup=wh_pop, 
              color= ~wh_pal(mchsmap$X..white),
              group="X..white",
              highlightOptions = highlightOptions(color = "gray", weight = 2, bringToFront = TRUE)) %>%
  
  addPolygons(stroke=TRUE, 
              weight=1, 
              smoothFactor = 0.2, 
              fillOpacity = 1, 
              popup=bl_pop, 
              color= ~bl_pal(mchsmap$X..black),
              group="X..black",
              highlightOptions = highlightOptions(color = "gray", weight = 2, bringToFront = TRUE))  %>%
  
  addPolygons(stroke=TRUE, 
              weight=1, 
              smoothFactor = 0.2, 
              fillOpacity = 1, 
              popup=his_pop, 
              color= ~his_pal(mchsmap$X..hispanic.latinx),
              group="X..hispanic.latinx",
              highlightOptions = highlightOptions(color = "gray", weight = 2, bringToFront = TRUE))  %>%
  
  addPolygons(stroke=TRUE, 
              weight=1, 
              smoothFactor = 0.2, 
              fillOpacity = 1, 
              popup=as_pop, 
              color= ~as_pal(mchsmap$X..asian),
              group="X..asian",
              highlightOptions = highlightOptions(color = "gray", weight = 2, bringToFront = TRUE))  %>%
  
  addLayersControl(baseGroups=c("X..white", "X..black", "X..hispanic.latinx", "X..asian"),
                   position = "bottomleft", 
                   options = layersControlOptions(collapsed = FALSE)) %>%
   htmlwidgets::onRender("
        function() {
            $('.leaflet-control-layers-list').prepend('Race'); }")
The End. Thank you!