From 4b78dfea729c836550fda9faf55f504983d29a3d Mon Sep 17 00:00:00 2001 From: Roger Bivand Date: Tue, 10 Sep 2024 10:38:11 +0200 Subject: [PATCH] update docs --- docs/articles/CO69.html | 6 +- docs/articles/nb_sf.html | 28 ++--- docs/articles/subgraphs.html | 175 +++++++++++++++++++------------- docs/pkgdown.yml | 2 +- docs/reference/bhicv.html | 16 ++- docs/reference/compon.html | 4 +- docs/reference/dnearneigh.html | 10 +- docs/reference/knearneigh.html | 2 +- docs/reference/localC.html | 4 +- docs/reference/mstree.html | 42 +++++++- docs/reference/nb2lines.html | 8 +- docs/reference/poly2nb.html | 4 +- docs/reference/read.gwt2nb.html | 17 +++- docs/reference/skater.html | 65 ++++++++++-- docs/search.json | 2 +- 15 files changed, 258 insertions(+), 127 deletions(-) diff --git a/docs/articles/CO69.html b/docs/articles/CO69.html index 8e7ee163..f89106aa 100644 --- a/docs/articles/CO69.html +++ b/docs/articles/CO69.html @@ -548,7 +548,7 @@

Measures of spatial autocorrelation Prop_stdR <- lapply(vars, function(x) moran.test(eire_ge1[[x]], listw=lw_std, randomisation=TRUE)) }) 
##    user  system elapsed 
-##   0.126   0.000   0.128
+## 0.12 0.00 0.12 
 res <- sapply(c("MoranN", "MoranR", "GearyN", "GearyR", "Prop_unstdN", "Prop_unstdR", "Prop_stdN", "Prop_stdR"), function(x) sapply(get(x), "[[", "statistic"))
 rownames(res) <- vars
@@ -780,7 +780,7 @@ 
Simulating measures of s
 Prop_stdSad  <- lapply(lm_objs, function(x) lm.morantest.sad(x, listw=lw_std))
 })
##    user  system elapsed 
-##   0.071   0.000   0.072
+## 0.066 0.000 0.066 
 res <- sapply(c("MoranSad", "Prop_unstdSad", "Prop_stdSad"), function(x) sapply(get(x), "[[", "statistic"))
 rownames(res) <- vars
@@ -843,7 +843,7 @@

Simulating measures of s Prop_stdEx <- lapply(lm_objs, function(x) lm.morantest.exact(x, listw=lw_std)) }) 
##    user  system elapsed 
-##   0.092   0.000   0.094
+## 0.082 0.000 0.082 
 res <- sapply(c("MoranEx", "Prop_unstdEx", "Prop_stdEx"), function(x) sapply(get(x), "[[", "statistic"))
 rownames(res) <- vars
diff --git a/docs/articles/nb_sf.html b/docs/articles/nb_sf.html index 8807b37a..7eda3152 100644 --- a/docs/articles/nb_sf.html +++ b/docs/articles/nb_sf.html @@ -230,7 +230,7 @@

Contiguity neighbours for pol eps <- sqrt(.Machine$double.eps) system.time(for(i in 1:reps) NY8_sf_1_nb <- poly2nb(NY8_sf, queen=TRUE, snap=eps))/reps 
##    user  system elapsed 
-##  0.1798  0.0100  0.1911
+## 0.1853 0.0094 0.1960

Using spatial indices to check intersection of polygons is much faster than the legacy method in poly2nb. From spdep 1.1-7, use is made of GEOS through sf to find candidate @@ -406,12 +406,12 @@

K-nearest neighbours
 system.time(for (i in 1:reps) suppressWarnings(NY88_nb_sf <- knn2nb(knearneigh(NY8_ct_sf, k=1))))/reps
##    user  system elapsed 
-##  0.0232  0.0017  0.0252
+## 0.0271 0.0013 0.0286

Legacy code may be used omitting the kd-tree:

 system.time(for (i in 1:reps) suppressWarnings(NY89_nb_sf <- knn2nb(knearneigh(NY8_ct_sf, k=1, use_kd_tree=FALSE))))/reps
##    user  system elapsed 
-##  0.0237  0.0015  0.0253
+## 0.0263 0.0019 0.0284

Distance neighbours @@ -430,7 +430,7 @@

Distance neighbours
 system.time(for (i in 1:reps) suppressWarnings(NY810_nb <- dnearneigh(NY8_ct_sf, d1=0, d2=0.75*max_1nn)))/reps

##    user  system elapsed 
-##  0.0616  0.0015  0.0636
+## 0.0619 0.0013 0.0634

By default, the function uses dbscan::frNN() to build a kd-tree in 2D or 3D which is then used to find distance neighbours. For small n, the argument use_kd_tree=FALSE may speed up @@ -440,7 +440,7 @@

Distance neighbours
 system.time(for (i in 1:reps) suppressWarnings(NY811_nb <- dnearneigh(NY8_ct_sf, d1=0, d2=0.75*max_1nn, use_kd_tree=FALSE)))/reps
##    user  system elapsed 
-##  0.0464  0.0013  0.0480
+## 0.0469 0.0008 0.0479
@@ -470,7 +470,7 @@

