简单的方法通过R中的属性来对SpatialPolygonsDataFrame进行子集(即删除多边形)



mapping maptools (4)

我想简单地从SpatialPolygonsDataFrame对象中删除一些基于@data数据框中相应属性值的多边形,这样我就可以绘制一个简化/子集的shapefile。 到目前为止,我还没有找到办法做到这一点。

例如,假设我想从这个世界shapefile中删除面积小于30000的所有多边形。我该如何去做这件事?

或者,同样,我怎样才能删除Antartica?

require(maptools)

getinfo.shape("TM_WORLD_BORDERS_SIMPL-0.3.shp") 
# Shapefile type: Polygon, (5), # of Shapes: 246
world.map <- readShapeSpatial("TM_WORLD_BORDERS_SIMPL-0.3.shp")

class(world.map)
# [1] "SpatialPolygonsDataFrame"
# attr(,"package")
# [1] "sp"

head(world.map@data)
#   FIPS ISO2 ISO3 UN                NAME   AREA  POP2005 REGION SUBREGION     LON     LAT
# 0   AC   AG  ATG 28 Antigua and Barbuda     44    83039     19        29 -61.783  17.078
# 1   AG   DZ  DZA 12             Algeria 238174 32854159      2        15   2.632  28.163
# 2   AJ   AZ  AZE 31          Azerbaijan   8260  8352021    142       145  47.395  40.430
# 3   AL   AL  ALB  8             Albania   2740  3153731    150        39  20.068  41.143
# 4   AM   AM  ARM 51             Armenia   2820  3017661    142       145  44.563  40.534
# 5   AO   AO  AGO 24              Angola 124670 16095214      2        17  17.544 -12.296

如果我做这样的事情,情节并不反映任何变化。

world.map@data = world.map@data[world.map@data$AREA > 30000,]
plot(world.map)

同样的结果,如果我这样做:

world.map@data = world.map@data[world.map@data$NAME != "Antarctica",]
plot(world.map)

任何帮助表示赞赏!

https://src-bin.com


Answer #1

作为第二个指针:这不适用于形状中带有“孔”的shapefile,因为它是按索引进行子集的。


Answer #2

只要提到这个subset也可以避免在条件中写入数据的名字。

world.map <- subset(world.map, AREA > 30000)
plot(world.map)

Answer #3

我用上面的技术来制作一个澳大利亚的地图:

australia.map < - world.map[world.map$NAME == "Australia",]
plot(australia.map)

事实证明,“澳大利亚”之后的逗号很重要。

这种方法的一个缺点是,它似乎保留所有其他国家的所有属性的列和行,并填充零。 我发现如果我写出一个.shp文件,然后使用readOGR(rgdal包)读取它,它会自动删除空的地理数据。 然后,我可以写只有我想要的数据的另一个形状文件。

writeOGR(australia.map,".","australia",driver="ESRI Shapefile")
australia.map < - readOGR(".","australia")
writeOGR(australia.map,".","australia_small",driver="ESRI Shapefile")

在我的系统上,至少是删除空数据的“读取”函数,所以我必须在读完文件后再写一次(如果我尝试重新使用文件名,则会出现错误)。 我相信有一个更简单的方法,但是这似乎对我的目的来说足够好了。


Answer #4

看起来像是覆盖数据,但不能删除多边形。 如果你想减少包括数据和多边形的数据集,请尝试例如

world.map <- world.map[world.map$AREA > 30000,]
plot(world.map)

[[编辑2016年4月19日]]这个解决方案曾经工作,但是@Bonnie报告了一个更新的R版本(虽然也许数据已经改变了吗?): world.map <- world.map[[email protected]$AREA > 30000, ] Upvote @ Bonnie的答案,如果有帮助。