Lab 3
Eve Farrat
2026-03-28
Load Packages
library(sf)## Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUElibrary(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.0 ✔ readr 2.2.0
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ ggplot2 4.0.2 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2
## ✔ purrr 1.2.1## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(ggplot2)
library(ggspatial)
library(patchwork)
library(leaflet)
library(htmlwidgets)
library(tmap)Load Final ZIP-Level Dataset from Lab 2
The final ZIP-level spatial dataset created in Lab 2 is used as the base dataset for visualization.
load("section_08/lab2_final_zip_dataset.RData")Create Derived Variables for Mapping
The following variables are created for choropleth mapping and comparison.
zip_map <- zip_final |>
mutate(
covid_rate = (positive / total) * 1000,
elderly_pct = (pop_65_plus / total_pop) * 100,
food_store_rate = (n_food_stores / total_pop) * 10000,
health_facility_rate = (n_health_facilities / total_pop) * 10000
)Static Map 1: COVID-19 Confirmed Case Rate
This choropleth map shows the spatial distribution of confirmed COVID-19 case rate by ZIP code in NYC.
covid_map <- ggplot(zip_map) +
geom_sf(aes(fill = covid_rate), color = "white", linewidth = 0.2) +
scale_fill_viridis_c(
name = "COVID rate\nper 1,000 tested",
option = "plasma",
na.value = "grey90"
) +
labs(
title = "COVID-19 Confirmed Case Rate by ZIP Code in NYC",
x = "Longitude",
y = "Latitude"
) +
annotation_scale(location = "bl", width_hint = 0.3) +
annotation_north_arrow(
location = "tr",
which_north = "true",
style = north_arrow_fancy_orienteering
) +
coord_sf(expand = F) +
theme_minimal() +
theme(
panel.grid.major = element_line(color = "grey80", linewidth = 0.3),
panel.border = element_rect(color = "black", fill = NA, linewidth = 0.8),
plot.title = element_text(face = "bold", size = 14),
legend.title = element_text(face = "bold")
)
covid_map
Static Map 2: Elderly Population Percentage
This choropleth map shows the percentage of elderly population by ZIP code in NYC.
elderly_map <- ggplot(zip_map) +
geom_sf(aes(fill = elderly_pct), color = "white", linewidth = 0.2) +
scale_fill_viridis_c(
name = "Population age 65+\n(percent)",
option = "magma",
na.value = "grey90"
) +
labs(
title = "Percent of Elderly Population by ZIP Code in NYC",
x = "Longitude",
y = "Latitude"
) +
annotation_scale(location = "bl", width_hint = 0.3) +
annotation_north_arrow(
location = "tr",
which_north = "true",
style = north_arrow_fancy_orienteering
) +
coord_sf(expand = F) +
theme_minimal() +
theme(
panel.grid.major = element_line(color = "grey80", linewidth = 0.3),
panel.border = element_rect(color = "black", fill = NA, linewidth = 0.8),
plot.title = element_text(face = "bold", size = 14),
legend.title = element_text(face = "bold")
)
elderly_map
Multi-Map Figure
The following side-by-side figure compares COVID-19 confirmed case rate and elderly population percentage using the same geographic scale.
multi_map <- covid_map + elderly_map + plot_layout(ncol = 2)
multi_map
Save the Multi-Map Figure
ggsave(
"section_08/lab3_multi_map_figure.png",
plot = multi_map,
width = 12,
height = 6,
dpi = 300
)Prepare Data for Interactive Mapping
ZIP centroids are created for overlay point layers representing food stores and health facilities.
zip_centroids <- st_centroid(zip_map)## Warning: st_centroid assumes attributes are constant over geometrieszip_map <- zip_map |>
mutate(
popup_covid = paste0(
"<b>ZIP:</b> ", zip, "<br>",
"<b>Positive cases:</b> ", positive, "<br>",
"<b>Total tested:</b> ", total, "<br>",
"<b>COVID rate per 1,000:</b> ", round(covid_rate, 1), "<br>",
"<b>Elderly %:</b> ", round(elderly_pct, 1)
)
)
zip_centroids <- zip_centroids |>
mutate(
popup_food = paste0(
"<b>ZIP:</b> ", zip, "<br>",
"<b>Food stores:</b> ", n_food_stores
),
popup_health = paste0(
"<b>ZIP:</b> ", zip, "<br>",
"<b>Health facilities:</b> ", n_health_facilities
)
)Interactive Map with Three Layers
This interactive map includes three layers: COVID-19 rate, food store counts, and health facility counts.
pal_covid <- colorNumeric("YlOrRd", domain = zip_map$covid_rate)
interactive_map <- leaflet() |>
addProviderTiles(providers$CartoDB.Positron, group = "Basemap") |>
addPolygons(
data = zip_map,
fillColor = ~pal_covid(covid_rate),
fillOpacity = 0.7,
color = "white",
weight = 1,
popup = ~popup_covid,
label = ~paste("ZIP:", zip, "- COVID rate:", round(covid_rate, 1)),
group = "COVID Rate"
) |>
addCircleMarkers(
data = zip_centroids,
radius = ~sqrt(n_food_stores) + 2,
stroke = T,
weight = 1,
fillOpacity = 0.8,
popup = ~popup_food,
label = ~paste("ZIP:", zip, "- Food stores:", n_food_stores),
group = "Food Stores"
) |>
addCircleMarkers(
data = zip_centroids,
radius = ~sqrt(n_health_facilities) + 2,
stroke = T,
weight = 1,
fillOpacity = 0.8,
popup = ~popup_health,
label = ~paste("ZIP:", zip, "- Health facilities:", n_health_facilities),
group = "Health Facilities"
) |>
addLegend(
position = "bottomright",
pal = pal_covid,
values = zip_map$covid_rate,
title = "COVID rate per 1,000",
group = "COVID Rate"
) |>
addLayersControl(
overlayGroups = c("COVID Rate", "Food Stores", "Health Facilities"),
options = layersControlOptions(collapsed = F)
) |>
addScaleBar(position = "bottomleft")
interactive_mapSave Interactive Map as HTML
saveWidget(
interactive_map,
file = "section_08/lab3_interactive_map.html",
selfcontained = T
)Summary
This lab produced two static choropleth maps, a side-by-side comparison figure, and an interactive web map with multiple layers for exploratory visualization of COVID-19 and related demographic and service factors at the ZIP code level in NYC.