Search code examples
rspatialr-maptoolsr-sp

get an empty SpatialPolygonsDataFrame via subset?

I'm looking to subset a SpatialPolygonsDataFrame by an attribute, but I want to allow it to return an empty SpatialPolygonsDataFrame.

If we are to treat objects of type SpatialPolygonsDataFrame like data.frames, as discussed here, we should be able to get and work with empty objects.

I'm interested because I want to incorporate this into a function that may try to subset by an attribute that will essentially pick no features. 

owd <- getwd() 
setwd(system.file("shapes", package = "maptools")) 
library(maptools) 
nc90 <- readShapeSpatial("co37_d90") 
setwd(owd)

nc90@data[nc90@data$AREA>0.15,]               # returns data.frame
bigctys <- nc90[nc90@data$AREA>0.15,]         # SpatialPolygonsDataFrame 
nc90@data[nc90@data$AREA>0.25,]               # returns empty data.frame
bigestctys <- nc90[nc90@data$AREA>0.25,]      # ERROR

Is there a way to make this work? If not, is there a way to initalize an empty SpatialPolygonsDataFrame object? The future actions I want to perform on such an object involve over plotting on an existing map, so I'd like the image to be produced anyways, even if blank.

Solution

  • Right now you can't. This is somewhat inconsistent, as for SpatialPointsDataFrame objects you can:

    library(sp)
demo(meuse, ask = FALSE)
x = meuse[F,]

    although with warnings; also, validObject(x) returns FALSE, so they are intended to be not allowed!

    It's a bit abstract what such objects should represent, but I can see the analogy with data.frame objects with zero rows: it is useful that they can exist.