Drafting automatic spatial reg sfdep/spdep

auto-spat-reg.R

library(sfdep)
library(dplyr)
library(spdep)

guerry_nb |> 
  reframe(across(where(is.numeric), \(.x) broom::tidy(global_moran_perm(.x, nb, wt)))) |> 
  tidyr::pivot_longer(everything()) |> 
  tidyr::unnest(value)

moran_all <- function(.data, vars, nb_col = "nb", wt_col = "wt") {
  .data |> 
    reframe(
      across(all_of(vars),
      \(.x) broom::tidy(global_moran_perm(.x, .data[[nb_col]], .data[[wt_col]])))
    ) |> 
    tidyr::pivot_longer(everything()) |> 
    tidyr::unnest(value)
}

# broom tidy wrapper
gm <- function(x, nb, wt) {
  broom::tidy(global_moran_perm(x, nb, wt))
}


guerry_nb |> 
  reframe(
    across(
      where(is.numeric),  gm, wt = wt, nb = nb)
    ) |> 
  tidyr::pivot_longer(everything()) |> 
  tidyr::unnest(value)

#' Calculates the morans i for each numeric varaible in a dataset
auto_moran <- function(.data, nb_col = "nb", wt_col = "wt", ...) {
  .data |> 
    reframe(
      across(where(is.numeric),
      \(.x) broom::tidy(global_moran_perm(.x, .data[[nb_col]], .data[[wt_col]], ...)))
    ) |> 
    tidyr::pivot_longer(everything(), names_to = "variable") |> 
    tidyr::unnest(value) |> 
    dplyr::select(-c("method")) |> 
    dplyr::rename(n_obs = parameter, p_value = p.value)
}

# function to report which model is most appropriate
# need to implement branching logic from anselin & rey 2014
which_sar_model <- function(mod, nb, wt) {
  rlang::check_installed("broom")
  lm_tests <- spdep::lm.RStests(mod, recreate_listw(nb, wt), test = "all")
  all_lm_tests <- dplyr::bind_rows(lapply(lm_tests, broom::tidy), .id = "test")
  # logic to print which test to perform 
  all_ps <- setNames(all_lm_tests$p.value, all_lm_tests$test)
  # insert the logic for returning which model 
  chosen_model <- names(all_ps[2])
  cli::cli_alert_info("Most appropriate model is {chosen_model}")
  attr(all_lm_tests, "chosen_model") <- chosen_model
  all_lm_tests
}

mod <- lm(crime_prop ~ literacy, data = guerry)

tmp <- which_sar_model(mod, guerry_nb$nb, guerry_nb$wt)

attr(tmp, "chosen_model")


auto_moran(guerry_nb, alternative = "greater")

guerry_nb |> 
  moran_all(vars = c("literacy", "infants"))


moran_all()