[R-sig-Geo] Subsetting SpatialPolygonsDataFrame and empty ID slots

[Agustin Lobo]

That's really helpful, thanks.

For "everyday" use, I've made:

"slotspol" <- function(pol, sl="ID")
sapply(slot(pol, "polygons"), function(x) slot(x, sl))

Which can be used for any slot name:
slotspol(Up)
slotspol(Up,"area")

Also, few questions:
1. I observe that pols$Site is identical to pols at data$Site. Is not 
pols at data$Site the correct syntax? I mean that we retrieve variable 
"Site" from slot "data" in object "pols".

2. In your code, I don't understand
sites_Lupuna <- length(grep("Lupuna", as.character(pols$Site))) > 0

because the result is just T (the result of length(grep("Lupuna", 
as.character(pols$Site)))  is 3)
and thus subsetting with sites_Lupuna yields an object that is identical
to the original:

 > pols_Lupuna <- pols[sites_Lupuna,]
 > length(pols_Lupuna at polygons)
[1] 140
 > length(pols at polygons)
[1] 140

Instead:
 > sites_Lupuna <- grep("Lupuna", as.character(pols$Site))
 > sites_Lupuna
[1] 1 2 3
 > pols_Lupuna <- pols[sites_Lupuna,]

seems to provide a correctly subset of the original pols object ("Lupuna 
1", "Lupuna 2", "Lupuna 3")

Am I wrong? I feel quite insecure on this.

Agus

Roger Bivand escribió:
> On Mon, 29 Oct 2007, Agustin Lobo wrote:
> 
>> Dear list,
>>
>> I have imported a shp file with:
>>
>>> pols <- readOGR("../AllTransectPolygons02",
>> layer="AllTransectPolygons02")
>>
>>> class(pols)
>> [1] "SpatialPolygonsDataFrame"
>> attr(,"package")
>> [1] "sp"
>>
>>> pols at data
>>     AREA PERIMETER               Site
>> 0   2500      1010           Lupuna 1
>> 1   2500      1010           Lupuna 3
>> 2   2500      1010           Lupuna 2
>> 3   2500      1010             Aleman
>> 4   2500      1010        Zancudoyacu
>> etc
>>
>> and now I want to select few polygons according
>> to a vector of Site names.
>>
>> Thinking on how doing this, I thought I could select
>> the corresponding IDs in the data table and then use
>> the ID slot in the polygons, but
>> I've found that the polygons themselves actually
>> lack the ID information:
>>
> 
> Please treat all the long named access thingies in sp as deprecated - 
> they will be flagged as such in the next release, and removed one 
> release after that. S4 is not like that.
> 
> library(sp)
> library(rgdal)
> example(readOGR)
> class(Up)
> getSlots(class(Up))
> class(slot(Up, "polygons")[[1]])
> getSlots(class(slot(Up, "polygons")[1]))
> 
> So:
> 
> sapply(slot(Up, "polygons"), function(x) slot(x, "ID"))
> 
>>> getPolygonsIDSlot(pols)
>> Error in getPolygonsIDSlot(pols) :
>>   no slot of name "ID" for this object of class 
>> "SpatialPolygonsDataFrame"
>>
>> How does sp relate the data table to the polygons if there is no ID?
>>
> 
> Of course there is an ID, but you can get at either by:
> 
> sapply(slot(pols, "polygons"), function(x) slot(x, "ID"))
> 
> or
> 
> rownames(as(pols, "data.frame"))
> 
> 
>> In any case, which would the most efficient way of subseting
>> the SpatialPolygonsDataFrame into another SpatialPolygonsDataFrame
>> with the selected polygons?
> 
> Returning to the original question, subsetting is simply by the "[" 
> operator, just like any other data.frame:
> 
> sites_Lupuna <- length(grep("Lupuna", as.character(pols$Site))) > 0
> summary(sites_Lupuna)
> pols_Lupuna <- pols[sites_Lupuna,]
> 
> Just think of Spatial*DataFrame objects as data.frame objects and things 
> will be much clearer.
> 
> Roger
> 
>> Should I convert to a maptools object
>> and use subset.polylist()? If so, how do I convert? I've tried:
>>
>> pols.maptools <- SpatialPolygons2PolySet(pols)
>>
> 
> PS. PolySet objects are really for drawing coastlines in this context.
> 
>> but the data table is lost.
>>
>> Thanks!
>>
>>
> 

-- 
Dr. Agustin Lobo
Institut de Ciencies de la Terra "Jaume Almera" (CSIC)
LLuis Sole Sabaris s/n
08028 Barcelona
Spain
Tel. 34 934095410
Fax. 34 934110012
email: Agustin.Lobo at ija.csic.es
http://www.ija.csic.es/gt/obster