A simple way to subset a SpatialPolygonsDataFrame (i.e. remove polygons) by attribute in R

Question

I would just like to remove some polygons from the SpatialPolygonsDataFrame object based on the corresponding attribute values ​​in the @data data frame so that I can build a simplified / multiple necked file. Until I found a way to do this.

For example, let's say I want to remove all polygons from this world shapefile , which has an area of ​​less than 30,000. How can I go about this?

Or similarly, how can I remove 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 

If I do something like this, the plot will not reflect any changes.

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

same result if i do this:

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

Any help is appreciated!

Answer 1

it looks like you are overwriting data but not deleting polygons. If you want to reduce the data set, including data and polygons, try, for example,

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

[[Edit April 19, 2016]] This solution was used to work, but @Bonnie reports differently for a newer version of R (although maybe the data has also changed?): world.map <- world.map[world.map@data$AREA > 30000, ] Upvote @Bonnie answer if this helped.

Answer 2

When I tried to do this in R 3.2.1, the tim riffe technique above did not work for me, although changing it fixed the problem a bit. I found that I had to specifically refer to the data slot and then specify the attribute for the subset, as shown below:

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

Adding this as an alternative answer in case others are facing the same problem.

Answer 3

I used the above technique to make a map of Australia only:

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

The comma after "Australia" is important, as it turned out.

One of the drawbacks of this method is that it saves all columns and rows of attributes for all other countries and simply fills them with zeros. I found that if I wrote a .shp file and then read it using readOGR (rgdal package), it would automatically delete the null geographic data. Then I could write another form file with only the data I want.

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

On my system, at least it is a “read” function that deletes null data, so I have to write the file after reading it once (and if I try to reuse the file name, I get an error message). I'm sure there is an easier way, but this seems to work well enough for my purposes.

Answer 4

Just note that subset also forces the job to avoid writing the data name in state.

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

Answer 5

As a second pointer: does not work for shapefiles with “holes” in the shapes, because this is a subset of the index.