I use the package
maptools to read in Shapefiles in R. I would like to know if there is a function that can built up a query/subset to select only some Polygons based on a specific attribute in the a.t. (for example one state out of a shapefile with all states in one country). And, can I save this subset as a spatial object to do further analyses in R?
I use the package
To summarize the very useful comments by @Spacedman, @JeffreyEvans and @AriB.Friedman, but also to address the questioner's second problem about saving the generated subset and to provide a reproducible example for future audience, let me add a brief example to cover the whole topic. In fact, I'm much more often working with raster data than point or polygon shapefiles, so this was also pretty useful for me.
As @Spacedman mentioned, you should rather use
rgdal instead of
maptools for the given reason.
Now, let's download a suitable ESRI shapefile containing all states of the US and extract it to a folder named 'data' in your current working directory. After that, you are ready to import the shapefile using
usa.states <- readOGR(dsn = "data", layer = "states")
You can now create a subset of this SpatialPolygonsDataFrame either via accessing the
@data slot (thanks @JeffreyEvans, I didn't know that...) or addressing a particular column solely via the
$ operator. Let's take New York for our testing purposes:
newyork <- usa.states[usa.states$STATE_NAME == "New York", ] # Or alternatively... :-) newyork <- usa.states[usa.states@data$STATE_NAME == "New York", ]
Have a look at the generated subset via
plot(newyork). Now that you know that everything worked, use
writeOGR() to save the shapefile to the desired folder ('data' in our case). If you're not sure whether 'ESRI Shapefile' or any other driver is available on your local system, check
ogrDrivers() prior to sending the
writeOGR() statement to the console.
writeOGR(oregon, dsn = "data", layer = "oregon", driver = "ESRI Shapefile", overwrite_layer = TRUE)
That's it. Hope that helps!
ok this was easier than i thought:
spatial.subset<-states.shp[5,] if you know which line in your attribute table should be the subset. If you don't know you could do something like
spatial.subset<-states.shp[states.shp@data$STATE=="California" ,] for a single state or
spatial.subset<-states.shp[c(states.shp@data$POPULATION>1000000),] for several states with a population greater than 1.Mio.