K-nearest neighbourssf_use_s2(TRUE) system.time(for (i in 1:reps) pts_ll1_nb <- knn2nb(knearneigh(pts_ll, k=6)))/reps

##    user  system elapsed 
-##  0.0335  0.0001  0.0339
+## 0.0331 0.0001 0.0332

For this smaller data set, the legacy approach without spatial indexing is adequate, but slows down as the number of observations increases:

@@ -480,7 +480,7 @@

K-nearest neighbours
 system.time(for (i in 1:reps) pts_ll2_nb <- knn2nb(knearneigh(pts_ll, k=6)))/reps
##    user  system elapsed 
-##  0.0286  0.0000  0.0287
+## 0.0287 0.0001 0.0289

The WGS84 ellipsoid Great Circle distances differ a very little from the s2 spherical distances, yielding output that here diverges for two tract centroids:

@@ -535,7 +535,7 @@

Distance neighbours=0.75*max_1nn_ll)))/(reps/5) } 
##    user  system elapsed 
-##   0.068   0.000   0.069
+## 0.068 0.000 0.068

Alternatively, spherical distances can be used with dwithin=FALSE and s2::s2_closest_edges(); although running in similar time, s2::s2_closest_edges() @@ -544,7 +544,7 @@

Distance neighbours
 system.time(for (i in 1:(reps/5)) suppressWarnings(pts_ll5_nb <- dnearneigh(pts_ll, d1=0, d2=0.75*max_1nn_ll, dwithin=FALSE)))/(reps/5)
##    user  system elapsed 
-##   0.055   0.000   0.055
+## 0.054 0.000 0.054 
 if (packageVersion("s2") > "1.0.7") all.equal(pts_ll3_nb, pts_ll5_nb, check.attributes=FALSE)
## [1] TRUE
@@ -557,7 +557,7 @@

Distance neighbours=0.75*max_1nn_ll, dwithin=FALSE)))/(reps/5) } 
##    user  system elapsed 
-##   0.046   0.000   0.046
+## 0.048 0.000 0.048

Using s2::s2_dwithin_matrix() requires a second pass, one for the lower bound, another for the upper bound, and a set difference operation to find neighbours in the distance band:

@@ -567,7 +567,7 @@

Distance neighbours=0.75*max_1nn_ll)))/(reps/5) } 
##    user  system elapsed 
-##   0.089   0.000   0.090
+## 0.0895 0.0005 0.0905 
 if (packageVersion("s2") > "1.0.7") all.equal(pts_ll3a_nb, pts_ll5a_nb, check.attributes=FALSE)
## [1] TRUE
@@ -576,7 +576,7 @@

Distance neighbours
 system.time(for (i in 1:reps) suppressWarnings(pts_ll6_nb <- dnearneigh(pts_ll, d1=0, d2=0.75*max_1nn_ll, use_s2=FALSE)))/reps
##    user  system elapsed 
-##  0.0382  0.0001  0.0386
+## 0.0376 0.0000 0.0377

Minor differences may occur between the legacy ellipsoid and s2 spherical approaches:

@@ -606,7 +606,7 @@ 
Distance neighbours
 system.time(for (i in 1:reps) suppressWarnings(pts_ll6a_nb <- dnearneigh(pts_ll, d1=5, d2=0.75*max_1nn_ll, use_s2=FALSE)))/reps
##    user  system elapsed 
-##  0.0299  0.0000  0.0303
+## 0.0314 0.0001 0.0318 
 if (packageVersion("s2") > "1.0.7") all.equal(pts_ll5a_nb, pts_ll6a_nb, check.attributes=FALSE)
##  [1] "Component 20: Numeric: lengths (6, 5) differ"       
@@ -658,7 +658,7 @@ 
Contiguity neighbours for spherical polygon supportsf_use_s2(TRUE)
 system.time(for (i in 1:reps) NY8_sf_1_nb_ll <- poly2nb(NY8_sf_ll, queen=TRUE, snap=eps))/reps
##    user  system elapsed 
-##  0.1613  0.0017  0.1643
+## 0.1660 0.0016 0.1682 
 all.equal(NY8_sf_1_nb, NY8_sf_1_nb_ll, check.attributes=FALSE)
## [1] TRUE
diff --git a/docs/articles/subgraphs.html b/docs/articles/subgraphs.html index 14b078b0..91d90574 100644 --- a/docs/articles/subgraphs.html +++ b/docs/articles/subgraphs.html @@ -156,26 +156,37 @@

No-neighbour observationsFreni-Sterrantino, Ventrucci, and Rue (2018), where Capraia and Giglio Isles are singleton nodes. Here we take Westminster constituencies for Wales used in the July 2024 UK general -election.

+election. If GDAL is at least version 3.7.0, the driver supports +compressed GeoPackage files, if not they must be decompressed first.

 +(GDAL37 <- as.numeric_version(unname(sf_extSoftVersion()["GDAL"])) >= "3.7.0") +
## [1] TRUE

