Packages

R Spatial Lab Assignment /# 3

Task 1:

# Load data from Week 8
load("../spatial_data.RData")

# Check object
ls()
## [1] "final_zip_sf"   "nyc_zipcode_sf" "nys_health_sf"  "nys_retail_sf"
class(final_zip_sf)
## [1] "sf"         "data.frame"
# Quick preview
head(final_zip_sf)
## Simple feature collection with 6 features and 16 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: 986490.1 ymin: 168910.5 xmax: 1043042 ymax: 189382.9
## Projected CRS: NAD83 / New York Long Island (ftUS)
##   ZIPCODE  PO_NAME POPULATION COUNTY Positive Total zcta_cum.perc_pos
## 1   11436  Jamaica      18681 Queens      342   567             60.32
## 2   11213 Brooklyn      62426  Kings      972  1653             58.80
## 3   11212 Brooklyn      83866  Kings     1086  1793             60.57
## 4   11225 Brooklyn      56527  Kings      814  1359             59.90
## 5   11218 Brooklyn      72280  Kings     1163  1967             59.13
## 6   11226 Brooklyn     106132  Kings     1336  2170             61.57
##   food_store_count health_facility_count totPop elderlyPop whitePop blackPop
## 1                3                     1  22377       2456     1192    13972
## 2              118                     1  66602       7662    16483    43625
## 3              187                     1  73069       9518     5278    59541
## 4              103                     1  60958       7592    15300    40048
## 5              135                     1  67426       8144    40315     5045
## 6              236                     3 103729      12278    16173    69058
##   asianPop hispanicPop elderlyPct                       geometry
## 1     2626        3226   10.97556 POLYGON ((1038098 188138.4,...
## 2     1836        6738   11.50416 POLYGON ((1001614 186926.4,...
## 3     1330       13900   13.02604 POLYGON ((1011174 183696.3,...
## 4     2272        5293   12.45448 POLYGON ((995908.4 183617.6...
## 5    14695       11169   12.07843 POLYGON ((991997.1 176307.5...
## 6     5207       18244   11.83661 POLYGON ((994821.5 177865.7...
names(final_zip_sf)
##  [1] "ZIPCODE"               "PO_NAME"               "POPULATION"           
##  [4] "COUNTY"                "Positive"              "Total"                
##  [7] "zcta_cum.perc_pos"     "food_store_count"      "health_facility_count"
## [10] "totPop"                "elderlyPop"            "whitePop"             
## [13] "blackPop"              "asianPop"              "hispanicPop"          
## [16] "elderlyPct"            "geometry"

Map 1

# Create COVID map
map1 <- ggplot(final_zip_sf) +
  geom_sf(aes(fill = Positive)) +
  scale_fill_distiller(palette = "YlOrRd", name = "COVID Cases") +
  labs(
    title = "COVID-19 Cases Across NYC ZIP Codes",
    x = "Longitude",
    y = "Latitude"
  )

# Show map 1
map1

Map 2

# Create elderly population map
map2 <- ggplot(final_zip_sf) +
  geom_sf(aes(fill = elderlyPct)) +
  scale_fill_distiller(palette = "PuBu", name = "Elderly %") +
  labs(
    title = "Elderly Population (%) Across NYC ZIP Codes",
    x = "Longitude",
    y = "Latitude"
  )

# Show map 2
map2

Task 2:

# Label ZIP codes with highest COVID cases
label_sf <- final_zip_sf %>%
  filter(Positive > quantile(Positive, 0.9, na.rm = TRUE))

# Comparison map 1
map1_compare <- ggplot(final_zip_sf) +
  geom_sf(aes(fill = Positive)) +
  geom_sf_text(data = label_sf, aes(label = ZIPCODE), size = 2) +
  scale_fill_distiller(palette = "YlOrRd", name = "COVID Cases") +
  coord_sf(crs = st_crs(4326)) +
  labs(
    title = "COVID-19 Cases Across NYC ZIP Codes",
    x = "Longitude",
    y = "Latitude"
  )

# Comparison map 2
map2_compare <- ggplot(final_zip_sf) +
  geom_sf(aes(fill = elderlyPct)) +
  scale_fill_distiller(palette = "PuBu", name = "Elderly %") +
  coord_sf(crs = st_crs(4326)) +
  labs(
    title = "Elderly Population (%) Across NYC ZIP Codes",
    x = "Longitude",
    y = "Latitude"
  )

# Arrange maps side by side
ggarrange(map1_compare, map2_compare, ncol = 2, nrow = 1)

Task 3:

pal <- colorQuantile("YlOrRd", domain = final_zip_sf$Positive, n = 5)

map <- leaflet(final_zip_sf %>% st_transform(4326)) %>%
  addTiles() %>%
  addPolygons(
    fillColor = ~pal(Positive),
    fillOpacity = 0.8,
    color = "white",
    weight = 1,
    popup = ~paste(
      "ZIP:", ZIPCODE,
      "<br>COVID Cases:", Positive,
      "<br>Elderly %:", round(elderlyPct, 2)
    )
  ) %>%
  addLegend(
    pal = pal,
    values = ~Positive,
    title = "COVID Cases"
  )

map

Save HTML

htmlwidgets::saveWidget(map, "covid_map.html")