Skip to content

Instantly share code, notes, and snippets.

@idshklein
Created August 11, 2022 14:26
Show Gist options
  • Save idshklein/749848c3114760c569c0c6cb05ee141d to your computer and use it in GitHub Desktop.
Save idshklein/749848c3114760c569c0c6cb05ee141d to your computer and use it in GitHub Desktop.
pacman::p_load(tidyverse,sf,sfdep,sfnetworks,waywiser,tidygraph)
nc = st_read(system.file("shape/nc.shp", package="sf"))
# waywiser- only polygons
waywiser::ww_build_neighbors(nc)
# sfnetworks - neighboring points connected by edges
st_intersects(nc,nc) %>%
imap_dfr(~data.frame(from=.y,to=.x)) %>%
filter(from!=to) %>%
as_tbl_graph() %>%
left_join(st_centroid(nc) %>% mutate(rn=row_number() %>% as.character()),by=c("name"="rn")) %>%
as_sfnetwork()
# both have 490 links/edges
@JosiahParry
Copy link

library(dplyr)
library(sfdep)
library(tidygraph)
library(sfnetworks)

nc <- sf::st_read(system.file("shape/nc.shp", package = "sf"))

# get neighbors from edges
sfn_nbs_from_edges <- function() {
  nbs <- igraph::neighborhood(.G())
  nb <- structure(lapply(nbs, as.integer), self.included = TRUE)
  remove_self(nb)
}

# get a variable nested as list (weights) from edges
sfn_var_from_edges <- function(var) {
  e_df <- .E()
  x <- rlang::ensym(var)
  res <- tapply(e_df[[x]], e_df[["from"]], FUN = c)
  names(res) <- NULL
  res
}

# cast nc as an sfnetwork
sfn <- nc |> 
  mutate(nb = st_contiguity(geometry)) |> 
  st_as_graph(nb)


sfn_lisa <- sfn |> 
  activate(edges) |> 
  # calculate edge length
  mutate(e_len = edge_length()) |> 
  # activate nodes to create nb and wt columns
  activate(nodes) |> 
  # create nb and wt columns
  mutate(wt = sfn_var_from_edges(e_len),
         nb = sfn_nbs_from_edges()) |> 
  # can't create data.frame columns as a tidygraph so 
  # cast to tibble
  as_tibble() |> 
  # calculate local_moran
  mutate(lisa = local_moran(BIR74, nb, wt)) 


sfn_lisa |> 
  as_tibble() |> 
  tidyr::unnest(lisa) |> 
  select(contains("ii"))


#> # A tibble: 100 × 6
#>          ii     eii      var_ii   z_ii    p_ii p_ii_sim
#>       <dbl>   <dbl>       <dbl>  <dbl>   <dbl>    <dbl>
#>  1   22693.  -143.  1176874592.  0.666 0.506      0.46 
#>  2   15694. -2462.  2121192046.  0.394 0.693      0.888
#>  3   -1429.   -94.6    7446793. -0.489 0.625      0.536
#>  4   44844. -7926.  3193846101.  0.934 0.350      0.012
#>  5   38915. -2650.  2495279893.  0.832 0.405      0.3  
#>  6   28204.   288.   785165819.  0.996 0.319      0.152
#>  7   46296. -2245.  1938775438.  1.10  0.270      0.016
#>  8   84264. -4087.  3904071505.  1.41  0.157      0.016
#>  9   30442. -2222.  4482811265.  0.488 0.626      0.728
#> 10 -112042.  1957.  1399196538. -3.05  0.00231    0.02 
#> # … with 90 more rows
#> # ℹ Use `print(n = ...)` to see more rows

Created on 2022-08-16 by the reprex package (v2.0.1)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment