Skip to contents

Designed to represent the relative resilience of communities within an area, SVI not only provides valuable information for local officials/policymakers in planning and managing public health emergencies, but also offers an important neighborhood metric to facilitate public health research, such as the correlation between health outcome and SVI (and the socioeconomic/demographic factors involved) of a certain area.

For example, to investigate the possible correlation between the COVID-19-related deaths and SVI in Philadelphia at the ZCTA level, we could use findSVI to get ZCTA-level SVI for Philadelphia (with geometry) and join the result with ZCTA-level COVID-19 data for visualization and correlation analysis.

As usual, we start with loading all the packages needed.

ZCTA-level SVI in Philadelphia for 2020

For retrieving census data with geometry, we’ll use get_census_data() and get_svi() to obtain SVI. We’ll keep the GEOID(ZCTA) and SVI-related columns in the resulting SVI table.

phl_zcta_2020_geo_data <- get_census_data(
  2020,
  geography = "zcta",
  state = "PA",
  county = "Philadelphia",
  geometry = TRUE
)

phl_zcta_2020_geo_svi <- get_svi(2020,
  data = phl_zcta_data_geo_2020) %>%
  select(GEOID, contains("RPL_theme"))

phl_zcta_2020_geo_svi
#> Simple feature collection with 48 features and 6 fields
#> Geometry type: MULTIPOLYGON
#> Dimension:     XY
#> Bounding box:  xmin: -75.28027 ymin: 39.84987 xmax: -74.95578 ymax: 40.13799
#> Geodetic CRS:  NAD83
#> # A tibble: 48 × 7
#>    GEOID RPL_theme1 RPL_theme2 RPL_theme3 RPL_theme4 RPL_themes
#>    <chr>      <dbl>      <dbl>      <dbl>      <dbl>      <dbl>
#>  1 19102     0.0444     0.0222     0.0889      0.956     0.178 
#>  2 19103     0.0222     0.222      0.178       0.933     0.267 
#>  3 19104     0.689      0.2        0.556       0.978     0.689 
#>  4 19106     0          0.111      0.0444      0.711     0.0444
#>  5 19107     0.333      0.0667     0.4         1         0.378 
#>  6 19109    NA         NA         NA          NA        NA     
#>  7 19111     0.578      0.867      0.511       0.578     0.711 
#>  8 19112    NA         NA         NA          NA        NA     
#>  9 19114     0.178      0.422      0.156       0.111     0.2   
#> 10 19115     0.378      0.822      0.289       0.822     0.533 
#> # ℹ 38 more rows
#> # ℹ 1 more variable: geometry <MULTIPOLYGON [°]>

ZCTA-level COVID-19 data in Philadelphia

Disease-related data at the ZCTA level is usually not easily accessible for privacy reasons. Here, we’ll use data from the COVID-19 Health Inequities in Cities Dashboard, a great resource released by Drexel University’s Urban Health Collaborative and the Big Cities Health Coalition (BCHC). In addition to data available for download, the dashboard provides informative visualizations of COVID-19 related outcomes and inequities over time and across BCHC cities.

After downloading the raw data, we can select Philadelphia city and the variables of interest (hospitalization per 100k).

#source:https://github.com/Drexel-UHC/covid_inequities_project
#bchc_raw <- read_csv("../../byZCTA_bchc.csv")

phl_hosp <- bchc_raw %>%
  filter(city == "Philadelphia") %>%
  mutate(GEOID = paste(zcta), .after = zcta) %>%
  select(GEOID, hosp_per_100k)

glimpse(phl_hosp)
#> Rows: 46
#> Columns: 2
#> $ GEOID         <chr> "19102", "19103", "19104", "19106", "19107", "19111", "1…
#> $ hosp_per_100k <dbl> 194.288, 846.618, 1280.807, 405.019, 1555.831, 1133.411,…

Joining data for visualzation

Once we have the ZCTA-level SVI and COVID-19 data ready, we can join them together by GEOID(ZCTA), keeping the spatial information.

phl_svi_covid <- phl_hosp %>%
  left_join(phl_zcta_2020_geo_svi, by = "GEOID") %>%
  #although geometry sticky, after wrangling, class() become df
  st_as_sf(sf_column_name = "geometry")

