s2-demo.R

library(dplyr)
library(sf)
library(s2) # remotes::install_github("r-spatial/s2")

counties_sf <- readRDS(url("https://biogeo.ucdavis.edu/data/gadm3.6/Rsf/gadm36_USA_2_sf.rds")) %>% 
  select(state = NAME_1, county = NAME_2) %>% 
  # USA atlas projection
  st_transform(2163)

counties_s2 <- counties_sf %>% 
  st_transform(4326) %>% 
  as_tibble() %>% 
  mutate(geometry = as_s2_geography(sf::st_as_binary(geometry)))

# neighbouring counties using intersects
system.time(
  neighbouring_county_int_indices_sf <- 
    st_intersects(counties_sf$geometry, counties_sf$geometry)
)

system.time(
  neighbouring_county_int_indices_s2 <- 
    s2_intersects_matrix(counties_s2$geometry, counties_s2$geometry, options = s2_options(model = 2))
)


attributes(neighbouring_county_int_indices_sf) <- NULL
all.equal(neighbouring_county_int_indices_s2, neighbouring_county_int_indices_sf)

# union states together
system.time(
  states_sf <- counties_sf %>%
    group_by(state) %>% 
    summarise()
)

system.time(
  states_s2 <- counties_s2 %>% 
    group_by(state) %>% 
    summarise(geometry = s2_union_agg(geometry))
)

library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
library(sf)
#> Linking to GEOS 3.7.2, GDAL 2.4.2, PROJ 5.2.0
library(s2) # remotes::install_github("r-spatial/s2")

counties_sf <- readRDS(url("https://biogeo.ucdavis.edu/data/gadm3.6/Rsf/gadm36_USA_2_sf.rds")) %>% 
  select(state = NAME_1, county = NAME_2) %>% 
  # USA atlas projection
  st_transform(2163)

counties_s2 <- counties_sf %>% 
  st_transform(4326) %>% 
  as_tibble() %>% 
  mutate(geometry = as_s2_geography(sf::st_as_binary(geometry)))

# neighbouring counties using intersects
system.time(
  neighbouring_county_int_indices_sf <- 
    st_intersects(counties_sf$geometry, counties_sf$geometry)
)
#>    user  system elapsed 
#>   5.082   0.209   5.297

system.time(
  neighbouring_county_int_indices_s2 <- 
    s2_intersects_matrix(counties_s2$geometry, counties_s2$geometry, options = s2_options(model = 2))
)
#>    user  system elapsed 
#>   2.503   0.315   2.824


attributes(neighbouring_county_int_indices_sf) <- NULL
all.equal(neighbouring_county_int_indices_s2, neighbouring_county_int_indices_sf)
#> [1] TRUE

# union states together
system.time(
  states_sf <- counties_sf %>%
    group_by(state) %>% 
    summarise()
)
#> `summarise()` ungrouping output (override with `.groups` argument)
#>    user  system elapsed 
#> 110.710   3.300 114.332

system.time(
  states_s2 <- counties_s2 %>% 
    group_by(state) %>% 
    summarise(geometry = s2_union_agg(geometry))
)
#> `summarise()` ungrouping output (override with `.groups` argument)
#>    user  system elapsed 
#>  13.333   0.608  13.992

^{Created on 2020-07-02 by the reprex package (v0.3.0)}