Finding rook geographic neighbors with sp package

General
,
1

I'm trying to find polygon neighbors that touch at more than just one point, referred to as rook neighbors in some places. I can do this using spdep::poly2nb(., queen=FALSE) but the sp package and its dependencies are going away. Below I show that sf:st_intersects(.) is not an adequate replacement. I want to be able to exclude the red connections. Any suggestions of what else to use?

library(tidyverse)
library(sf)
#> Linking to GEOS 3.11.2, GDAL 3.6.2, PROJ 9.2.0; sf_use_s2() is TRUE
library(spdep)
#> Warning: package 'spdep' was built under R version 4.3.1
#> Loading required package: spData
#> Warning: package 'spData' was built under R version 4.3.1
#> The legacy packages maptools, rgdal, and rgeos, underpinning this package
#> will retire shortly. Please refer to R-spatial evolution reports on
#> https://r-spatial.org/r/2023/05/15/evolution4.html for details.
#> This package is now running under evolution status 0
#> To access larger datasets in this package, install the spDataLarge
#> package with: `install.packages('spDataLarge',
#> repos='https://nowosad.github.io/drat/', type='source')`

cnty_sf_in <- tigris::counties(state = "AR") # Arkansas
#> Retrieving data for the year 2021

cnty_sf <- cnty_sf_in %>%
  filter(STATEFP == "05") %>%
  select(GEOID, INTPTLAT, INTPTLON, geometry) %>%
  mutate(across(c(INTPTLAT, INTPTLON), as.numeric))

nb_queen <- poly2nb(cnty_sf, queen = TRUE)
nb_rook <- poly2nb(cnty_sf, queen = FALSE)
nb_inter <- st_intersects(cnty_sf)

nb_inter_mod <- nb_inter %>% as.list()
comp1 <- rep(NA, length(nb_inter))
comp2 <- rep(NA, length(nb_inter))

for (i in seq_along(nb_inter)) {
  nb_inter_mod[[i]] <- setdiff(nb_inter[[i]], i)
  comp1[i] <- setequal(nb_inter_mod[[i]], nb_queen[[i]])
  comp2[i] <- setequal(nb_inter_mod[[i]], nb_rook[[i]])
}

rm(i)

mean(comp1) # nb_inter same as nb_queen
#> [1] 1
mean(comp2) # not the same
#> [1] 0.9466667

cnty_centroid_sf <- cnty_sf %>%
  as_tibble() %>%
  select(-geometry) %>%
  st_as_sf(coords = c("INTPTLON", "INTPTLAT"))

cnty_centroid <- cnty_sf %>%
  as_tibble() %>%
  select(-geometry)

create_seg <- function(nb) {
  nb_sel <- function(i) {
    GEOIDA <- cnty_centroid %>%
      slice(i) %>%
      pull(GEOID)
    GEOIDB <- cnty_centroid %>%
      slice(nb[[i]]) %>%
      pull(GEOID)
    tibble(GEOIDA = GEOIDA, GEOIDB = GEOIDB)
  }

  seq_along(nb) %>%
    map(nb_sel) %>%
    list_rbind() %>%
    left_join(select(cnty_centroid, GEOIDA = GEOID, INTPTLATA = INTPTLAT, INTPTLONA = INTPTLON), by = "GEOIDA") %>%
    left_join(select(cnty_centroid, GEOIDB = GEOID, INTPTLATB = INTPTLAT, INTPTLONB = INTPTLON), by = "GEOIDB")
}

seg_queen <- create_seg(nb_queen)
seg_rook <- create_seg(nb_rook)
seg_inter <- create_seg(nb_inter_mod)

seg_diff <- mutate(seg_inter, Intersection = TRUE) %>%
  full_join(
    seg_rook %>% select(-starts_with("INT")) %>% mutate(Rook = TRUE),
    by = c("GEOIDA", "GEOIDB")
  ) %>%
  full_join(
    seg_queen %>% select(-starts_with("INT")) %>% mutate(Queen = TRUE),
    by = c("GEOIDA", "GEOIDB")
  ) %>%
  replace_na(list(Intersection = FALSE, Rook = FALSE, Queen = TRUE)) %>%
  mutate(
    MatchType = str_c(if_else(Intersection, "Int", ""),
      if_else(Queen, "Queen", ""),
      if_else(Rook, "Rook", ""),
      sep = ", "
    )
  )