phl_svi_covid %>% head(10)
#> Simple feature collection with 10 features and 64 fields
#> Geometry type: MULTIPOLYGON
#> Dimension:     XY
#> Bounding box:  xmin: -75.24074 ymin: 39.92578 xmax: -74.97371 ymax: 40.13799
#> Geodetic CRS:  NAD83
#> # A tibble: 10 × 65
#>    GEOID hosp_per_100k NAME                        geometry E_TOTPOP  E_HU  E_HH
#>    <chr>         <dbl> <chr>             <MULTIPOLYGON [°]>    <dbl> <dbl> <dbl>
#>  1 19102          194. ZCTA5 191… (((-75.16854 39.94663, -…     5335  3967  3365
#>  2 19103          847. ZCTA5 191… (((-75.18166 39.95151, -…    25113 18038 15765
#>  3 19104         1281. ZCTA5 191… (((-75.21367 39.96003, -…    52480 20416 16508
#>  4 19106          405. ZCTA5 191… (((-75.15476 39.94573, -…    13536  8597  7625
#>  5 19107         1556. ZCTA5 191… (((-75.16506 39.95361, -…    14689  8920  7823
#>  6 19111         1133. ZCTA5 191… (((-75.10653 40.04938, -…    67159 26250 24355
#>  7 19114         1336. ZCTA5 191… (((-75.03509 40.07462, -…    31448 13909 13032
#>  8 19115         1502. ZCTA5 191… (((-75.07456 40.08912, -…    35504 14976 14391
#>  9 19116         1226. ZCTA5 191… (((-75.04426 40.11571, -…    33976 13372 12850
#> 10 19118         1038. ZCTA5 191… (((-75.24022 40.0836, -7…    10826  4900  4651
#> # ℹ 58 more variables: E_POV150 <dbl>, E_UNEMP <dbl>, E_HBURD <dbl>,
#> #   E_NOHSDP <dbl>, E_UNINSUR <dbl>, E_AGE65 <dbl>, E_AGE17 <dbl>,
#> #   E_DISABL <dbl>, E_SNGPNT <dbl>, E_LIMENG <dbl>, E_MINRTY <dbl>,
#> #   E_MUNIT <dbl>, E_MOBILE <dbl>, E_CROWD <dbl>, E_NOVEH <dbl>,
#> #   E_GROUPQ <dbl>, EP_POV150 <dbl>, EP_UNEMP <dbl>, EP_HBURD <dbl>,
#> #   EP_NOHSDP <dbl>, EP_UNINSUR <dbl>, EP_AGE65 <dbl>, EP_AGE17 <dbl>,
#> #   EP_DISABL <dbl>, EP_SNGPNT <dbl>, EP_LIMENG <dbl>, EP_MINRTY <dbl>, …

Maps

To visualize their spatial pattern directly, we can plot the SVI and COVID-19 data in Philadelphia neighborhoods (ZCTAs) directly on maps:


covid_hosp <- phl_svi_covid %>%
  select(GEOID, geometry, hosp_per_100k) %>%
  drop_na() %>% 
  tm_shape(projection = sf::st_crs(26915))+
  tm_polygons("hosp_per_100k",
    style = "quantile",
    title = "Hospitalizations/100k")+
  tm_layout(title = "COVID-19 Hospitalizations",
    title.size = 1,
    title.position = c("left", "TOP"),
    legend.position = c("RIGHT", "bottom"),
    legend.title.size = 0.9,
    legend.width = 2)

svi <- phl_svi_covid %>%
  select(GEOID, geometry, RPL_themes) %>%
  drop_na() %>%
  tm_shape(projection = sf::st_crs(26915))+
  tm_polygons("RPL_themes",
    style = "quantile",
    title = "Overall SVI (2020)",
    palette = "Blues")+
  tm_layout(title = "Social Vulnerability",
    title.size = 1,
    title.position = c("left", "TOP"),
    legend.position = c("RIGHT", "bottom"),
    legend.title.size = 0.9,
    legend.width = 2)

plots <- list(svi, covid_hosp)

current.mode <- tmap_mode("plot")
#> tmap mode set to plotting
#> tmap mode set to plotting

tmap_arrange(
  plots, 
  nrow = 1,
  width = c(0.5, 0.5)
  ) 
SVI and COVID-19 hospitalizations in Philadelphia Neighborhoods (ZCTAs). COVID-19 data are cumulative till 8/2022.

SVI and COVID-19 hospitalizations in Philadelphia Neighborhoods (ZCTAs). COVID-19 data are cumulative till 8/2022.

Correlation

To look at possible correlation between SVI and COVID-19 hospitalizations, we’ll use scatter plot to visualize their relationships:

phl_svi_covid %>%
  ggplot(aes(x = RPL_themes, y = hosp_per_100k)) +
  geom_point()+
  labs(
    title = "Social Vulnerability and COVID-19 Hospitalizations \nin Philadelphia Neighborhoods",
    caption = "COVID-19 data are cumulative till 8/2022", 
    x = "Overall SVI (2020)",
    y = "Hospitalizations per 100k") +
  theme_bw()+
  theme(
    text = element_text(size = 13),
    plot.title = element_text(size = 14, face = "bold"),
    axis.title = element_text(size = 12, face = "bold"),
    legend.title = element_text(size = 13),
    plot.caption = element_text(size = 9, color = "#4D4948")
  )

With a correlation coefficient of 0.7030584, overall SVI is shown to be strongly associated with COVID-19 hospitalizations in Philadelphia neighborhoods, where higher hospitalization rate is found in neighborhoods with higher social vulnerability.

This is consistent with the story on the COVID-19 Health Inequities in Cities Dashboard, where they found most socially vulnerable neighborhoods have 129.4% higher hospitalizations per 100k compared to the least socially vulnerable neighborhoods.

Reference

Diez Roux, A., Kolker, J., Barber, S., Bilal, U., Mullachery, P., Schnake-Mahl, A., McCulley, E., Vaidya, V., Ran, L., Rollins, H., Furukawa, A., Koh, C., Sharaf, A., Dureja, K. (2021). COVID-19 Health Inequities In Cities Dashboard. Drexel University: Urban Health Collaborative. http://www.covid-inequities.info/.