The boundaries are taken from the Ordnance Survey Boundary-Line site, https://osdatahub.os.uk/downloads/open/BoundaryLine, choosing the 2024 Westminster constituencies (https://www.os.uk/opendata/licence), simplified using a tolerance of 50m to reduce object size, and merged with selected voting outcomes for constituencies in Great Britain https://electionresults.parliament.uk/countries/1, (https://www.nationalarchives.gov.uk/doc/open-government-licence/version/3/). Here, the subset for Wales is useful as we will see:

-
-w50m <- st_read(system.file("etc/shapes/GB_2024_Wales_50m.gpkg.zip", package="spdep"))
+
+file <- "etc/shapes/GB_2024_Wales_50m.gpkg.zip"
+zipfile <- system.file(file, package="spdep")
+if (GDAL37) {
+    w50m <- st_read(zipfile)
+} else {
+    td <- tempdir()
+    bn <- sub(".zip", "", basename(file), fixed=TRUE)
+    target <- unzip(zipfile, files=bn, exdir=td)
+    w50m <- st_read(target)
+}
## Reading layer `GB_2024_Wales_50m' from data source 
-##   `/tmp/Rtmprq0AG5/temp_libpath6568f436245/spdep/etc/shapes/GB_2024_Wales_50m.gpkg.zip' 
+##   `/tmp/RtmpPpr0yz/temp_libpathafa452560f839/spdep/etc/shapes/GB_2024_Wales_50m.gpkg.zip' 
 ##   using driver `GPKG'
 ## Simple feature collection with 32 features and 19 fields
 ## Geometry type: MULTIPOLYGON
 ## Dimension:     XY
 ## Bounding box:  xmin: 146597.1 ymin: 164536.5 xmax: 355287 ymax: 395993.5
 ## Projected CRS: OSGB36 / British National Grid
-
+
 (w50m |> poly2nb(row.names=as.character(w50m$Constituency)) -> nb_W_50m)

 
## Warning in poly2nb(w50m, row.names = as.character(w50m$Constituency)): some observations have no neighbours;
 ## if this seems unexpected, try increasing the snap argument.

@@ -191,7 +202,7 @@ 
No-neighbour observations## 2 disjoint connected subgraphs

The two subgraphs are the singleton Ynys Môn and all the other 31 constituencies:

-
+
 attr(nb_W_50m, "ncomp")$comp.id |>table() |> table()

 
## 
 ##  1 31 
@@ -201,7 +212,7 @@ 
No-neighbour observationscard, filling the polygon. The right map shows the location
 of the island, known in English as Anglesey, north-west of the Welsh
 mainland, and with no neighbour links:

-
+
 ynys_mon <- w50m$Constituency == "Ynys Môn"
 pts <- st_point_on_surface(st_geometry(w50m))
 opar <- par(mfrow=c(1, 2))
@@ -211,13 +222,13 @@ 
No-neighbour observationsplot(st_geometry(w50m), border="grey75")
 plot(nb_W_50m, pts, add=TRUE)

 

-
+
 par(opar)

 
From the maps, we can see that the island is close to two
 constituencies across the Afon Menai (Menai Strait in English), the
 three simplified polygons being less than 280m apart, measured between
 polygon boundaries:

-
+
 dym <- c(st_distance(w50m[ynys_mon,], w50m))
 sort(dym)[1:12]

 
## Units: [m]
@@ -225,14 +236,14 @@ 
No-neighbour observations##  [7]  58146.4320  65550.2491  67696.3323  93741.9873 113007.3659 137858.1826

 
Using a snap distance of 280m, we can join the island to
 its two obvious proximate neighbours:

-
+
 (nb_W_50m_snap <- poly2nb(w50m, row.names=as.character(w50m$Constituency), snap=280))

 
## Neighbour list object:
 ## Number of regions: 32 
 ## Number of nonzero links: 140 
 ## Percentage nonzero weights: 13.67188 
 ## Average number of links: 4.375

-
+
 plot(st_geometry(w50m), border="grey75")
 plot(nb_W_50m_snap, pts, add=TRUE)

 

@@ -240,7 +251,7 @@ 
No-neighbour observationssqrt(.Machine$double.eps) 1.4901161^{-8} in all cases)
 helps. The symmetric links added are to:

-
+
 attr(nb_W_50m_snap, "region.id")[nb_W_50m_snap[[which(ynys_mon)]]]

 
## [1] "Bangor Aberconwy"   "Dwyfor Meirionnydd"

 
This is not always going to be the case, but here the strait is
@@ -254,12 +265,12 @@ 
No-neighbour observationsWe can also use the distances to pick out those neighbour candidates
 that meet our criterion of 280m, taking care not to lose the ordering
 needed to identify the correct observations:

-
+
 (meet_criterion <- sum(dym <= units::set_units(280, "m")))

 
## [1] 3

 
These candidates are the island itself, and the two neighbours across
 the Menai Strait:

-
+
 (cands <- attr(nb_W_50m, "region.id")[order(dym)[1:meet_criterion]])

 
## [1] "Ynys Môn"           "Bangor Aberconwy"   "Dwyfor Meirionnydd"

 
The addlinks1 function can be used to add both the links
@@ -267,14 +278,14 @@ 
No-neighbour observations
-
+
 (nb_W_50m_add <- addlinks1(nb_W_50m, from = cands[1], to = cands[2:meet_criterion]))

 
