Data Sources: Saporito and Casey rely on 2010 5-year ACS estimates at the census block group level for the demographic data. The census block group serves as an approximation to a neighborhood. The race groups considered in the analysis are: non-Hispanic white, non-Hispanic Black, and Hispanic. They limit their analyses to census block groups in census places with more than 25,000 residents.
To measure green space, the authors draw upon satellite imagery from Landsat 7. They use the normalized difference vegetation index (NDVI) as an indicator of live green vegetation. NDVI is calculated from the visible and near-infrared light reflected by vegetation. Nearly all satellite Vegetation Indices employ this difference formula to quantify the density of plant growth on the Earth. Specifically, NDVI takes the subtracts the amount of visible radiation from the amount of near-infrared radiation and divides the difference by the sum of near-infrared radiation and visible radiation. Each pixel always takes on a number that ranges from -1 to +1. A zero means no vegetation and a value close to +1 (0.8 - 0.9) indicates the highest possible density of green vegetation.
I use data from the 2010 U.S. Census as well as satellite imagery from the USGS Landsat 7 Collection in the summer of 2010. The demographic variables are the same as those used in Saporito and Casey’s analysis. I also use both NDVI and EVI - Enhanced Vegetation Index - as measures of green space. EVI is a higher resolution measure than NDVI and will capture smaller patches of vegetation, such as a small cluster of trees or a pocket park like those found in cities.
## DATA IMPORTS-----------------------------------------------------------------
evi_all_decades_cbg <- readRDS(paste0(wd,
green_space_path,
"evi_all_years_on_2010cbg.rds"))
ndvi_2010_cbg <- read_csv(paste0(wd,
green_space_path,
"NDVI_2010_by_cbg10.csv"))
race_all_decades_cbg <- readRDS(paste0(wd,
census_data_path,
"/all_decades_race_on_cbg10.rds"))
counties <- get_decennial(geography = "county",
variables = "P001001",
year = 2010)
places <- get_decennial(geography = "place",
variables = "P001001",
year = 2010)
cbg_higher_geog_ids <- readRDS(paste0(shp_repo,
"/cbg_full_geog_identifers.rds"))
top100_cbsa <- read_csv(paste0(shp_repo,
"/top100_cbsa_2010bounds.csv"))
top100_places <- read_csv(paste0(shp_repo,
"/top100_places_2010bounds.csv"))#DATA PREPARATION---------------------------------------------------------------
## MERGE ALL GEOIDs
top100_cbsa <- top100_cbsa %>%
mutate(CBSAFP = as.character(CBSAFP)) %>%
select(-CBSA_NM) %>%
rename(top100_cbsa="top_100", cbsa_rank_pop_size="rank_pop_size")
top100_places <- top100_places %>%
select(-PLACE_NM) %>%
rename(top100_place="top_100", place_rank_pop_size="rank_pop_size")
counties <- counties %>%
select(PLACEFP_DERIVED="GEOID", PLACE_DERIVED_NM="NAME")
cbg_all_geoids <- cbg10 %>%
st_set_geometry(NULL) %>%
select(GEOID10, GISJOIN) %>%
left_join(., cbg_higher_geog_ids, by=c("GEOID10", "GISJOIN")) %>%
mutate(PLACEFP_DERIVED = if_else(is.na(PLACEID),
str_sub(GEOID10, end=5),
PLACEID)) %>%
left_join(., counties, by = "PLACEFP_DERIVED") %>%
mutate(PLACE_DERIVED_NM=coalesce(PLACE_DERIVED_NM, PLACE_NM)) %>%
relocate(PLACEFP_DERIVED, .after=PLACEID) %>%
relocate(PLACE_DERIVED_NM, .after=PLACE_NM) %>%
left_join(., top100_places, by=c("PLACEID"="PLACEFP")) %>%
left_join(., top100_cbsa, by=c("CBSAFP10"="CBSAFP"))
cbg_area <- cbg10 %>%
select(GEOID10) %>%
mutate(area = st_area(.),
area = as.numeric(set_units(area, mi^2))
) %>%
st_set_geometry(NULL)
ndvi_2010_cbg <- ndvi_2010_cbg %>%
select(GEOID10, GISJOIN, ndvi10 = "mean")
## MERGE ALL DEMOGRAPHIC & GREEN SPACE DATA
# PERCENT
race_pct_all_decades_cbg <- race_all_decades_cbg %>%
mutate_at(vars(WHITE_NH_10:HISPANIC_10), funs("percent" = (./TOT_POP_10) * 100)) %>%
mutate_at(vars(WHITE_NH_00:HISPANIC_00), funs("percent" = (./TOT_POP_00) * 100)) %>%
select(CBG_10, TOT_POP_10, TOT_POP_00, WHITE_NH_10_percent:HISPANIC_00_percent)
cbg_geog_evi_race_all_years <- cbg_all_geoids %>%
left_join(., evi_all_decades_cbg, by="GISJOIN") %>%
left_join(., race_pct_all_decades_cbg, by=c("GEOID10"="CBG_10")) %>%
left_join(., ndvi_2010_cbg, by=c("GEOID10"))
## DECADE SUBSETS
cbg_2010_data_subset <- cbg_geog_evi_race_all_years %>%
left_join(., cbg_area, by="GEOID10") %>%
mutate(POP_DENS_10 = TOT_POP_10/area) %>%
select(where(is_numeric) & contains("10") & !contains("top100"))
cbg_2000_data_subset <- cbg_geog_evi_race_all_years %>%
select(where(is_numeric) & contains("00") & !contains("top100"))Saporito and Casey do not include correlations in their paper, but they mention that the strongest correlation is between a place’s share of Hispanic residents and NDVI.
## CORRELATIONS
##2010
corr10 <- round(cor(cbg_2010_data_subset, use="complete.obs"), 1)
pmat10 <- cor_pmat(cbg_2010_data_subset, use="complete.obs")
ggcorrplot(corr10,
hc.order=T,
type="lower",
lab=T,
#p.mat = pmat10,
ggtheme = ggplot2::theme_gray,
colors = c("#6D9EC1", "white", "#E46726")) * NDVI and EVI have almost a perfect correlation (correlations are rounded in the graphic above). * NDVI/EVI is most positively correlated with the share of non-Hispanic whites in a neighborhood. * NDVI/EVI is most strongly negatively correlated with the share of Hispanic residents in a neighborhoods. * Increasing population density is associated with lower NDVI/EVI in a neighborhood. This negative association is second in magnitude to the Hispanic-NDVI/EVI association. * The share of Asian residents is also moderately negatively correlated with EVI/NDVI levels. * The percent of Black residents at the neighborhood level is very weakly negatively associated with EVI/NDVI levels.
Following Saporito and Casey’s analytic approach of assessing differences in green space within the same city, I construct their verdancy exposure for each racial group by city. This race-group verdancy exposure measure involves several steps. First, I multiplied the NDVI and EVI values by a factor of 10,000 for ease of comparison. Then, I computed race weighted vegetation exposure by multiplying EVI by each race group’s total count (code for each step is shown below).
#REPLICATE SAPORITO & CASEY DESCRIPTIVES----------------------------------------
# 1) WEIGHT EVI AND NDVI BY RACE
cbg_for_replication <- cbg_all_geoids %>%
left_join(., race_all_decades_cbg, by=c("GEOID10"="CBG_10")) %>%
left_join(., evi_all_decades_cbg, by="GISJOIN") %>%
left_join(., ndvi_2010_cbg, by="GEOID10") %>%
select(-TOT_POP_00:-evi00, -evi20)
cbg_for_rep_wtd_evi <- cbg_for_replication %>%
mutate(evi10 = evi10 * 10000) %>%
mutate_each(funs(.*evi10), ends_with("_10"))
cbg_for_rep_wtd_ndvi <- cbg_for_replication %>%
mutate(ndvi10 = ndvi10 * 10000) %>%
mutate_each(funs(.*ndvi10), ends_with("_10"))Next, I aggregated the race-weighted green space values from census block group to places. All census block groups that are not a city or census designated place were treated as county-non-incorporated places. Thus all census block groups were included in the analysis.
# 2) SUM RACE-WEIGHTED EVI TO PLACES
race_wtd_evi_summed_to_place <- cbg_for_rep_wtd_evi %>%
select(PLACEFP_DERIVED, PLACE_DERIVED_NM, WHITE_NH_10:HISPANIC_10) %>%
group_by(PLACEFP_DERIVED, PLACE_DERIVED_NM) %>%
summarize_all(sum, na.rm=T) %>%
rename_at(vars(-PLACEFP_DERIVED, -PLACE_DERIVED_NM),function(x) paste0(x,"_EVI"))
race_wtd_ndvi_summed_to_place <- cbg_for_rep_wtd_ndvi %>%
select(PLACEFP_DERIVED, PLACE_DERIVED_NM, WHITE_NH_10:HISPANIC_10) %>%
group_by(PLACEFP_DERIVED, PLACE_DERIVED_NM) %>%
summarize_all(sum, na.rm=T) %>%
rename_at(vars(-PLACEFP_DERIVED, -PLACE_DERIVED_NM),function(x) paste0(x,"_NDVI"))Finally, the race-weighted green-space measures were divided by the total number of residents belonging to the given race group to give an average racial-green space exposure. Racial groups’ green-space exposure levels were compared by subtracting two groups weighted verdancy scores.
Tables 1 and 2 below include summary statistics for the verdancy exposure differences.
# 3) DIVIDE RACE-WEIGHTED EVI BY TOTAL # RACE IN PLACE
race_counts_per_place <- cbg_for_replication %>%
select(PLACEFP_DERIVED, PLACE_DERIVED_NM, WHITE_NH_10:HISPANIC_10) %>%
group_by(PLACEFP_DERIVED, PLACE_DERIVED_NM) %>%
summarize_all(sum, na.rm=T)
green_exposure_race_diffs_by_place <- race_wtd_evi_summed_to_place %>%
left_join(., race_counts_per_place, by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM")) %>%
mutate(WHITE_NH_EVI = WHITE_NH_10_EVI/WHITE_NH_10,
BLACK_NH_EVI = BLACK_NH_10_EVI/BLACK_NH_10,
ASIAN_NH_EVI = ASIAN_NH_10_EVI/ASIAN_NH_10,
HISPANIC_EVI = HISPANIC_10_EVI/HISPANIC_10) %>%
mutate(BLACK_WHITE_EVI_DIFF = (BLACK_NH_EVI - WHITE_NH_EVI),
ASIAN_WHITE_EVI_DIFF = (ASIAN_NH_EVI - WHITE_NH_EVI),
HISPANIC_WHITE_EVI_DIFF = (HISPANIC_EVI - WHITE_NH_EVI),
BLACK_ASIAN_EVI_DIFF = (BLACK_NH_EVI - ASIAN_NH_EVI),
BLACK_HISPANIC_EVI_DIFF = (BLACK_NH_EVI - HISPANIC_EVI),
HISPANIC_ASIAN_EVI_DIFF = (HISPANIC_EVI - ASIAN_NH_EVI)
) %>%
select(PLACEFP_DERIVED, PLACE_DERIVED_NM, BLACK_WHITE_EVI_DIFF:HISPANIC_ASIAN_EVI_DIFF)
NDVI_exposure_race_diffs_by_place <- race_wtd_ndvi_summed_to_place %>%
left_join(., race_counts_per_place, by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM")) %>%
mutate(WHITE_NH_NDVI = WHITE_NH_10_NDVI/WHITE_NH_10,
BLACK_NH_NDVI = BLACK_NH_10_NDVI/BLACK_NH_10,
ASIAN_NH_NDVI = ASIAN_NH_10_NDVI/ASIAN_NH_10,
HISPANIC_NDVI = HISPANIC_10_NDVI/HISPANIC_10) %>%
mutate(BLACK_WHITE_NDVI_DIFF = (BLACK_NH_NDVI - WHITE_NH_NDVI),
ASIAN_WHITE_NDVI_DIFF = (ASIAN_NH_NDVI - WHITE_NH_NDVI),
HISPANIC_WHITE_NDVI_DIFF = (HISPANIC_NDVI - WHITE_NH_NDVI),
BLACK_ASIAN_NDVI_DIFF = (BLACK_NH_NDVI - ASIAN_NH_NDVI),
BLACK_HISPANIC_NDVI_DIFF = (BLACK_NH_NDVI - HISPANIC_NDVI),
HISPANIC_ASIAN_NDVI_DIFF = (HISPANIC_NDVI - ASIAN_NH_NDVI)
) %>%
select(PLACEFP_DERIVED, PLACE_DERIVED_NM, BLACK_WHITE_NDVI_DIFF:HISPANIC_ASIAN_NDVI_DIFF)
## EVI DIFFERENCES
green_exposure_race_diffs_by_place %>%
ungroup() %>%
select(BLACK_WHITE_EVI_DIFF:HISPANIC_WHITE_EVI_DIFF) %>%
datasummary((`Black-White`=BLACK_WHITE_EVI_DIFF) +
(`Asian-White`=ASIAN_WHITE_EVI_DIFF) +
(`Hispanic-White`=HISPANIC_WHITE_EVI_DIFF) ~
Mean + SD + Median + Min + Max,
fmt=2,
data = .,
title = 'Table 1: Summary of Race-Weighted EVI Differences',
output = 'markdown'
)| Mean | SD | Median | Min | Max | |
|---|---|---|---|---|---|
| Black-White | -102.33 | 268.83 | -30.74 | -2630.96 | 1648.96 |
| Asian-White | -71.17 | 207.71 | -22.14 | -1759.23 | 1846.60 |
| Hispanic-White | -82.78 | 201.95 | -24.23 | -1958.75 | 1614.99 |
Like Saporito and Casey’s study, white people tend to live in neighborhoods with a higher average verdancy score than the Black, Hispanic, and Asian people. The average differences are quite modest. However, there is substantial variation across places, especially for the Black-white comparison. For instance, in La Grange - a suburb of Chicago, Illinois - on average, whites live in neighborhoods with verdancy exposure that is 1,333 greater than Black residents. By contrast, in the city of Baltimore, black residents lived in neighborhoods with average verdancy scores that are 1,093 greater than whites.
## NDVI DIFFERENCES
NDVI_exposure_race_diffs_by_place %>%
ungroup() %>%
select(BLACK_WHITE_NDVI_DIFF:HISPANIC_WHITE_NDVI_DIFF) %>%
datasummary((`Black-White`=BLACK_WHITE_NDVI_DIFF) +
(`Asian-White`=ASIAN_WHITE_NDVI_DIFF) +
(`Hispanic-White`=HISPANIC_WHITE_NDVI_DIFF) ~
Mean + SD + Median + Min + Max,
fmt=2,
data = .,
title = 'Table 2: Summary of Race-Weighted NDVI Differences',
output = 'markdown'
)| Mean | SD | Median | Min | Max | |
|---|---|---|---|---|---|
| Black-White | -102.72 | 258.77 | -33.47 | -2750.47 | 1894.68 |
| Asian-White | -71.76 | 202.81 | -23.59 | -2426.23 | 1445.38 |
| Hispanic-White | -85.33 | 194.63 | -28.05 | -1893.84 | 1368.49 |
The differences based on NDVI, shown in Table 2, look very similar to those based on EVI in Table 1.
I have attempted to replicate the regression analyses featured in Saporito and Casey’s paper. I constructed segregation measures using the disimilarity index as they did for Hispanic-white and Black-white comparisons. I also included most of the control variables they used: population (in 1,000s), population density (1,000 people per square mile), percent Black, and percent Hispanic. I did not include the poverty segregation measure as a covariate because the 2010 census does not include poverty/socioeconomic status.
## CONSTRUCT COVARIATES
evi_place_pop_percent_covars <- cbg_for_replication %>%
group_by(PLACEFP_DERIVED, PLACE_DERIVED_NM) %>%
summarize(TOT_POP_10 = sum(TOT_POP_10, na.rm=F)) %>%
left_join(., race_counts_per_place, by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM")) %>%
select(PLACEFP_DERIVED, PLACE_DERIVED_NM, TOT_POP_10, BLACK_NH_10, HISPANIC_10) %>%
mutate_each(funs(./1000), ends_with("_10")) %>%
rename(TOT_POP_1000S="TOT_POP_10") %>%
mutate_at(c("BLACK_NH_10", "HISPANIC_10"), funs("percent" = (./TOT_POP_1000S) * 100))
evi_place_pop_percent_popdens_covars <- cbg_for_replication %>%
select(GEOID10:PLACE_DERIVED_NM) %>%
left_join(., cbg_area, by = "GEOID10") %>%
group_by(PLACEFP_DERIVED, PLACE_DERIVED_NM) %>%
summarize(area_sqmi = sum(area)) %>%
left_join(., evi_place_pop_percent_covars, by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM")) %>%
mutate(pop_dens_1000s=TOT_POP_1000S/area_sqmi)# CREATE DISIM. MEASURES--------------------------------------------------------
dissim_by_place <- cbg_for_replication %>%
select(PLACEFP_DERIVED, PLACE_DERIVED_NM, WHITE_NH_10:HISPANIC_10) %>%
rename_at(vars(-PLACEFP_DERIVED, -PLACE_DERIVED_NM),function(x) paste0(x,"_CBG")) %>%
left_join(., race_counts_per_place, by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM")) %>%
mutate(diss_bw=abs(BLACK_NH_10_CBG/BLACK_NH_10 - WHITE_NH_10_CBG/WHITE_NH_10),
diss_hw=abs(HISPANIC_10_CBG/HISPANIC_10 - WHITE_NH_10_CBG/WHITE_NH_10)) %>%
group_by(PLACEFP_DERIVED, PLACE_DERIVED_NM) %>%
summarise(dissim_bw= .5*sum(diss_bw, na.rm=T),
dissim_hw= .5*sum(diss_hw, na.rm=T)
) %>%
mutate_at(vars(dissim_bw:dissim_hw),
.funs = funs(. * 100))
full_place_data <- dissim_by_place %>%
left_join(.,
green_exposure_race_diffs_by_place,
by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM")) %>%
left_join(.,
evi_place_pop_percent_popdens_covars,
by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM"))
full_place_NDVI <- dissim_by_place %>%
left_join(.,
NDVI_exposure_race_diffs_by_place,
by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM")) %>%
left_join(.,
evi_place_pop_percent_popdens_covars,
by = c("PLACEFP_DERIVED", "PLACE_DERIVED_NM"))The regression analyses mirror Saporito and Casey’s results (see Table 2 in their paper). There are strong correlations between racial segregation and disparities in exposure to green space.
## Regressions Predicting Black-White Green Space Differences
m1_bw_evi <- lm(BLACK_WHITE_EVI_DIFF ~ dissim_bw, data=full_place_data)
m2_bw_evi <- lm(BLACK_WHITE_EVI_DIFF ~ dissim_bw +
TOT_POP_1000S +
BLACK_NH_10_percent +
HISPANIC_10_percent +
pop_dens_1000s,
data=full_place_data
)
m3_bw_ndvi <- lm(BLACK_WHITE_NDVI_DIFF ~ dissim_bw, data=full_place_NDVI)
m4_bw_ndvi <- lm(BLACK_WHITE_NDVI_DIFF ~ dissim_bw +
TOT_POP_1000S +
BLACK_NH_10_percent +
HISPANIC_10_percent +
pop_dens_1000s,
data=full_place_NDVI)
stargazer(m1_bw_evi, m2_bw_evi, m3_bw_ndvi, m4_bw_ndvi, type="html",
title="Table 3: Regression of Black-White Difference in Exposure to Green Space",
dep.var.caption = "Green Space Exposure",
covariate.labels=c("Black-White Dissimilarity",
"Total Population (1,000s)",
"Percent Black",
"Percent Hispanic",
"People per Square Mile (1,000s)"),
dep.var.labels=c("EVI", "NDVI"),
single.row=T)| Green Space Exposure | ||||
| EVI | NDVI | |||
| (1) | (2) | (3) | (4) | |
| Black-White Dissimilarity | -3.719*** (0.141) | -3.633*** (0.152) | -3.682*** (0.135) | -3.578*** (0.146) |
| Total Population (1,000s) | -0.135*** (0.020) | -0.134*** (0.020) | ||
| Percent Black | 1.177*** (0.168) | 1.021*** (0.161) | ||
| Percent Hispanic | 1.233*** (0.151) | 1.078*** (0.146) | ||
| People per Square Mile (1,000s) | -1.384 (0.848) | -1.201 (0.815) | ||
| Constant | 8.939* (5.097) | -13.125** (6.373) | 7.458 (4.894) | -12.021** (6.125) |
| Observations | 8,121 | 8,121 | 8,121 | 8,121 |
| R2 | 0.079 | 0.095 | 0.084 | 0.098 |
| Adjusted R2 | 0.079 | 0.094 | 0.083 | 0.097 |
| Residual Std. Error | 258.007 (df = 8119) | 255.821 (df = 8115) | 247.734 (df = 8119) | 245.853 (df = 8115) |
| F Statistic | 696.341*** (df = 1; 8119) | 170.327*** (df = 5; 8115) | 740.578*** (df = 1; 8119) | 176.128*** (df = 5; 8115) |
| Note: | p<0.1; p<0.05; p<0.01 | |||
For every one point increase in Black-white segregation, the Black-white difference in exposure to green space favors whites by about 3.7 points. This strength of the relationship holds across both EVI and NDVI as well as when the control variables are added.
## Regressions Predicting Black-White Green Space Differences
m1_hw_evi <- lm(HISPANIC_WHITE_EVI_DIFF ~ dissim_hw, data=full_place_data)
m2_hw_evi <- lm(HISPANIC_WHITE_EVI_DIFF ~ dissim_hw +
TOT_POP_1000S +
BLACK_NH_10_percent +
HISPANIC_10_percent +
pop_dens_1000s,
data=full_place_data
)
m3_hw_ndvi <- lm(HISPANIC_WHITE_NDVI_DIFF ~ dissim_hw, data=full_place_NDVI)
m4_hw_ndvi <- lm(HISPANIC_WHITE_NDVI_DIFF ~ dissim_hw +
TOT_POP_1000S +
BLACK_NH_10_percent +
HISPANIC_10_percent +
pop_dens_1000s,
data=full_place_NDVI)
stargazer(m1_hw_evi, m2_hw_evi, m3_hw_ndvi, m4_hw_ndvi, type="html",
title="Table 4: Regression of Hispanic-White Difference in Exposure to Green Space",
dep.var.caption = "Green Space Exposure",
covariate.labels=c("Hispanic-White Dissimilarity",
"Total Population (1,000s)",
"Percent Black",
"Percent Hispanic",
"People per Square Mile (1,000s)"),
dep.var.labels=c("EVI","NDVI"),
single.row=T)| Green Space Exposure | ||||
| EVI | NDVI | |||
| (1) | (2) | (3) | (4) | |
| Hispanic-White Dissimilarity | -4.918*** (0.137) | -4.929*** (0.147) | -5.047*** (0.130) | -5.014*** (0.140) |
| Total Population (1,000s) | -0.139*** (0.015) | -0.141*** (0.014) | ||
| Percent Black | 0.902*** (0.122) | 0.820*** (0.117) | ||
| Percent Hispanic | 1.044*** (0.108) | 0.849*** (0.103) | ||
| People per Square Mile (1,000s) | -3.703*** (0.594) | -3.491*** (0.567) | ||
| Constant | 19.101*** (3.504) | 12.210*** (4.198) | 19.232*** (3.341) | 14.377*** (4.006) |
| Observations | 8,229 | 8,229 | 8,229 | 8,229 |
| R2 | 0.136 | 0.159 | 0.155 | 0.176 |
| Adjusted R2 | 0.136 | 0.159 | 0.154 | 0.176 |
| Residual Std. Error | 187.699 (df = 8227) | 185.209 (df = 8223) | 178.966 (df = 8227) | 176.716 (df = 8223) |
| F Statistic | 1,297.969*** (df = 1; 8227) | 311.976*** (df = 5; 8223) | 1,503.946*** (df = 1; 8227) | 351.452*** (df = 5; 8223) |
| Note: | p<0.1; p<0.05; p<0.01 | |||
For every one point increase in Hispanic-white segregation, the Hispanic-white difference in exposure to green space favors whites by 4.9 points. As above, the disparity holds across verdancy measures and with the control variables added. Like Saporito and Casey’s regressions, the disparity is greater in magnitude for Hispanics than Blacks.
Takeaways for next analyses: