Local Spatial Autocorrelation
GEOG 4/597: Advanced Spatial Quantitative
Analysis, Winter 2023
Jackson Voelkel | Portland State
University
Global Example
summary(uhi$uhi_c)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 28.30 30.36 30.74 30.66 31.06 31.92hist(uhi$uhi_c,breaks=20, main="Histogram of Evening UHI by Block Group", xlab = "Evening UHI")
uhi %>% ggplot(aes(x=uhi_c)) +
geom_histogram(color ='white', lwd=0.2) +
xlab("Average Temperature (Evening)") +
ylab("Frequency") +
theme_minimal()
Local Getis-Ord G Statistic
suppressPackageStartupMessages(require(spdep))
suppressPackageStartupMessages(require(sf))
suppressPackageStartupMessages(require(tidycensus))
suppressPackageStartupMessages(require(dplyr))
tract <- tidycensus::get_acs(
geography = 'tract',
variables = c("Households_" = 'B22010_001', "SNAP_" = 'B22010_002'),
state="CA",
county='Los Angeles',
geometry = TRUE,
output = 'wide') %>%
st_transform(4326) %>%
mutate(
pct_snap = SNAP_E/Households_E
# pct_snap = ifelse(is.nan(pct_snap), 0, pct_snap)
# pct_snap = ifelse(is.na(pct_snap), 0, pct_snap)
) %>%
filter(
as.vector(sf::st_area(.)) > 0,
st_coordinates(st_centroid(.))[,2] > 33.667,
!is.na(pct_snap),
!is.nan(pct_snap) # Lets remove the NA values
)
suppressMessages(tmap_mode('plot'))
tm_shape(tract) +
tm_polygons(
col='pct_snap',
palette = 'Blues',
border.col = 'white',
lwd=0.3,
style = 'quantile',
n=5,
contrast=c(0,1),
title = "Share of Households (0 to 1)") +
tm_shape(tract %>% st_union()) + # adding a nice border!
tm_polygons(col = NA, alpha = 0, border.col = 'black') +
tm_layout(main.title = "Households on Food Assistance (SNAP)")
tract <- tract %>% filter(st_coordinates(st_centroid(.))[,2] < 34.3)
tm_shape(tract) +
tm_polygons(
col='pct_snap',
palette = "Blues",
border.col = 'white',
lwd=0.3,
style = 'quantile',
contrast=c(0,1),
n=5,
title = "Share of Households (0 to 1)") +
tm_shape(tract %>% st_union()) + # adding a nice border!
tm_polygons(col = NA, alpha = 0, border.col = 'black') +
tm_layout(main.title = "Households on Food Assistance (SNAP)")
Assessing Significance
# Apply our colors on the fly to the 'col' column
tract %>% mutate(G=make_g_bin(tract_local_g)) %>%
tm_shape() +
tm_polygons(
col = "G",
palette = g_palette,
border.col = 'grey80',
lwd=0.4,
title = "Local Getis-Ord G") +
tm_layout(main.title = "Households on Food Assistance (SNAP)")
# DO NOT USE THIS IN YOUR LAB! THIS IS JUST AN EXAMPLE!
# You'll need to recreate the one in a few steps!
suppressMessages(tmap_mode('plot'))
tract %>% mutate(I=make_i_bin(tract$lisa_p)) %>%
tm_shape() +
tm_polygons(
col = "I",
palette = c('grey','cornsilk','lightgreen'),
border.col = 'grey60',
lwd=0.4,
title = "Local Moran's I") +
tm_layout(main.title = "Households on Food Assistance (SNAP)")
Creating a lag vector
tract$lag_snap <- spdep::lag.listw(x = tract_W, var = tract$pct_snap)
# Does this look familiar?
moran_plot = function(var, w_var, xlab, ylab, main_title = NULL, normalize = T){
if(missing(main_title)){
main_title = "Moran's I Plot"
}
if(normalize){
var = scale(var) %>% as.vector
w_var = scale(w_var) %>% as.vector()
}
plot_dat = data.frame(var,w_var,stringsAsFactors = F) %>%
mutate(quad = case_when(
var > 0 & w_var > 0 ~ "HH",
var > 0 & w_var < 0 ~ "HL",
var < 0 & w_var < 0 ~ "LL",
var < 0 & w_var > 0 ~ "LH",
TRUE ~ 'Unknown'))
plot_dat %>% ggplot() +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
geom_point(aes(x=var,y=w_var,color=quad)) +
xlab(xlab) +
ylab(ylab) +
geom_smooth(
mapping = aes(x=var,y=w_var),
formula = 'y ~ x',
method = 'lm',
color="black",
size=0.5,
se=F) +
scale_color_manual(
name = "",
values = c('Unknown' = 'grey','HH'="firebrick",'HL'="brown1",'LL'="navy",'LH'="cornflowerblue")) +
ggtitle(main_title)
}
moran_plot(
var = tract$pct_snap,
w_var = tract$lag_snap,
xlab = "Households on SNAP (%)",
ylab = "Lagged Households on SNAP (%)")
tract %>% tm_shape() +
tm_polygons(
col = "lag_snap",
palette = '-RdBu',
border.col = 'grey80',
lwd=0.4,
style = 'cont',
title = "Lagged Households on SNAP (%)") +
tm_layout(main.title = "Households on Food Assistance (SNAP)")
tract %>% mutate(quad = make_i_quads(v = pct_snap, w_v = lag_snap)) %>%
tm_shape() +
tm_polygons(
col = "quad",
palette = i_palette,
border.col = 'grey80',
lwd=0.4,
title = "Local Moran's I")
tract %>% mutate(quad = make_i_quads(v = pct_snap, w_v = lag_snap,p = lisa_p)) %>%
tm_shape() +
tm_polygons(
col = "quad",
palette = i_palette,
border.col = 'grey70',
lwd=0.4,
title = "Local Moran's I") +
tm_layout(main.title = "Local Moran's I: SNAP Benefit Usage in LA County")
Bivariate Moran’s I Plot
moran_plot(
var = oak$mhi,
w_var = lag_edu,
ylab = "Lagged Adults 25+ with at Least a Bachelors Degree (%)",
xlab = "Median Household Income ($)",
main_title = "Bivariate Moran's I Plot")
We are observing positive spatial autocorrelation between “education” and lagged income!
Lecture: Spatial Regression
