rspdep

Difficulty with poly2nb() and nb2mat() in spdep package


I would appreciate some advice on making a spatial weights matrix using the spdep() package.

(1) If I have spatial data download from the rnaturalearth package, and a list of countries in another dataset, what is an efficient way to subset the data. To be specific I have:

library(rnaturalearth)
ne_sf <- ne_countries(returnclass = "sf")

And I have a list of countries so that dput(in_data_set) returns:

c("ALB", "ARE", "ARG", "ARM", "AUS", "AUT", "AZE", "BEL", "BGD", 
"BGR", "BIH", "BLR", "BLZ", "BOL", "BRA", "BWA", "CAF", "CAN", 
"CHE", "CHL", "CHN", "CMR", "COD", "COG", "COL", "CRI", "CYP", 
"DEU", "DJI", "DNK", "DZA", "ECU", "ERI", "ESP", "EST", "ETH", 
"FIN", "FJI", "FRA", "GAB", "GBR", "GEO", "GNB", "GNQ", "GRC", 
"GUY", "HRV", "HUN", "IDN", "IND", "IRL", "ISL", "ISR", "ITA", 
"JAM", "JOR", "JPN", "KAZ", "KEN", "KHM", "KWT", "LBN", "LBR", 
"LBY", "LKA", "LSO", "LTU", "LUX", "LVA", "MAR", "MDG", "MEX", 
"MLI", "MNE", "MNG", "MOZ", "MWI", "MYS", "NAM", "NIC", "NLD", 
"NOR", "NPL", "NZL", "OMN", "PAK", "PAN", "PER", "PHL", "POL", 
"PRT", "QAT", "ROU", "SAU", "SEN", "SLE", "SLV", "SRB", "SVN", 
"SWE", "TCD", "THA", "TJK", "TTO", "TUN", "UKR", "URY", "USA", 
"VNM", "ZAF", "ZMB")

So the question is how could I use this vector of three letter country codes to subset the data. one of the objects in ne_sf is iso_a3 against which the country codes could be matched.

  1. When I skip ahead and try to create a (non-subsetted) spatial matrix I get an error as follows (which I assume will not disappear when I subset the data):

    w_ne_sf <- poly2nb(ne_sf) %>% nb2mat(style = "W")

    Error in wk_handle.wk_wkb(wkb, s2_geography_writer(oriented = oriented, : Loop 7 is not valid: Edge 172 crosses edge 174


Solution

  • The answer to the first problem is easy:

    library(rnaturalearth)
    library(sf)
    library(dplyr)
    ne_sf <- ne_countries(returnclass = "sf")
    
    incl <- c("ALB", "ARE", "ARG", "ARM", "AUS", "AUT", "AZE", "BEL", "BGD", 
              "BGR", "BIH", "BLR", "BLZ", "BOL", "BRA", "BWA", "CAF", "CAN", 
              "CHE", "CHL", "CHN", "CMR", "COD", "COG", "COL", "CRI", "CYP", 
              "DEU", "DJI", "DNK", "DZA", "ECU", "ERI", "ESP", "EST", "ETH", 
              "FIN", "FJI", "FRA", "GAB", "GBR", "GEO", "GNB", "GNQ", "GRC", 
              "GUY", "HRV", "HUN", "IDN", "IND", "IRL", "ISL", "ISR", "ITA", 
              "JAM", "JOR", "JPN", "KAZ", "KEN", "KHM", "KWT", "LBN", "LBR", 
              "LBY", "LKA", "LSO", "LTU", "LUX", "LVA", "MAR", "MDG", "MEX", 
              "MLI", "MNE", "MNG", "MOZ", "MWI", "MYS", "NAM", "NIC", "NLD", 
              "NOR", "NPL", "NZL", "OMN", "PAK", "PAN", "PER", "PHL", "POL", 
              "PRT", "QAT", "ROU", "SAU", "SEN", "SLE", "SLV", "SRB", "SVN", 
              "SWE", "TCD", "THA", "TJK", "TTO", "TUN", "UKR", "URY", "USA", 
              "VNM", "ZAF", "ZMB")
    
    ne_sub <- ne_sf %>% 
      filter(iso_a3 %in% incl)
    nrow(ne_sf)
    #> [1] 177
    nrow(ne_sub)
    #> [1] 111
    

    Created on 2023-01-19 by the reprex package (v2.0.1)

    The answer to the second is more challenging. This GitHub issue suggests that the spherical geometry may not work well here and that you could rebuild it if that made sense.