## Neighbour list object:
 ## Number of regions: 32 
 ## Number of nonzero links: 140 
 ## Percentage nonzero weights: 13.67188 
 ## Average number of links: 4.375

-
+
 all.equal(nb_W_50m_add, nb_W_50m_snap, check.attributes=FALSE)

 
## [1] TRUE

 
Since these constituency observations have areal support, it is not
@@ -283,11 +294,11 @@ 
No-neighbour observations
-
+
 k2 <- knn2nb(knearneigh(pts, k=2), row.names=as.character(w50m$Constituency), sym=TRUE)

 
## Warning in knn2nb(knearneigh(pts, k = 2), row.names =
 ## as.character(w50m$Constituency), : neighbour object has 2 sub-graphs

-
+
 attr(k2, "region.id")[k2[[which(ynys_mon)]]]

 
## [1] "Bangor Aberconwy" "Clwyd North"

 
Here, Clwyd North, east of Bangor Aberconwy, is given as a neighbour
@@ -302,7 +313,7 @@ 
Subgraphs
+
 (k6 <- knn2nb(knearneigh(pts, k=6), row.names=as.character(w50m$Constituency), sym=TRUE))

 
## Warning in knn2nb(knearneigh(pts, k = 6), row.names =
 ## as.character(w50m$Constituency), : neighbour object has 2 sub-graphs

@@ -312,13 +323,13 @@ 
Subgraphs## Percentage nonzero weights: 23.24219 
 ## Average number of links: 7.4375 
 ## 2 disjoint connected subgraphs

-
+
 plot(st_geometry(w50m), border="grey75")
 plot(k6, pts, add=TRUE)

 

 We can show the block structure by displaying the binary spatial weights
 matrix:

-
+
 o <- order(attr(k6, "ncomp")$comp.id)
 image(t(nb2mat(k6, style="B")[o, rev(o)]), axes=FALSE, asp=1)

 

@@ -334,7 +345,7 @@ 
Subgraphsigraph::components can handle directed graphs without such
 conversion (see https://github.com/r-spatial/spdep/issues/160 for
 details).

-
+
 (k6a <- knn2nb(knearneigh(pts, k=6), row.names=as.character(w50m$Constituency)))

 
## Warning in knn2nb(knearneigh(pts, k = 6), row.names =
 ## as.character(w50m$Constituency)): neighbour object has 2 sub-graphs

@@ -348,17 +359,26 @@ 
Subgraphs
-sc50m <- st_read(system.file("etc/shapes/GB_2024_southcoast_50m.gpkg.zip", package="spdep"))

+
+file <- "etc/shapes/GB_2024_southcoast_50m.gpkg.zip"
+zipfile <- system.file(file, package="spdep")
+if (GDAL37) {
+    sc50m <- st_read(zipfile)
+} else {
+    td <- tempdir()
+    bn <- sub(".zip", "", basename(file), fixed=TRUE)
+    target <- unzip(zipfile, files=bn, exdir=td)
+    sc50m <- st_read(target)
+}

 
## Reading layer `GB_2024_southcoast_50m' from data source 
-##   `/tmp/Rtmprq0AG5/temp_libpath6568f436245/spdep/etc/shapes/GB_2024_southcoast_50m.gpkg.zip' 
+##   `/tmp/RtmpPpr0yz/temp_libpathafa452560f839/spdep/etc/shapes/GB_2024_southcoast_50m.gpkg.zip' 
 ##   using driver `GPKG'
 ## Simple feature collection with 119 features and 19 fields
 ## Geometry type: MULTIPOLYGON
 ## Dimension:     XY
 ## Bounding box:  xmin: 82643.12 ymin: 5342.9 xmax: 640301.6 ymax: 187226.2
 ## Projected CRS: OSGB36 / British National Grid

-
+
 (nb_sc_50m <- poly2nb(sc50m, row.names=as.character(sc50m$Constituency)))

 
## Warning in poly2nb(sc50m, row.names = as.character(sc50m$Constituency)): neighbour object has 2 sub-graphs;
 ## if this sub-graph count seems unexpected, try increasing the snap argument.

@@ -370,17 +390,17 @@ 
Subgraphs## 2 disjoint connected subgraphs

 
The second subgraph only has two members, who are each others’ only
 neighbours, known as a cyclical component.

-
+
 nc <- attr(nb_sc_50m, "ncomp")$comp.id
 table(nc)

 
## nc
 ##   1   2 
 ## 117   2

 
Both constituencies are on the Isle of Wight:

-
+
 (sub2 <- attr(nb_sc_50m, "region.id")[nc == 2L])

 
## [1] "Isle of Wight East" "Isle of Wight West"

-
+
 pts <- st_point_on_surface(st_geometry(sc50m))
 plot(st_geometry(sc50m), border="grey75")
 plot(st_geometry(sc50m)[nc == 2L], border="orange", lwd=2, add=TRUE)