ggplot() +
  geom_sf(data = cnty_sf, aes(geometry = geometry)) +
  geom_segment(data = seg_diff, aes(x = INTPTLONA, y = INTPTLATA, xend = INTPTLONB, yend = INTPTLATB, colour = MatchType), lwd = 1.5)

Created on 2023-06-22 with reprex v2.0.2

Session info 
sessioninfo::session_info()
#> ─ Session info ───────────────────────────────────────────────────────────────
#>  setting  value
#>  version  R version 4.3.0 (2023-04-21 ucrt)
#>  os       Windows 10 x64 (build 19045)
#>  system   x86_64, mingw32
#>  ui       RTerm
#>  language (EN)
#>  collate  English_United States.utf8
#>  ctype    English_United States.utf8
#>  tz       America/New_York
#>  date     2023-06-22
#>  pandoc   3.1.1 @ C:/Program Files/RStudio/resources/app/bin/quarto/bin/tools/ (via rmarkdown)
#> 
#> ─ Packages ───────────────────────────────────────────────────────────────────
#>  package     * version  date (UTC) lib source
#>  boot          1.3-28.1 2022-11-22 [1] CRAN (R 4.3.0)
#>  class         7.3-21   2023-01-23 [1] CRAN (R 4.3.0)
#>  classInt      0.4-9    2023-02-28 [1] CRAN (R 4.3.0)
#>  cli           3.6.1    2023-03-23 [1] CRAN (R 4.3.0)
#>  colorspace    2.1-0    2023-01-23 [1] CRAN (R 4.3.0)
#>  curl          5.0.0    2023-01-12 [1] CRAN (R 4.3.0)
#>  DBI           1.1.3    2022-06-18 [1] CRAN (R 4.3.0)
#>  deldir        1.0-9    2023-05-17 [1] CRAN (R 4.3.0)
#>  digest        0.6.31   2022-12-11 [1] CRAN (R 4.3.0)
#>  dplyr       * 1.1.2    2023-04-20 [1] CRAN (R 4.3.0)
#>  e1071         1.7-13   2023-02-01 [1] CRAN (R 4.3.0)
#>  evaluate      0.20     2023-01-17 [1] CRAN (R 4.3.0)
#>  fansi         1.0.4    2023-01-22 [1] CRAN (R 4.3.0)
#>  farver        2.1.1    2022-07-06 [1] CRAN (R 4.3.0)
#>  fastmap       1.1.1    2023-02-24 [1] CRAN (R 4.3.0)
#>  forcats     * 1.0.0    2023-01-29 [1] CRAN (R 4.3.0)
#>  fs            1.6.2    2023-04-25 [1] CRAN (R 4.3.0)
#>  generics      0.1.3    2022-07-05 [1] CRAN (R 4.3.0)
#>  ggplot2     * 3.4.2    2023-04-03 [1] CRAN (R 4.3.0)
#>  glue          1.6.2    2022-02-24 [1] CRAN (R 4.3.0)
#>  gtable        0.3.3    2023-03-21 [1] CRAN (R 4.3.0)
#>  highr         0.10     2022-12-22 [1] CRAN (R 4.3.0)
#>  hms           1.1.3    2023-03-21 [1] CRAN (R 4.3.0)
#>  htmltools     0.5.5    2023-03-23 [1] CRAN (R 4.3.0)
#>  httr          1.4.6    2023-05-08 [1] CRAN (R 4.3.0)
#>  KernSmooth    2.23-20  2021-05-03 [1] CRAN (R 4.3.0)
#>  knitr         1.42     2023-01-25 [1] CRAN (R 4.3.0)
#>  lattice       0.21-8   2023-04-05 [1] CRAN (R 4.3.0)
#>  lifecycle     1.0.3    2022-10-07 [1] CRAN (R 4.3.0)
#>  lubridate   * 1.9.2    2023-02-10 [1] CRAN (R 4.3.0)
#>  magrittr      2.0.3    2022-03-30 [1] CRAN (R 4.3.0)
#>  mime          0.12     2021-09-28 [1] CRAN (R 4.3.0)
#>  munsell       0.5.0    2018-06-12 [1] CRAN (R 4.3.0)
#>  pillar        1.9.0    2023-03-22 [1] CRAN (R 4.3.0)
#>  pkgconfig     2.0.3    2019-09-22 [1] CRAN (R 4.3.0)
#>  proxy         0.4-27   2022-06-09 [1] CRAN (R 4.3.0)
#>  purrr       * 1.0.1    2023-01-10 [1] CRAN (R 4.3.0)
#>  R.cache       0.16.0   2022-07-21 [1] CRAN (R 4.3.1)
#>  R.methodsS3   1.8.2    2022-06-13 [1] CRAN (R 4.3.0)
#>  R.oo          1.25.0   2022-06-12 [1] CRAN (R 4.3.0)
#>  R.utils       2.12.2   2022-11-11 [1] CRAN (R 4.3.1)
#>  R6            2.5.1    2021-08-19 [1] CRAN (R 4.3.0)
#>  rappdirs      0.3.3    2021-01-31 [1] CRAN (R 4.3.0)
#>  Rcpp          1.0.10   2023-01-22 [1] CRAN (R 4.3.0)
#>  readr       * 2.1.4    2023-02-10 [1] CRAN (R 4.3.0)
#>  reprex        2.0.2    2022-08-17 [1] CRAN (R 4.3.0)
#>  rlang         1.1.1    2023-04-28 [1] CRAN (R 4.3.0)
#>  rmarkdown     2.21     2023-03-26 [1] CRAN (R 4.3.0)
#>  rstudioapi    0.14     2022-08-22 [1] CRAN (R 4.3.0)
#>  s2            1.1.4    2023-05-17 [1] CRAN (R 4.3.0)
#>  scales        1.2.1    2022-08-20 [1] CRAN (R 4.3.0)
#>  sessioninfo   1.2.2    2021-12-06 [1] CRAN (R 4.3.0)
#>  sf          * 1.0-13   2023-05-24 [1] CRAN (R 4.3.0)
#>  sp            1.6-1    2023-05-31 [1] CRAN (R 4.3.0)
#>  spData      * 2.2.2    2023-03-01 [1] CRAN (R 4.3.1)
#>  spdep       * 1.2-8    2023-02-28 [1] CRAN (R 4.3.1)
#>  stringi       1.7.12   2023-01-11 [1] CRAN (R 4.3.0)
#>  stringr     * 1.5.0    2022-12-02 [1] CRAN (R 4.3.0)
#>  styler        1.10.1   2023-06-05 [1] CRAN (R 4.3.1)
#>  tibble      * 3.2.1    2023-03-20 [1] CRAN (R 4.3.0)
#>  tidyr       * 1.3.0    2023-01-24 [1] CRAN (R 4.3.0)
#>  tidyselect    1.2.0    2022-10-10 [1] CRAN (R 4.3.0)
#>  tidyverse   * 2.0.0    2023-02-22 [1] CRAN (R 4.3.0)
#>  tigris        2.0.3    2023-05-19 [1] CRAN (R 4.3.0)
#>  timechange    0.2.0    2023-01-11 [1] CRAN (R 4.3.0)
#>  tzdb          0.4.0    2023-05-12 [1] CRAN (R 4.3.0)
#>  units         0.8-2    2023-04-27 [1] CRAN (R 4.3.0)
#>  utf8          1.2.3    2023-01-31 [1] CRAN (R 4.3.0)
#>  uuid          1.1-0    2022-04-19 [1] CRAN (R 4.3.0)
#>  vctrs         0.6.2    2023-04-19 [1] CRAN (R 4.3.0)
#>  withr         2.5.0    2022-03-03 [1] CRAN (R 4.3.0)
#>  wk            0.7.3    2023-05-06 [1] CRAN (R 4.3.0)
#>  xfun          0.39     2023-04-20 [1] CRAN (R 4.3.0)
#>  xml2          1.3.4    2023-04-27 [1] CRAN (R 4.3.0)
#>  yaml          2.3.7    2023-01-23 [1] CRAN (R 4.3.0)
#> 
#>  [1] C:/Program Files/R/R-4.3.0/library
#> 
#> ──────────────────────────────────────────────────────────────────────────────
2

I believe there is some misunderstanding going on - not that I blame you. The spatial package ecosystem is convoluted.

The best piece of information I have is from twitter, where the package maintainer (Roger Bivand) declared a plan to maintain and actively develop spdep in the forseeable future, but based on sf (not sp) data objects.

See this thread:

My understanding is that rgdal, rgeos and maptools are going away, sp gets kind of frozen in limbo - but spdep will live on...

1 Like
3

Thanks, that's awesome! I am using sf objects

4

As I could tell from your code :slight_smile: You should be fine, but you can be excused for feeling confused - spdep on load loads sp, which emits a big fat warning, which is completely unrelated to your use case...

5

This topic was automatically closed 7 days after the last reply. New replies are no longer allowed.

If you have a query related to it or one of the replies, start a new topic and refer back with a link.