Zip Code Data

COVID, Retail, and Health

Quarto markdown is different from R markdown in terms of chunk options. See chunk options at Quarto website.

popData <- popData %>%
  mutate(BIPOC_percentage = round((blackPop / totPop) * 100, 1))


plotbipoc <- ggplot(popData) +
  geom_sf(aes(fill = BIPOC_percentage)) +
  scale_fill_distiller(palette = "OrRd", direction = 1, name = "Percentage BIPOC") +
  coord_sf(xlim = c(-74.4, -73.65), default_crs = sf::st_crs(4326) ) +
  labs(x = "Longitude", y = "Latitude",
       title = "BIPOC Population Percentage in Each NYC ZIP Code") +
  theme_minimal() +
  theme(legend.position = "right") 

# Ensure COVID_DEATH_COUNT is numeric
nyc_covid_data_sf_merged$COVID_DEATH_COUNT <- as.numeric(nyc_covid_data_sf_merged$COVID_DEATH_COUNT)

# Compute the quartile breaks based on COVID_DEATH_COUNT
breaks_qt <- quantile(nyc_covid_data_sf_merged$COVID_DEATH_COUNT, probs = 0:4/4, na.rm = TRUE)

# Create a new column that categorizes COVID_DEATH_COUNT into quartiles
nyc_covid_data_sf_merged <- mutate(nyc_covid_data_sf_merged, 
                                   covid_death_cat = cut(COVID_DEATH_COUNT, 
                                                         breaks = breaks_qt, 
                                                         include.lowest = TRUE, 
                                                         labels = c("Q1 (Low)", "Q2", "Q3", "Q4 (High)"),
                                                         dig.lab = 4, 
                                                         digits = 1))

# Plot using ggplot2
plotcaserate <- ggplot(nyc_covid_data_sf_merged) + 
  geom_sf(aes(fill = covid_death_cat)) +
  scale_fill_brewer(palette = "OrRd", direction = 1, name = 'COVID-19 Death Count') +  # Color scale for quartiles
  coord_sf(xlim = c(-74.4, -73.65), default_crs = sf::st_crs(4326)) +
  labs(x = 'Longitude', y = 'Latitude', 
       title = 'COVID-19 Death Count by ZIP Code') +
  theme_minimal() +
  theme(legend.position = "right")




plotelder <- ggplot(popData) +
  geom_sf(aes(fill = elderly_pct)) +
  scale_fill_distiller(palette = "OrRd", direction = 1, name = "% Elderly") +
  coord_sf(xlim = c(-74.4, -73.65), default_crs = sf::st_crs(4326) ) +
  labs(x = "Longitude", y = "Latitude",
       title = "Elderly Population Percentage in Each NYC ZIP Code") +
  theme_minimal() +
  theme(legend.position = "right") 


plotcaserate

plotbipoc

plotelder

#  Code for Elderly Pop


popData <- popData %>%
  mutate(elderly_pct = round((elderlyPop / totPop) * 100, 1))

install.packages("gridExtra")
## 
## The downloaded binary packages are in
##  /var/folders/dm/sj7012g577qdnv76gj2_szd80000gp/T//RtmpxCyu9L/downloaded_packages
library(gridExtra)
## 
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
## 
##     combine
require(ggmap)


install.packages("patchwork")
## 
## The downloaded binary packages are in
##  /var/folders/dm/sj7012g577qdnv76gj2_szd80000gp/T//RtmpxCyu9L/downloaded_packages
require(patchwork)
## Loading required package: patchwork
## 
## Attaching package: 'patchwork'
## The following object is masked from 'package:raster':
## 
##     area
plotelder + plotcaserate + plot_layout(ncol = 2, nrow = 1, widths = 2, heights = 2)

library(leaflet) 

st_crs(nyc_covid_data_sf_merged)
## Coordinate Reference System:
##   User input: NAD83 / New York Long Island (ftUS) 
##   wkt:
## PROJCRS["NAD83 / New York Long Island (ftUS)",
##     BASEGEOGCRS["NAD83",
##         DATUM["North American Datum 1983",
##             ELLIPSOID["GRS 1980",6378137,298.257222101,
##                 LENGTHUNIT["metre",1]]],
##         PRIMEM["Greenwich",0,
##             ANGLEUNIT["degree",0.0174532925199433]],
##         ID["EPSG",4269]],
##     CONVERSION["SPCS83 New York Long Island zone (US survey foot)",
##         METHOD["Lambert Conic Conformal (2SP)",
##             ID["EPSG",9802]],
##         PARAMETER["Latitude of false origin",40.1666666666667,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8821]],
##         PARAMETER["Longitude of false origin",-74,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8822]],
##         PARAMETER["Latitude of 1st standard parallel",41.0333333333333,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8823]],
##         PARAMETER["Latitude of 2nd standard parallel",40.6666666666667,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8824]],
##         PARAMETER["Easting at false origin",984250,
##             LENGTHUNIT["US survey foot",0.304800609601219],
##             ID["EPSG",8826]],
##         PARAMETER["Northing at false origin",0,
##             LENGTHUNIT["US survey foot",0.304800609601219],
##             ID["EPSG",8827]]],
##     CS[Cartesian,2],
##         AXIS["easting (X)",east,
##             ORDER[1],
##             LENGTHUNIT["US survey foot",0.304800609601219]],
##         AXIS["northing (Y)",north,
##             ORDER[2],
##             LENGTHUNIT["US survey foot",0.304800609601219]],
##     USAGE[
##         SCOPE["Engineering survey, topographic mapping."],
##         AREA["United States (USA) - New York - counties of Bronx; Kings; Nassau; New York; Queens; Richmond; Suffolk."],
##         BBOX[40.47,-74.26,41.3,-71.8]],
##     ID["EPSG",2263]]
nyc_covid_data_sf_merged <- st_transform(nyc_covid_data_sf_merged, 4326)

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


p_tip <- paste0("Neighborhood: ", covid_data$NEIGHBORHOOD_NAME);

p_popup <- paste0("<strong>COVID Death Density: </strong>", 
                  nyc_covid_data_sf_merged$COVID_DEATH_COUNT%>%round(3)%>%format(nsmall = 3), 
                  " /sqkm",
                  " <br/>",
                  "<strong> Number of Deaths: </strong>",
                  nyc_covid_data_sf_merged$COVID_DEATH_RATE,
                  sep="")



htmlMap <- leaflet(nyc_covid_data_sf_merged) %>%
  addPolygons(
    stroke = FALSE, 
    fillColor = ~pal_fun(COVID_DEATH_COUNT),
    fillOpacity = 0.8, 
    smoothFactor = 0.5,
    label = p_tip,
    popup = p_popup) %>%
    addTiles() %>%
     addLegend("bottomright",  # location
            pal=pal_fun,    # palette function
            opacity = 0.8, # default is .5. Make it the same as the map
            values=~COVID_DEATH_COUNT,  # value to be passed to palette function
            title = 'COVID Deaths <br> by Zip Code <br> - Quantile')
htmlMap