@@ -391,31 +411,31 @@ 
SubgraphsBivand, Hauke, and Kossowski
 (2013). With row-standardised weights, the eigenvalues of this
 component are:

-
+
 1/range(eigen(cbind(c(0, 1), c(1, 0)))$values)

 
## [1] -1  1

-
+
 1/range(eigen(nb2mat(subset(nb_sc_50m, nc == 2L), style="W"))$values)

 
## [1] -1  1

 
This “takes over” the lower domain boundary, which for the whole data
 set is now the same:

-
+
 1/range(eigen(nb2mat(nb_sc_50m, style="W"))$values)

 
## [1] -1  1

 
compared to the lower domain boundary for the remainder of the study
 area:

-
+
 1/range(eigen(nb2mat(subset(nb_sc_50m, nc == 1L), style="W"))$values)

 
## [1] -1.094637  1.000000

 
This subgraph may be added to the remainder as shown above:

-
+
 iowe <- match(sub2[1], attr(nb_sc_50m, "region.id"))
 diowe <- c(st_distance(sc50m[iowe,], sc50m))
 sort(diowe)[1:12]

 
## Units: [m]
 ##  [1]     0.000     0.000  1886.833  3509.366  6693.575  6943.672  7678.999
 ##  [8]  8576.454 10579.530 12163.332 16875.920 17161.786

-
+
 ioww <- match(sub2[2], attr(nb_sc_50m, "region.id"))
 dioww <- c(st_distance(sc50m[ioww,], sc50m))
 sort(dioww)[1:12]

@@ -425,17 +445,17 @@ 
Subgraphs
+
 (meet_criterion <- sum(diowe <= units::set_units(5000, "m")))

 
## [1] 4

 
Obviously the contiguous neighbour is among them with zero distance,
 and needs to be dropped, although addlinks1 would drop the
 duplicate:

-
+
 (cands <- attr(nb_sc_50m, "region.id")[order(diowe)[1:meet_criterion]])

 
## [1] "Isle of Wight East" "Isle of Wight West" "Portsmouth South"  
 ## [4] "Gosport"

-
+
 (nb_sc_50m_iowe <- addlinks1(nb_sc_50m, from = cands[1], to = cands[3:meet_criterion]))

 
## Neighbour list object:
 ## Number of regions: 119 
@@ -445,25 +465,25 @@ 
Subgraphs
+
 (meet_criterion <- sum(dioww <= units::set_units(5000, "m")))

 
## [1] 4

 
It, does, but we need to beware of the sorting order of the zero
 self-distance and contiguous neighbour distance, where from
 is now in the second position:

-
+
 (cands <- attr(nb_sc_50m, "region.id")[order(dioww)[1:meet_criterion]])

 
## [1] "Isle of Wight East" "Isle of Wight West" "New Forest West"   
 ## [4] "New Forest East"

 
This then yields links to the north-west:

-
+
 (nb_sc_50m_iow <- addlinks1(nb_sc_50m_iowe, from = cands[2], to = cands[3:meet_criterion]))

 
## Neighbour list object:
 ## Number of regions: 119 
 ## Number of nonzero links: 538 
 ## Percentage nonzero weights: 3.799167 
 ## Average number of links: 4.521008

-
+
 pts <- st_point_on_surface(st_geometry(sc50m))
 plot(st_geometry(sc50m), border="grey75")
 plot(st_geometry(sc50m)[nc == 2L], border="orange", lwd=2, add=TRUE)
@@ -480,13 +500,13 @@ 
Per-session control of functi
 argument to many methods and functions was NULL. If the
 value was NULL, zero.policy would be set from
 get.ZeroPolicyOption:

-
+
 
## [1] FALSE

 
On loading spdep, the internal option is set to
 FALSE, so functions and methods using
 zero.policy need to choose how to handle islands:

-
+
 try(nb2listw(nb_W_50m))

 
## Error in nb2listw(nb_W_50m) : 
 ##   Empty neighbour sets found (zero.policy: FALSE)

@@ -495,12 +515,12 @@ 
Per-session control of functi
 the user wishes to override the default,
 set.ZeroPolicyOption may be used to set a different
 per-session default:

-
 
+
 
## [1] TRUE

-
+
 (lw <- nb2listw(nb_W_50m))

 
## Characteristics of weights list object:
 ## Neighbour list object:
@@ -516,10 +536,10 @@ 
Per-session control of functi
 ## Weights constants summary:
 ##    n  nn S0       S1      S2
 ## W 31 961 31 15.36355 129.051

-
+
 attr(lw, "zero.policy")

 
## [1] TRUE

-
+
 
When a listw object is created with
 zero.policy set to TRUE, this choice is added
@@ -532,13 +552,13 @@ 
Per-session control of functi
 
Other internal options have been introduced to suppress no-neighbour
 and subgraph warnings when creating nb objects. The default
 values are as follows:

-
+
 
## [1] TRUE

-
+
 
## [1] TRUE

-
+
 
## [1] 100000

 
get.NoNeighbourOption controls the issuing of warnings
@@ -559,7 +579,7 @@ 
Per-session control of functi
 creation of training and test data sets. Here the Welsh constituency
 boundaries will be used to show the behaviour of the internal
 options:

