在 R 中查找多边形的最近邻
Find nearest neighbours of polygons in R
我有一个坐标数据框，我已将其转换为 R 中的 sf 对象，如下所示：
> head(df1)
Cell_ID Spot_ID X Y
1 0 600000000 193.722 175.733
2 0 600000001 192.895 176.727
3 0 600000002 193.828 177.462
4 8 600000003 178.173 178.220
5 7 600000004 187.065 178.285
6 0 600000005 190.754 178.186
> df1_sf <- st_as_sf(df1,
coords = c('X', 'Y')) %>%
group_by(Cell_ID) %>%
summarise() %>%
ungroup() %>%
st_convex_hull()
>plot(st_geometry(df1_sf), border = "red")
然后我可以绘制我所有的多边形，它看起来像这样：
现在我想获取每个多边形的邻居的 ID。为此，我正在做
n = st_set_geometry(st_intersection(df1_sf,df1_sf), NULL)
head(n)
# A tibble: 6 x 2
Cell_ID Cell_ID.1
<int> <int>
1 0 0
2 7 0
3 51 0
4 1 1
5 4 1
6 5 1
但这是一项平庸的工作，因为它需要一个交叉点，而如果它们是最近的交叉点，我也对它们很感兴趣（关闭但不像下面的图片那样接触，Cell_ID 1 将具有邻居单元格 3-6 但也会检测到单元格 7，因为它位于给定的半径内）。
谁能帮我解决这个问题？
谢谢！！
从你的问题来看，你似乎对通用的最近邻类型的方法更感兴趣。如果这过于简单化，请纠正我。
无需考虑每个多边形及其边界，您可以简单地获取中心坐标并使用任何 knn 类型的算法将 k nearest neighbours 分类为给定坐标。
由于我无法访问您的数据，因此我创建了一些虚拟坐标。
使用包 RANN 和函数 nn2 see here.
install.packages('RANN')
library(RANN)
# Make dummy coordinates
df <-
data.frame( X = runif(100)
, Y = runif(100)
)
# Find closest 5 points between df and itself
closest <- nn2(data = df, query = df , k = 5)
closest$nn.idx # Index of Closest neigbours
closest$nn.dists # Euclidean distance of Closest neigbours
# Note the first colum is a reference to itself, so real 5 nearest neighbours (not including itself) would mean you select k = 6.
> head(closest$nn.idx) # Euclidean distance of Closest neigbours
[,1] [,2] [,3] [,4] [,5]
[1,] 1 82 31 86 49
[2,] 2 22 41 34 91
[3,] 3 96 20 55 32
[4,] 4 65 53 77 14
[5,] 5 38 48 59 30
[6,] 6 36 43 97 61
> head(closest$nn.dists) # Euclidean distance of Closest neigbours
[,1] [,2] [,3] [,4] [,5]
[1,] 0 0.04971692 0.06305752 0.08597908 0.09485483
[2,] 0 0.03668956 0.05248395 0.09570358 0.10489092
[3,] 0 0.07257007 0.10263107 0.11204297 0.13275642
[4,] 0 0.07209561 0.07227328 0.07259919 0.07326718
[5,] 0 0.02842711 0.06003873 0.08930219 0.12286905
[6,] 0 0.08018734 0.09312385 0.10844622 0.11368332
您也可以使用 searchtype = "radius" 和 radius.
根据问题中提到的半径方法执行此操作
为了说明在每个多边形周围使用缓冲区的出色建议
（每个多边形的数学膨胀）这是一个快速而肮脏的 spatstat
解决方案。
首先加载包并制作一些示例数据：
library(spatstat)
dat <- tiles(dirichlet(cells))
ii <- seq(2, 42, by=2)
dat[ii] <- lapply(dat[ii], erosion, r = .01)
dat <- lapply(seq_along(dat), function(i) cbind(Cell_ID = i, as.data.frame(dat[[i]])))
dat <- Reduce(rbind, dat)
df1 <- cbind(Spot_ID = 1:nrow(dat), dat)
head(df1)
#> Spot_ID Cell_ID x y
#> 1 1 1 0.4067780 0.0819020
#> 2 2 1 0.3216680 0.1129640
#> 3 3 1 0.1967080 0.0000000
#> 4 4 1 0.4438430 0.0000000
#> 5 5 2 0.5630909 0.1146781
#> 6 6 2 0.4916145 0.1649979
拆分每个 Cell_ID，找到凸包并绘制数据：
dat <- split(df1[,c("x", "y")], df1$Cell_ID)
dat <- lapply(dat, convexhull)
plot(owin(), main = "")
for(i in seq_along(dat)){
plot(dat[[i]], add = TRUE, border = "red")
}
扩大每个多边形：
bigdat <- lapply(dat, dilation, r = 0.0125)
天真地 for-loop 分配哪些扩张的多边形重叠（即完整
n^2 成对交集):
neigh <- list()
for(i in seq_along(bigdat)){
overlap <- sapply(bigdat[-i], function(x) !is.empty(intersect.owin(x, bigdat[[i]])))
neigh[[i]] <- which(overlap)
}
绘制具有邻居数量的扩张多边形（邻居的 ID 在
列表 neigh):
plot(owin(), main = "")
for(i in seq_along(bigdat)){
plot(bigdat[[i]], add = TRUE, border = "red")
}
text.ppp(cells, labels = sapply(neigh, length))
基于曲面细分的替代解决方案
是否要求使用convexhull作为cell的定义
地区？我很想简单地用质心代表每个细胞
样本点，然后使用 Dirichlet/Voronoi 镶嵌作为
地区。它们到处都有 well-defined 个邻居，唯一的问题是
如何定义单元格集合的边界区域。
拆分每个 Cell_ID，找到质心，细分并绘制数据：
dat <- split(df1[,c("x", "y")], df1$Cell_ID)
dat <- t(sapply(dat, colMeans))
X <- as.ppp(dat, W = ripras)
D <- dirichlet(X)
plot(D)
查找邻居 ID 的额外代码：
eps <- sqrt(.Machine$double.eps) # Epsilon for numerical comparison below
tilelist <- tiles(D)
v_list <- lapply(tilelist, vertices.owin)
v_list <- lapply(v_list, function(v){ppp(v$x, v$y, window = Window(X), check = FALSE)})
neigh <- list()
dd <- safedeldir(X)
for(i in seq_len(npoints(X))){
## All neighbours from deldir (infinite border tiles)
all_neigh <- c(dd$delsgs$ind1[dd$delsgs$ind2==i],
dd$delsgs$ind2[dd$delsgs$ind1==i])
## The remainder keeps only neighbour tiles that share a vertex with tile i:
true_neigh <- sapply(v_list[all_neigh], function(x){min(nncross.ppp(v_list[[i]], x))}) < eps
neigh[[i]] <- sort(all_neigh[true_neigh])
}
plot(D, main = "Tessellation with Cell_ID")
text(X)
neigh[[1]] # Neighbours of tile 1
#> [1] 2 7 8
neigh[[10]] # Neighbours of tile 10
#> [1] 3 4 5 9 15 16 20
我有一个坐标数据框，我已将其转换为 R 中的 sf 对象，如下所示：
> head(df1)
Cell_ID Spot_ID X Y
1 0 600000000 193.722 175.733
2 0 600000001 192.895 176.727
3 0 600000002 193.828 177.462
4 8 600000003 178.173 178.220
5 7 600000004 187.065 178.285
6 0 600000005 190.754 178.186
> df1_sf <- st_as_sf(df1,
coords = c('X', 'Y')) %>%
group_by(Cell_ID) %>%
summarise() %>%
ungroup() %>%
st_convex_hull()
>plot(st_geometry(df1_sf), border = "red")
然后我可以绘制我所有的多边形，它看起来像这样：
现在我想获取每个多边形的邻居的 ID。为此，我正在做
n = st_set_geometry(st_intersection(df1_sf,df1_sf), NULL)
head(n)
# A tibble: 6 x 2
Cell_ID Cell_ID.1
<int> <int>
1 0 0
2 7 0
3 51 0
4 1 1
5 4 1
6 5 1
但这是一项平庸的工作，因为它需要一个交叉点，而如果它们是最近的交叉点，我也对它们很感兴趣（关闭但不像下面的图片那样接触，Cell_ID 1 将具有邻居单元格 3-6 但也会检测到单元格 7，因为它位于给定的半径内）。 谁能帮我解决这个问题？
谢谢！！
从你的问题来看，你似乎对通用的最近邻类型的方法更感兴趣。如果这过于简单化，请纠正我。
无需考虑每个多边形及其边界，您可以简单地获取中心坐标并使用任何 knn 类型的算法将 k nearest neighbours 分类为给定坐标。
由于我无法访问您的数据，因此我创建了一些虚拟坐标。
使用包 RANN 和函数 nn2 see here.
install.packages('RANN')
library(RANN)
# Make dummy coordinates
df <-
data.frame( X = runif(100)
, Y = runif(100)
)
# Find closest 5 points between df and itself
closest <- nn2(data = df, query = df , k = 5)
closest$nn.idx # Index of Closest neigbours
closest$nn.dists # Euclidean distance of Closest neigbours
# Note the first colum is a reference to itself, so real 5 nearest neighbours (not including itself) would mean you select k = 6.
> head(closest$nn.idx) # Euclidean distance of Closest neigbours
[,1] [,2] [,3] [,4] [,5]
[1,] 1 82 31 86 49
[2,] 2 22 41 34 91
[3,] 3 96 20 55 32
[4,] 4 65 53 77 14
[5,] 5 38 48 59 30
[6,] 6 36 43 97 61
> head(closest$nn.dists) # Euclidean distance of Closest neigbours
[,1] [,2] [,3] [,4] [,5]
[1,] 0 0.04971692 0.06305752 0.08597908 0.09485483
[2,] 0 0.03668956 0.05248395 0.09570358 0.10489092
[3,] 0 0.07257007 0.10263107 0.11204297 0.13275642
[4,] 0 0.07209561 0.07227328 0.07259919 0.07326718
[5,] 0 0.02842711 0.06003873 0.08930219 0.12286905
[6,] 0 0.08018734 0.09312385 0.10844622 0.11368332
您也可以使用 searchtype = "radius" 和 radius.
为了说明在每个多边形周围使用缓冲区的出色建议
（每个多边形的数学膨胀）这是一个快速而肮脏的 spatstat
解决方案。
首先加载包并制作一些示例数据：
library(spatstat)
dat <- tiles(dirichlet(cells))
ii <- seq(2, 42, by=2)
dat[ii] <- lapply(dat[ii], erosion, r = .01)
dat <- lapply(seq_along(dat), function(i) cbind(Cell_ID = i, as.data.frame(dat[[i]])))
dat <- Reduce(rbind, dat)
df1 <- cbind(Spot_ID = 1:nrow(dat), dat)
head(df1)
#> Spot_ID Cell_ID x y
#> 1 1 1 0.4067780 0.0819020
#> 2 2 1 0.3216680 0.1129640
#> 3 3 1 0.1967080 0.0000000
#> 4 4 1 0.4438430 0.0000000
#> 5 5 2 0.5630909 0.1146781
#> 6 6 2 0.4916145 0.1649979
拆分每个 Cell_ID，找到凸包并绘制数据：
dat <- split(df1[,c("x", "y")], df1$Cell_ID)
dat <- lapply(dat, convexhull)
plot(owin(), main = "")
for(i in seq_along(dat)){
plot(dat[[i]], add = TRUE, border = "red")
}
扩大每个多边形：
bigdat <- lapply(dat, dilation, r = 0.0125)
天真地 for-loop 分配哪些扩张的多边形重叠（即完整 n^2 成对交集):
neigh <- list()
for(i in seq_along(bigdat)){
overlap <- sapply(bigdat[-i], function(x) !is.empty(intersect.owin(x, bigdat[[i]])))
neigh[[i]] <- which(overlap)
}
绘制具有邻居数量的扩张多边形（邻居的 ID 在
列表 neigh):
plot(owin(), main = "")
for(i in seq_along(bigdat)){
plot(bigdat[[i]], add = TRUE, border = "red")
}
text.ppp(cells, labels = sapply(neigh, length))
基于曲面细分的替代解决方案
是否要求使用convexhull作为cell的定义 地区？我很想简单地用质心代表每个细胞 样本点，然后使用 Dirichlet/Voronoi 镶嵌作为 地区。它们到处都有 well-defined 个邻居，唯一的问题是 如何定义单元格集合的边界区域。
拆分每个 Cell_ID，找到质心，细分并绘制数据：
dat <- split(df1[,c("x", "y")], df1$Cell_ID)
dat <- t(sapply(dat, colMeans))
X <- as.ppp(dat, W = ripras)
D <- dirichlet(X)
plot(D)
查找邻居 ID 的额外代码：
eps <- sqrt(.Machine$double.eps) # Epsilon for numerical comparison below
tilelist <- tiles(D)
v_list <- lapply(tilelist, vertices.owin)
v_list <- lapply(v_list, function(v){ppp(v$x, v$y, window = Window(X), check = FALSE)})
neigh <- list()
dd <- safedeldir(X)
for(i in seq_len(npoints(X))){
## All neighbours from deldir (infinite border tiles)
all_neigh <- c(dd$delsgs$ind1[dd$delsgs$ind2==i],
dd$delsgs$ind2[dd$delsgs$ind1==i])
## The remainder keeps only neighbour tiles that share a vertex with tile i:
true_neigh <- sapply(v_list[all_neigh], function(x){min(nncross.ppp(v_list[[i]], x))}) < eps
neigh[[i]] <- sort(all_neigh[true_neigh])
}
plot(D, main = "Tessellation with Cell_ID")
text(X)
neigh[[1]] # Neighbours of tile 1
#> [1] 2 7 8
neigh[[10]] # Neighbours of tile 10
#> [1] 3 4 5 9 15 16 20