Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

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? Best

share|improve this question
4  
Comment: don't use maptools for reading, use readOGR from package:rgdal - it will preserve any coordinate reference info from the shape file. –  Spacedman May 12 '13 at 22:31
    
thanks for the hint! –  Dspanes May 13 '13 at 8:12
add comment

2 Answers

up vote 3 down vote accepted

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.

library(rgdal)

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 readOGR().

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!

share|improve this answer
    
thanks for this summary! –  Dspanes May 14 '13 at 9:59
add comment

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.

share|improve this answer
1  
This doesn't look like it's an attribute query. Perhaps you could elaborate on how it answers your question? –  whuber May 12 '13 at 16:28
    
@whuber you are right. Please note the edits in the answer to make it more comprehensive. –  Dspanes May 12 '13 at 20:40
1  
Using "which" could have unexpected behavior and is quite unnecessary. You should also be accessing the data slot in the bracket index. This will get you what you want. states.shp[states.shp@data$STATE=="California" ,] –  Jeffrey Evans May 13 '13 at 3:41
2  
@JeffreyEvans Why is accessing the @data slot necessary? $ is overloaded for Spatial*DataFrame objects to return the named column of the @data slot, and thus reaching inside the object's structure would seem to be disfavored over using the appropriate method.... –  Ari B. Friedman May 13 '13 at 11:46
1  
Because calling the slot is much more version proof. Until very recently you could not call a column in the data slot without explicitly accessing the slot object. –  Jeffrey Evans May 13 '13 at 17:06
show 1 more comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.