-
+
 set.NoNeighbourOption(FALSE)
 (w50m |> poly2nb(row.names=as.character(w50m$Constituency)) -> nb_W_50mz)

 
## Warning in poly2nb(w50m, row.names = as.character(w50m$Constituency)): neighbour object has 2 sub-graphs;
@@ -573,7 +593,7 @@ 
Per-session control of functi
 ## Ynys Môn
 ## 2 disjoint connected subgraphs

 
Turning both off removes the warnings:

-
+
 set.SubgraphOption(FALSE)
 (w50m |> poly2nb(row.names=as.character(w50m$Constituency)) -> nb_W_50my)

 
## Neighbour list object:
@@ -585,12 +605,12 @@ 
Per-session control of functi
 ## Ynys Môn

 
When get.SubgraphOption is FALSE, the attribute
 containing the output of n.comp.nb is not added:

-
+
 str(attr(nb_W_50my, "ncomp"))

 
##  NULL

 
The reduction of the ceiling to below node count 32 plus edge count
 136 also supresses the calculation of graph components:

-
+
 set.SubgraphOption(TRUE)
 set.SubgraphCeiling(100L)
 (w50m |> poly2nb(row.names=as.character(w50m$Constituency)) -> nb_W_50mx)

@@ -601,11 +621,11 @@ 
Per-session control of functi
 ## Average number of links: 4.25 
 ## 1 region with no links:
 ## Ynys Môn

-
+
 str(attr(nb_W_50mx, "ncomp"))

 
##  NULL

 
Restoring the remaining default values:

-
@@ -622,10 +642,19 @@ 
Unintentional disconnected graphs
-
-tokyo <- st_read(system.file("etc/shapes/tokyo.gpkg.zip", package="spdep"))

+
+file <- "etc/shapes/tokyo.gpkg.zip"
+zipfile <- system.file(file, package="spdep")
+if (GDAL37) {
+    tokyo <- st_read(zipfile)
+} else {
+    td <- tempdir()
+    bn <- sub(".zip", "", basename(file), fixed=TRUE)
+    target <- unzip(zipfile, files=bn, exdir=td)
+    tokyo <- st_read(target)
+}

 
## Reading layer `tokyo' from data source 
-##   `/tmp/Rtmprq0AG5/temp_libpath6568f436245/spdep/etc/shapes/tokyo.gpkg.zip' 
+##   `/tmp/RtmpPpr0yz/temp_libpathafa452560f839/spdep/etc/shapes/tokyo.gpkg.zip' 
 ##   using driver `GPKG'
 ## Simple feature collection with 262 features and 3 fields
 ## Geometry type: MULTIPOLYGON
@@ -633,15 +662,15 @@ 
Unintentional disconnected graphs## Bounding box:  xmin: 266206.6 ymin: -90932.11 xmax: 411400.3 ymax: 37142.75
 ## Projected CRS: Tokyo / Japan Plane Rectangular CS VI

 
After correcting invalid polygons:

-
+
 
## [1] TRUE

-
+
 tokyo <- st_make_valid(tokyo)

 
applying poly2nb with the legacy default snap value
 produced numerous singleton observations as well as many
 multiple-observation subgraphs:

-
+
 (nb_t0 <- poly2nb(tokyo, snap=sqrt(.Machine$double.eps)))

 
## Warning in poly2nb(tokyo, snap = sqrt(.Machine$double.eps)): some observations have no neighbours;
 ## if this seems unexpected, try increasing the snap argument.

@@ -658,38 +687,38 @@ 
Unintentional disconnected graphsThe legacy default snap value when the coordinates are
 measured in metres was 15 nanometres, which effectively assumed that the
 coordinates making up polygon boundaries were identical:

-
+
 units::set_units(units::set_units(attr(nb_t0, "snap"), "m"), "nm")

 
## 14.90116 [nm]

 
Stepping out a little to 2mm, the lack of contact ceased to be a
 problem.

-
+
 (nb_t1 <- poly2nb(tokyo, snap=0.002))

 
## Neighbour list object:
 ## Number of regions: 262 
 ## Number of nonzero links: 1390 
 ## Percentage nonzero weights: 2.02494 
 ## Average number of links: 5.305344

-
+
 units::set_units(units::set_units(attr(nb_t1, "snap"), "m"), "mm")

 
## 2 [mm]

 
On that basis, the default was changed from spdep 1.3-6
 to 10mm for projected polygons, and the snap value used was returned as
 an attribute of the neighbour object:

-
+
 (nb_t2 <- poly2nb(tokyo))

 
## Neighbour list object:
 ## Number of regions: 262 
 ## Number of nonzero links: 1390 
 ## Percentage nonzero weights: 2.02494 
 ## Average number of links: 5.305344

-
+
 units::set_units(units::set_units(attr(nb_t2, "snap"), "m"), "mm")

 
## 10 [mm]

 
Where the polygons are represented by geographical (spherical)
 coordinates, the new default from spdep 1.3-6 is set to a
 value mimicking 10mm:

-
+
 (nb_t3 <- poly2nb(st_transform(tokyo, "OGC:CRS84")))

 
