SVI and health outcome
SVI and COVID-19 hospitalizations in Philadelphia
Source:vignettes/articles/svi-covid.Rmd
svi-covid.Rmd
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)
)
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/.