R Spatial Lab Assignment 9

Introduction

This assignment explores spatial data visualization using ZIP-code polygons and COVID-related demographic data. I created two static choropleth maps with ggplot2, a side by side comparison map, and one interactive leaflet web map with multiple layers.

Load packages

library(sf)
library(dplyr)
library(ggplot2)
library(classInt)
library(RColorBrewer)
library(leaflet)
library(ggpubr)
library(htmlwidgets)

Load and inspect data

acsCovidDat <- st_read("data/acsPopByZip.gpkg", layer = "NYC_COVID_Pop_by_Zip") %>%
  st_transform(4326)

Reading layer `NYC_COVID_Pop_by_Zip' from data source 
  `C:\Users\danyu\OneDrive - Hunter - CUNY\week9\data\acsPopByZip.gpkg' 
  using driver `GPKG'
Simple feature collection with 180 features and 13 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 913129 ymin: 120020.9 xmax: 1067113 ymax: 272710.9
Projected CRS: NAD83 / New York Long Island (ftUS)

Static Map 1: Confirmed COVID-19 Cases by ZIP Code

breaks_qt <- classIntervals(
  c(min(acsCovidDat$Positiv, na.rm = TRUE) - 0.00001, acsCovidDat$Positiv),
  n = 7,
  style = "quantile"
)

acsCovidDat <- mutate(
  acsCovidDat,
  Positiv_cat = cut(Positiv, breaks_qt$brks, dig.lab = 4, digits = 1)
)

g1 <- ggplot(acsCovidDat) +
  geom_sf(aes(fill = Positiv_cat), color = "white", linewidth = 0.2) +
  scale_fill_brewer(palette = "OrRd", name = "COVID Cases") +
  labs(
    x = "Longitude",
    y = "Latitude",
    title = "Confirmed COVID-19 Cases by ZIP Code",
    caption = "Data Source: ACS / ZIP-level COVID data"
  ) +
  theme_bw() +
  theme(
    panel.grid.major = element_line(color = "gray80"),
    panel.border = element_rect(color = "black", fill = NA),
    legend.position = "right"
  )

g1

Static Map 2: Elderly Population by ZIP Code

breaks_elderly <- classIntervals(
  c(min(acsCovidDat$eldrlyP, na.rm = TRUE) - 0.00001, acsCovidDat$eldrlyP),
  n = 7,
  style = "quantile"
)

acsCovidDat <- mutate(
  acsCovidDat,
  elderly_cat = cut(eldrlyP, breaks_elderly$brks, dig.lab = 4, digits = 1)
)

g2 <- ggplot(acsCovidDat) +
  geom_sf(aes(fill = elderly_cat), color = "white", linewidth = 0.2) +
  scale_fill_brewer(palette = "YlGnBu", name = "Elderly Population") +
  labs(
    x = "Longitude",
    y = "Latitude",
    title = "Elderly Population by ZIP Code",
    caption = "Data Source: ACS / ZIP-level COVID data"
  ) +
  theme_bw() +
  theme(
    panel.grid.major = element_line(color = "gray80"),
    panel.border = element_rect(color = "black", fill = NA),
    legend.position = "right"
  )

g2

Side-by-side comparison

combined_map <- ggarrange(g1, g2, ncol = 2, nrow = 1)
combined_map

Interactive Leaflet Map

pal_fun <- colorQuantile("YlOrRd", NULL, n = 5)

pal_elderly <- colorQuantile(
  "YlGnBu",
  domain = c(min(acsCovidDat$eldrlyP, na.rm = TRUE),
             max(acsCovidDat$eldrlyP, na.rm = TRUE)),
  alpha = TRUE
)

pal_tests <- colorQuantile(
  "Blues",
  domain = c(min(acsCovidDat$Total, na.rm = TRUE),
             max(acsCovidDat$Total, na.rm = TRUE)),
  alpha = TRUE
)

covid_popup <- paste0(
  "<strong>Zip Code: </strong>", acsCovidDat$ZIPCODE,
  "<br><strong>Place Name: </strong>", acsCovidDat$PO_NAME,
  "<br><strong>Confirmed COVID-19 Cases: </strong>", acsCovidDat$Positiv
)

polyHighlightOption <- leaflet::highlightOptions(opacity = 1.0, fillColor = "black")
polyLabelOption <- leaflet::labelOptions(opacity = 1.0, textsize = "11px")

htmlMap <- leaflet(acsCovidDat) %>%
  addPolygons(
    stroke = TRUE,
    fillOpacity = 0.8,
    fillColor = ~pal_elderly(eldrlyP),
    color = "grey",
    weight = 1,
    highlightOptions = polyHighlightOption,
    label = paste("Elderly population:", acsCovidDat$eldrlyP),
    popup = paste0(
      "<strong>ZIP:</strong> ", acsCovidDat$ZIPCODE,
      "<br><strong>Elderly population:</strong> ", acsCovidDat$eldrlyP
    ),
    group = "Elderly Population"
  ) %>%
  addPolygons(
    stroke = FALSE,
    fillColor = ~pal_fun(Positiv),
    fillOpacity = 0.8,
    smoothFactor = 0.5,
    popup = covid_popup,
    group = "NYC-COVID",
    label = paste0(
      acsCovidDat$Positiv, " at ", acsCovidDat$PO_NAME, ", ", acsCovidDat$ZIPCODE
    ),
    highlightOptions = polyHighlightOption,
    labelOptions = polyLabelOption
  ) %>%
  addPolygons(
    stroke = FALSE,
    fillColor = ~pal_tests(Total),
    fillOpacity = 0.8,
    smoothFactor = 0.5,
    popup = paste0(
      "<strong>ZIP:</strong> ", acsCovidDat$ZIPCODE,
      "<br><strong>Total tests:</strong> ", acsCovidDat$Total
    ),
    group = "Total Tests",
    label = paste("Total tests:", acsCovidDat$Total),
    highlightOptions = polyHighlightOption,
    labelOptions = polyLabelOption
  ) %>%
  addTiles(group = "OSM") %>%
  addProviderTiles("CartoDB.DarkMatter", group = "Carto") %>%
  addLegend(
    "bottomright",
    pal = pal_fun,
    values = acsCovidDat$Positiv,
    title = "Confirmed COVID-19 Cases",
    group = "NYC-COVID"
  ) %>%
  addLayersControl(
    baseGroups = c("OSM", "Carto"),
    overlayGroups = c("NYC-COVID", "Elderly Population", "Total Tests")
  )

htmlMap

Conclusion

This assignment produced two static choropleth maps, a side byside comparison figure, and one interactive Leaflet map with multiple layers. These visualizations help compare COVID cases, total tests, and elderly population patterns across ZIP codes.