## Neighbour list object:
 ## Number of regions: 262 
@@ -698,12 +727,12 @@ 
Unintentional disconnected graphs## Average number of links: 5.099237

 
The default snap value used in poly2nb when
 the polygons are expressed in decimal degrees is:

-
+
 attr(nb_t3, "snap")

 
## [1] 9e-08

 
This was set based on the apparent “size” of 10mm in decimal
 degrees:

-
+
 (180 * 0.01) / (pi * 6378137)

 
## [1] 8.983153e-08

 

diff --git a/docs/pkgdown.yml b/docs/pkgdown.yml
index 11094647..a91fa60a 100644
--- a/docs/pkgdown.yml
+++ b/docs/pkgdown.yml
@@ -7,7 +7,7 @@ articles:
   nb: nb.html
   sids: sids.html
   subgraphs: subgraphs.html
-last_built: 2024-09-09T14:11Z
+last_built: 2024-09-10T08:34Z
 urls:
   reference: https://r-spatial.github.io/spdep/reference
   article: https://r-spatial.github.io/spdep/articles
diff --git a/docs/reference/bhicv.html b/docs/reference/bhicv.html
index 7b1c898b..e0475f4f 100644
--- a/docs/reference/bhicv.html
+++ b/docs/reference/bhicv.html
@@ -77,12 +77,20 @@ 
Format<
 
     

     
Examples

-    
if (as.numeric_version(unname(sf_extSoftVersion()["GDAL"])) >= "3.7.0") {
-bh <- st_read(system.file("etc/shapes/bhicv.gpkg.zip",
-      package="spdep")[1])
+    
(GDAL37 <- as.numeric_version(unname(sf_extSoftVersion()["GDAL"])) >= "3.7.0")
+#> [1] TRUE
+file <- "etc/shapes/bhicv.gpkg.zip"
+zipfile <- system.file(file, package="spdep")
+if (GDAL37) {
+    bh <- st_read(zipfile)
+} else {
+    td <- tempdir()
+    bn <- sub(".zip", "", basename(file), fixed=TRUE)
+    target <- unzip(zipfile, files=bn, exdir=td)
+    bh <- st_read(target)
 }
 #> Reading layer `bhicv' from data source 
-#>   `/tmp/Rtmprq0AG5/temp_libpath6568f436245/spdep/etc/shapes/bhicv.gpkg.zip' 
+#>   `/tmp/RtmpPpr0yz/temp_libpathafa452560f839/spdep/etc/shapes/bhicv.gpkg.zip' 
 #>   using driver `GPKG'
 #> Simple feature collection with 98 features and 8 fields
 #> Geometry type: POLYGON
diff --git a/docs/reference/compon.html b/docs/reference/compon.html
index 0191cf42..28873812 100644
--- a/docs/reference/compon.html
+++ b/docs/reference/compon.html
@@ -189,12 +189,12 @@ 
Examplessystem.time(udir <- n.comp.nb(make.sym.nb(k6)))
 }
 #>    user  system elapsed 
-#>   2.196   0.000   2.223 
+#>   2.104   0.000   2.111 
 if (run) {
 system.time(dir <- n.comp.nb(k6))
 }
 #>    user  system elapsed 
-#>   0.734   0.003   0.745 
+#>   0.701   0.003   0.707 
 if (run) {
 udir$nc
 }
diff --git a/docs/reference/dnearneigh.html b/docs/reference/dnearneigh.html
index 43b7c25c..c226d66f 100644
--- a/docs/reference/dnearneigh.html
+++ b/docs/reference/dnearneigh.html
@@ -222,14 +222,14 @@ 
Examples#> Warning: neighbour object has 11 sub-graphs
 system.time(o <- nbdists(gck1b, xy1))
 #>    user  system elapsed 
-#>   0.006   0.000   0.006 
+#>   0.007   0.000   0.007 
 (all.linked <- max(unlist(o)))
 #> [1] 522.4464
 # use s2 brute-force dwithin_matrix approach for s2 <= 1.0.7
 system.time(gc.nb.dwithin <- dnearneigh(xy1, 0, all.linked, use_s2=TRUE, dwithin=TRUE))
 #> Warning: neighbour object has 3 sub-graphs
 #>    user  system elapsed 
-#>   0.013   0.000   0.013 
+#>   0.012   0.000   0.011 
 summary(gc.nb, xy1, scale=0.5)
 #> Neighbour list object:
 #> Number of regions: 48 
@@ -251,13 +251,13 @@ 
Examples}
 #> Warning: neighbour object has 3 sub-graphs
 #>    user  system elapsed 
-#>   0.011   0.000   0.010 
+#>    0.01    0.00    0.01 
 if (packageVersion("s2") > "1.0.7") {
 system.time(gc.nb.dwithin <- dnearneigh(xy1, 0, all.linked, use_s2=TRUE, dwithin=TRUE))
 }
 #> Warning: neighbour object has 3 sub-graphs
 #>    user  system elapsed 
-#>   0.013   0.000   0.012 
+#>   0.012   0.000   0.012 
 if (packageVersion("s2") > "1.0.7") {
 summary(gc.nb.dwithin, xy1, scale=0.5)
 }
