Skip to content

Instantly share code, notes, and snippets.

@JosiahParry
Created September 7, 2024 14:08
Show Gist options
  • Save JosiahParry/ad3d591738b3e284c590f8b35f10bd5a to your computer and use it in GitHub Desktop.
Save JosiahParry/ad3d591738b3e284c590f8b35f10bd5a to your computer and use it in GitHub Desktop.
spatial lag categorical development
# remotes::install_github("simonpcouch/forested")
library(dplyr)
library(sfdep)
library(spdep)
trees <- forested::forested |>
sf::st_as_sf(coords = c("lon", "lat"), crs = 4326)
k <- ceiling(nrow(trees)^(1/3))
nb <- st_knn(trees, k)
wt <- st_kernel_weights(nb, trees, "gaussian", adaptive = TRUE)
listw <- nb2listw(nb, wt, "B")
glimpse(trees)
table(trees$land_type)
# sum of the weights for each class
lag_categorical_weighted_sum <- function(f, listw) {
f_levs <- levels(f)
n_cats <- length(f_levs)
n_row <- length(f)
cols <- setNames(
replicate(n_cats, double(n_row), FALSE),
cat_names
)
idx <- rep.int(1:n_row, lengths(listw$weights))
fj <- f[unlist(listw$neighbours)]
wj <- unlist(listw$weights)
res <- hardhat::weighted_table(
from = as.factor(idx),
fj,
weights = wj
)
# row standardize it
res <- as.data.frame(res / rowSums(res))
# clean column names
colnames(res) <- heck::to_snek_case(colnames(res))
# add tbl and data.frame class
class(res) <- c("tbl", "data.frame")
# return
res
}
lag_categorical_weighted_sum(trees$land_type, listw)
data.frame(f = f[nb[[1]]], wt = wt[[1]]) |>
group_by(f) |> summarise(v = sum(wt))
lag_categorical_most_freq <- function(f, listw, weighted = TRUE) {
}
cols <- list(barren = NULL, non_tree_vegetation = NULL, tree = NULL)
f <- rnorm(100)
vctrs::new_data_frame(
cols,
n = length(f),
class = "tbl"
)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment