在 R 中合并和分解多边形时的问题

Question

有时，当使用 sp 和 rgeos 包在 R 中使用 SpatialPolygons 时，我在子集和分解这些多边形时会遇到问题。 此示例显示了我遇到的问题类型，尽管我也遇到过其他数据集。

require(maptools)
require(rgeos)

data(wrld_simpl)

countries <- c("Argentina","Bolivia","Brazil","Chile","Colombia","Ecuador","Guyana","Paraguay","Peru","Suriname","Uruguay","Venezuela")
SAmerica <- subset(wrld_simpl,wrld_simpl@data$NAME == countries[1])

for (i in 2:length(countries)) {
    x <- subset(wrld_simpl, wrld_simpl@data$NAME == countries[i])
    SAmerica <- gUnion(SAmerica, x)
}

plot(SAmerica)
gIsValid(SAmerica)

shoreline <- gUnaryUnion(SAmerica)
plot(shoreline)

多边形相互接壤的区域有可见的线段，这个联合的多边形在 rgeos 包中是无效的。

我很好奇 (1) 为什么会发生这种情况？以及 (2) 解决此类问题的最佳方法是什么。在这里，我绝对可以为国家矢量文件找到更好的数据，但有时我可能会使用更多独特的数据。我该如何清理它？

Answer 1

问题是 wrld_simpl 是一个糟糕的数据集。应该换成更好的数据，没有拓扑错误。

wrld_simpl 数据集中检测到的错误的流形 GIS 视图

您可以使用更好的数据尝试您的方法，结果会更好。

countries <- c("Bolivia","Brazil","Chile","Colombia","Ecuador","Guyana",
               "Paraguay","Peru","Suriname","Uruguay","Venezuela")

writeOGR(wrld_simpl, dsn = 'd:', layer = 'wrld_simpl', driver = 'ESRI Shapefile' )
# Correct topology (normalize topology with Manifol GIS
# reducing precision to 0.00001

从Dropbox下载拓扑归一化数据

sam_topo <- readOGR(dsn = 'd:/', layer = 'Samerica')

#plot(sam_topo, axes = T)
sam2 <- sam_topo[sam_topo$NAME == "Argentina", ]

for (i in countries) {
  x <- sam_topo[sam_topo$NAME == i,]
  sam2 <- gUnion(sam2, x)
}

gIsValid(sam2)
shoreline <- gUnaryUnion(sam2)
plot(shoreline)

Answer 2

这是使用raster 包中的getData(...) 执行此操作的简单方法。这避免了使用for 循环和外部软件。

library(rgeos)
library(raster)
world     <- getData("countries")
countries <- c("Argentina","Bolivia","Brazil","Chile","Colombia","Ecuador","Guyana","Paraguay","Peru","Suriname","Uruguay","Venezuela")
SAmerica  <- gUnaryUnion(world[world$ENGLISH %in% countries,])
plot(SAmerica)
gIsValid(SAmerica)
# [1] TRUE

getData(...) 加载的世界数据库在属性表中也有一个CONTINENTS 字段，在这种特定情况下，您可以使用以下方法获得完全相同的结果：

SAmerica  <- gUnaryUnion(world[world$CONTINENT=="South America",])