In JosiahParry/sfweight: What the Package Does (One Line, Title Case)

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.path = "man/figures/README-",
  out.width = "100%"
)

NOTE: this package is under active and experimental development. Functions are likely to change.

sfweight is an opinionated translation of the wonderful spdep package. The goal is to provide a streamlined method of doing spatial statistics that works with sf objects, data frames, and the tidyverse. spdep is more flexible with the types of input objects and a bit more idiosyncratic in the syntax that is used.

Installation

You can install the development version from GitHub with

remotes::install_github("Josiahparry/sfweight")

Motivating examples

Spatial OLS

We can fit a spatial Durbin model by calculating spatially lagged predictors.

library(sfweight)
library(tidyverse)

acs_lagged <- acs %>% 
  mutate(nb = st_contiguity(geometry),
         wts = st_weights(nb),
         trans_lag = st_lag(by_pub_trans, nb, wts),
         bach_lag = st_lag(bach, nb, wts))


durbin_lm <- lm(med_house_income ~ trans_lag + by_pub_trans + bach_lag + bach, 
   data = acs_lagged)

broom::tidy(durbin_lm)

Local Autocorrelation

We can create a Moran plot by creating a spatially lagged variable. Additionally the function categorize_lisa() will categorize high-high, high-low, etc., groupings of these variables.

acs_lagged %>% 
  mutate(inc_lag = st_lag(med_house_income, nb, wts),
         lisa_group = categorize_lisa(med_house_income, inc_lag)) %>% 
  ggplot(aes(med_house_income, inc_lag, color = lisa_group)) +
  geom_vline(aes(xintercept = mean(med_house_income)), lty = 2, alpha = 1/3) +
  geom_hline(aes(yintercept = mean(inc_lag)), lty = 2, alpha = 1/3) + 
  geom_point() +
  labs(title = "Moran Plot",
       y = "Med. HH Income Spatial Lag",
       x = "Median Household Income") +
  theme_minimal() +
  scale_x_continuous(labels = scales::dollar) + 
  scale_y_continuous(labels = scales::dollar)

We can also calculate the Local Moran's I for each observation using the function local_moran() this will create a dataframe column containing the I, expected I, variance, Z-value, and P-value for each observation. You can extract this using tidyr::unnest().

acs %>% 
  mutate(nb = st_contiguity(geometry),
         wt = st_weights(nb),
    lisa = local_moran(med_house_income, nb, wt)) %>% 
  unnest(lisa) %>% 
  ggplot(aes(fill = lisa_category)) + 
  geom_sf(color = "black", lwd = 0.25) +
  scale_fill_manual(values = c("HH" = "#528672","LL" = "#525586", "Insignificant" = NA))

Basic usage & contiguities

str(acs)

We can get neighbors based on Queen contiguities with st_contiguity().

nbs <- st_contiguity(acs)

nbs[1:5]

If needed, we can also identify the cardinalities from the neighbors list as well.

st_cardinalties(nbs)

We can get the weights from the neighbor contiguities as well. By default, st_weights() uses row standardization.

wts <- st_weights(nbs)

wts[1:5]

We can also calculate the spatial lag with the weights and neighbors.

inc_lag <- st_lag(acs$med_house_income, nbs, wts)

inc_lag[1:5]

K-Nearest Neighbor Distances

If we have point data we can also identify the k-nearest neighbors with st_knn(). For an example we can use the airbnb dataset that's imported with sfweight.

airbnb

airbnb_knn <- st_knn(airbnb)

airbnb_knn[1:5]

Other weights

Point based weights implemented based on Luc Anselin and Grant Morrison's notes.

Inverse distance band

airbnb_idw <- st_inverse_weights(airbnb$geometry, airbnb_knn)

airbnb_idw[1]

Kernel based weights

Available kernels are:

• uniform
• triangular
• epanechnikov
• quartic
• gaussian

airbnb_gauss <- st_kernel_weight(airbnb$geometry, airbnb_knn, "gaussian")

airbnb_gauss[1]

Higher order neighbors

acs %>% 
  transmute(nb = st_contiguity(geometry),
            nb_2 = st_nb_lag(nb, 2),
            nb_cumul_2 = st_nb_lag_cumul(nb, 2))

JosiahParry/sfweight documentation built on April 4, 2022, 1:52 a.m.