Tutorial: Spatial Autocorrelation — Completed
July 21, 2026
knitr::opts_chunk$set(echo = TRUE, message = FALSE, warning = FALSE)Tutorial
Spatial autocorrelation describes whether nearby locations have similar values. Moran’s I summarizes the overall pattern, while Local Indicators of Spatial Autocorrelation (LISA) identify where significant clusters and spatial outliers occur.
The analysis below uses 2023 American Community Survey data for North Carolina census tracts and chronic-absence data for North Carolina public schools.
Reading in the data
library(tmap)
library(tidyverse)
library(spdep)
library(sf)
# Both datasets use NAD83 / North Carolina (EPSG:2264), measured in feet.
tracts_sf <- st_read(
"https://drive.google.com/uc?export=download&id=1mUxk0R_lqcfuJppw_A-1vy-pcB9XdCnS",
quiet = TRUE
)
schools_sf <- st_read(
"https://drive.google.com/uc?export=download&id=1_X3s6sTw5zeIuXa6wDQ61r8Px7_tI7Nr",
quiet = TRUE
)
set.seed(734)
tmap_mode("plot")Calculating global autocorrelation
Moran’s I ranges roughly from -1 to 1. Positive values show that nearby locations tend to have similar values, values near zero suggest a spatially random pattern, and negative values show dispersion.
Global autocorrelation in median household income
tm_shape(tracts_sf) +
tm_polygons(
fill = "med_hh_inc",
fill.scale = tm_scale_intervals(style = "fisher"),
col_alpha = 0
)
Q1: Just by looking at the values on the map, do you see evidence of a spatial pattern?
Yes. Higher-income tracts appear grouped around the Charlotte, Triangle, and Triad metropolitan areas, while many lower-income tracts are grouped in more rural parts of the state. The values do not look randomly scattered.
# Queen neighbors share either a boundary or a corner.
nb <- poly2nb(tracts_sf, queen = TRUE)
nbw <- nb2listw(nb, style = "W", zero.policy = TRUE)
gmoran <- moran.mc(
tracts_sf$med_hh_inc,
nbw,
nsim = 99,
zero.policy = TRUE
)
gmoranMonte-Carlo simulation of Moran I
statistic = 0.58868, observed rank = 100, p-value = 0.01
alternative hypothesis: greaterThe Moran’s I of about 0.59 shows strong positive spatial autocorrelation. Nearby tracts tend to have similar household incomes, and the permutation p-value of 0.01 shows that this pattern is unlikely to be random.
Global autocorrelation in school absence rate
tm_shape(schools_sf) +
tm_dots(
fill = "pre_covid_pct",
fill.scale = tm_scale_intervals(style = "fisher"),
col_alpha = 0
)
Q2: Just by looking at the values on the map, do you see evidence of a spatial pattern?
There is some pattern, but it is much less obvious than the income pattern. Schools with similar absence rates form a few pockets, especially around metropolitan areas and parts of eastern and northern North Carolina, but there is also a lot of mixing.
# A 26,400-foot distance band is five miles in this projected CRS.
nb_school <- dnearneigh(schools_sf, d1 = 0, d2 = 26400)
nbw_school <- nb2listw(nb_school, style = "W", zero.policy = TRUE)
gmoran_school <- moran.mc(
schools_sf$pre_covid_pct,
nbw_school,
nsim = 99,
zero.policy = TRUE
)
gmoran_schoolMonte-Carlo simulation of Moran I
statistic = 0.21811, observed rank = 100, p-value = 0.01
alternative hypothesis: greaterQ3: What does our Moran’s I value tell us about the spatial structure of our data?
The Moran’s I of about 0.22 shows weak positive spatial autocorrelation. Nearby schools tend to have somewhat similar pre-COVID absence rates, but the clustering is not especially strong. The p-value of 0.01 still shows that the pattern is statistically significant.
Local Indicators of Spatial Autocorrelation (LISA)
Median household income
locali <- localmoran_perm(
tracts_sf$med_hh_inc,
nbw,
zero.policy = TRUE
) |>
as_tibble() |>
set_names(c(
"local_i", "exp_i", "var_i", "z_i", "p_i", "p_i_sim",
"p_i_sim_folded", "skewness", "kurtosis"
))
tracts_sf <- tracts_sf |>
bind_cols(locali) |>
mutate(
incz = as.numeric(scale(med_hh_inc)),
inczlag = lag.listw(nbw, incz, zero.policy = TRUE),
lisa_cluster = case_when(
p_i >= 0.05 ~ "Not significant",
incz > 0 & local_i > 0 ~ "High-high",
incz > 0 & local_i < 0 ~ "High-low",
incz < 0 & local_i > 0 ~ "Low-low",
incz < 0 & local_i < 0 ~ "Low-high"
)
)
tm_shape(tracts_sf) +
tm_polygons(fill = "lisa_cluster", col_alpha = 0)
Q4: What is the spatial pattern of clustering in the median household income variable?
The LISA map shows clear high-high income clusters around Charlotte, the Triangle, and the Triad, with a few smaller metropolitan clusters elsewhere. Low-low clusters are more common in rural northeastern and south-central North Carolina and in a few western areas. High-low and low-high spatial outliers are much less common and usually sit along the edges of the larger clusters.
Pre-COVID school chronic absence
locali_school <- localmoran_perm(
schools_sf$pre_covid_pct,
nbw_school,
zero.policy = TRUE
) |>
as_tibble() |>
set_names(c(
"local_i", "exp_i", "var_i", "z_i", "p_i", "p_i_sim",
"p_i_sim_folded", "skewness", "kurtosis"
))
schools_sf <- schools_sf |>
bind_cols(locali_school) |>
mutate(
absz = as.numeric(scale(pre_covid_pct)),
abszlag = lag.listw(nbw_school, absz, zero.policy = TRUE),
lisa_cluster = case_when(
p_i >= 0.05 ~ "Not significant",
absz > 0 & local_i > 0 ~ "High-high",
absz > 0 & local_i < 0 ~ "High-low",
absz < 0 & local_i > 0 ~ "Low-low",
absz < 0 & local_i < 0 ~ "Low-high"
)
)
tm_shape(schools_sf) +
tm_dots(fill = "lisa_cluster", col_alpha = 0)
Q5: What is the spatial pattern of clustering in the school absence variable?
The school pattern is more localized than the income pattern. Low-low absence clusters appear around the Charlotte, Triangle, and Triad areas. High-high clusters are smaller and more scattered, including pockets in northern, eastern, and coastal parts of the state. Most schools are not part of a statistically significant local cluster, which fits the weaker global Moran’s I.
Mini Challenge 1: Average commuting time
tm_shape(tracts_sf) +
tm_polygons(
fill = "av_commute",
fill.scale = tm_scale_intervals(style = "fisher"),
col_alpha = 0
)
gmoran_commute <- moran.mc(
tracts_sf$av_commute,
nbw,
nsim = 99,
zero.policy = TRUE
)
gmoran_commuteMonte-Carlo simulation of Moran I
statistic = 0.46659, observed rank = 100, p-value = 0.01
alternative hypothesis: greaterAverage commuting time has moderate-to-strong positive spatial autocorrelation. The Moran’s I of about 0.47 means neighboring tracts often have similar commuting times, and the p-value of 0.01 shows that the clustering is statistically significant.
locali_commute <- localmoran_perm(
tracts_sf$av_commute,
nbw,
zero.policy = TRUE
) |>
as_tibble() |>
set_names(c(
"local_i", "exp_i", "var_i", "z_i", "p_i", "p_i_sim",
"p_i_sim_folded", "skewness", "kurtosis"
))
tracts_commute <- tracts_sf |>
select(GEOID, av_commute, geometry) |>
bind_cols(locali_commute) |>
mutate(
commute_z = as.numeric(scale(av_commute)),
lisa_cluster = case_when(
p_i_sim_folded >= 0.05 ~ "Not significant",
commute_z > 0 & local_i > 0 ~ "High-high",
commute_z > 0 & local_i < 0 ~ "High-low",
commute_z < 0 & local_i > 0 ~ "Low-low",
commute_z < 0 & local_i < 0 ~ "Low-high"
)
)
tm_shape(tracts_commute) +
tm_polygons(fill = "lisa_cluster", col_alpha = 0)
The LISA map shows high-high commuting-time clusters across parts of northern and eastern North Carolina, along with several outer suburban and rural commuting belts. Low-low clusters are concentrated in and near major employment centers such as Charlotte, the Triangle, the Triad, and Asheville. A small number of high-low and low-high tracts act as local spatial outliers.
Mini Challenge 2: Chronic school absences after COVID
tm_shape(schools_sf) +
tm_dots(
fill = "after_covid_pct",
fill.scale = tm_scale_intervals(style = "fisher"),
col_alpha = 0
)
gmoran_post <- moran.mc(
schools_sf$after_covid_pct,
nbw_school,
nsim = 99,
zero.policy = TRUE
)
gmoran_postMonte-Carlo simulation of Moran I
statistic = 0.33041, observed rank = 100, p-value = 0.01
alternative hypothesis: greaterThe post-COVID absence rate has moderate positive spatial autocorrelation. Its Moran’s I of about 0.33 is higher than the pre-COVID value of 0.22, suggesting that school absence rates became more spatially clustered after COVID. The p-value of 0.01 shows that this pattern is statistically significant.
locali_post <- localmoran_perm(
schools_sf$after_covid_pct,
nbw_school,
zero.policy = TRUE
) |>
as_tibble() |>
set_names(c(
"local_i", "exp_i", "var_i", "z_i", "p_i", "p_i_sim",
"p_i_sim_folded", "skewness", "kurtosis"
))
schools_post <- schools_sf |>
select(agency_code, school_nam, after_covid_pct, geometry) |>
bind_cols(locali_post) |>
mutate(
post_z = as.numeric(scale(after_covid_pct)),
lisa_cluster = case_when(
p_i_sim_folded >= 0.05 ~ "Not significant",
post_z > 0 & local_i > 0 ~ "High-high",
post_z > 0 & local_i < 0 ~ "High-low",
post_z < 0 & local_i > 0 ~ "Low-low",
post_z < 0 & local_i < 0 ~ "Low-high"
)
)
tm_shape(schools_post) +
tm_dots(fill = "lisa_cluster", col_alpha = 0)
High-high post-COVID absence clusters appear in several northern, eastern, and south-central parts of the state. Low-low clusters are especially visible around Charlotte, the Triangle, and the Triad, with a few smaller western clusters. Many schools remain statistically insignificant, but the significant clusters are more noticeable than they were before COVID, which is consistent with the larger global Moran’s I.
Overall takeaway
All four variables show statistically significant positive spatial autocorrelation. Median household income is the most strongly clustered of the variables examined, followed by average commuting time, post-COVID school absence, and pre-COVID school absence. The LISA results add important detail by showing that the statewide Moran’s I values are produced by specific local clusters rather than one uniform pattern across North Carolina.