@@ -300,7 +300,7 @@ 
Examplessystem.time(gc.nb.legacy <- dnearneigh(xy1, 0, all.linked, use_s2=FALSE))
 #> Warning: neighbour object has 3 sub-graphs
 #>    user  system elapsed 
-#>   0.007   0.000   0.007 
+#>   0.007   0.000   0.006 
 summary(gc.nb, xy1, scale=0.5)
 #> Neighbour list object:
 #> Number of regions: 48 
diff --git a/docs/reference/knearneigh.html b/docs/reference/knearneigh.html
index fae88d25..03b8c654 100644
--- a/docs/reference/knearneigh.html
+++ b/docs/reference/knearneigh.html
@@ -166,7 +166,7 @@ 
Examplessf_use_s2(TRUE)
 system.time(gck4a.nb <- knn2nb(knearneigh(xy1, k=4)))
 #>    user  system elapsed 
-#>   0.013   0.000   0.013 
+#>   0.011   0.000   0.012 
 summary(gck4a.nb, xy1, scale=0.5)
 #> Neighbour list object:
 #> Number of regions: 48 
diff --git a/docs/reference/localC.html b/docs/reference/localC.html
index a265ab2b..d7b45e51 100644
--- a/docs/reference/localC.html
+++ b/docs/reference/localC.html
@@ -414,7 +414,7 @@ 
Examples    alternative="two.sided"))
 }
 #>    user  system elapsed 
-#>   0.370   0.009   0.386 
+#>   0.345   0.008   0.354 
 if (run) {
   if (require(parallel, quietly=TRUE)) {
     ncpus <- max(2L, detectCores(logical=FALSE), na.rm = TRUE)-1L
@@ -668,7 +668,7 @@ 
Examples    nsim=9999, alternative="two.sided"))
 }
 #>    user  system elapsed 
-#>   0.631   0.003   0.639 
+#>   0.611   0.004   0.617 
 if (run) {
   all.equal(lisa_values(lmc_vars), c(lmc_vars_spdep))
 }
diff --git a/docs/reference/mstree.html b/docs/reference/mstree.html
index 92a5809d..08ac4976 100644
--- a/docs/reference/mstree.html
+++ b/docs/reference/mstree.html
@@ -102,9 +102,26 @@ 
Author<
     

     
Examples

     
### loading data
-if (as.numeric_version(unname(sf_extSoftVersion()["GDAL"])) >= "3.7.0") {
-bh <- st_read(system.file("etc/shapes/bhicv.gpkg.zip",
-      package="spdep")[1], quiet=TRUE)
+(GDAL37 <- as.numeric_version(unname(sf_extSoftVersion()["GDAL"])) >= "3.7.0")
+#> [1] TRUE
+file <- "etc/shapes/bhicv.gpkg.zip"
+zipfile <- system.file(file, package="spdep")
+if (GDAL37) {
+    bh <- st_read(zipfile)
+} else {
+    td <- tempdir()
+    bn <- sub(".zip", "", basename(file), fixed=TRUE)
+    target <- unzip(zipfile, files=bn, exdir=td)
+    bh <- st_read(target)
+}
+#> Reading layer `bhicv' from data source 
+#>   `/tmp/RtmpPpr0yz/temp_libpathafa452560f839/spdep/etc/shapes/bhicv.gpkg.zip' 
+#>   using driver `GPKG'
+#> Simple feature collection with 98 features and 8 fields
+#> Geometry type: POLYGON
+#> Dimension:     XY
+#> Bounding box:  xmin: -45.02175 ymin: -20.93007 xmax: -42.50321 ymax: -18.08342
+#> Geodetic CRS:  Corrego Alegre 1970-72
 ### data padronized
 dpad <- data.frame(scale(as.data.frame(bh)[,5:8]))
 
@@ -119,17 +136,34 @@ 
Examples
 ### find a minimum spanning tree
 system.time(mst.bh <- mstree(nb.w,5))
+#>    user  system elapsed 
+#>   0.002   0.000   0.002 
 dim(mst.bh)
+#> [1] 97  3
 head(mst.bh)
+#>      [,1] [,2]      [,3]
+#> [1,]    5   12 1.2951120
+#> [2,]   12   13 0.6141101
+#> [3,]   13   11 0.7913745
+#> [4,]   13    6 0.9775650
+#> [5,]   11   31 0.9965625
+#> [6,]   31   39 0.6915158
 tail(mst.bh)
+#>       [,1] [,2]      [,3]
+#> [92,]   89   90 2.5743702
+#> [93,]   26   56 2.6235317
+#> [94,]   86   87 2.6471303
+#> [95,]   87   72 0.7874461
+#> [96,]   49   36 2.8743677
+#> [97,]   24   25 3.4675168
 ### the mstree plot
 par(mar=c(0,0,0,0))
 plot(st_geometry(bh), border=gray(.5))
 plot(mst.bh, st_coordinates(st_centroid(bh)), col=2, 
      cex.lab=.6, cex.circles=0.035, fg="blue", add=TRUE)
-}
 #> Warning: st_centroid assumes attributes are constant over geometries
 
+