Work with individual Polygons of a SpatialPolygonsDataFrame (R)

Question

I'm trying to work with individual polygons of a SpatialPolygonsDataFrame (SPDF). I've seen these two questions: Using sapply on SPDF and a for-loop on a SPDF. The latter I thought may be helpful, but I don't understand the %>%-syntax and don't know how to apply it to conditional code.

I want to iterate over an SPDF, shrink (gBuffer) polygons, as long as they don't get smaller than 1 ha (10'000 m²), then maybe apply other steps to those polygons. Something like this:

# get shp file with multiple polygons, creates SpatialPolygonsDataFrame
shpA = shapefile(file.path(PATH_DATA, 'sample_test_A.shp'))

# iterate over shp and work with individual polygons
for (i in 1:length(shpA)) {
  # THIS DOESN'T GIVE ME A SPATIAL POLYGON I CAN WORK WITH:
  pol_temp = shpA[[1]] # or shpA@polygons[[i]]
  
  # create buffered polygon and check if it is > 1ha
  pol_buffered = gBuffer(shpA, width=-30)
  if(area(pol_buffered) > 10000){
    # replace polygon
    shpA[[1]] == pol_buffered
  }
  
  # maybe do some other stuff, like sampling the polygon
  ...
}

Answer 1

The way you subset polygons using a bracket index is not with a double bracket try shpA[1,]. You also do not need to do this in a for loop.

Even though it is not polygon data, we can use the meuse data to illustrate

library(sp)
library(rgeos)

data(meuse)
  coordinates(meuse) <- ~x+y

Create a variable buffer for the meuse points

b <- gBuffer(meuse, byid = TRUE, width = sample(10:500, nrow(meuse)))

Calculate areas for each polygon buffer

( a <- gArea(b, byid=TRUE) )

Now, to subset we can use which to create an index representing our query. If passed to a bracket index this can be used to subset the data. The if statement is helpful in avoiding the creation of an empty dataset due to a zero length the index resulting from no features meeting the query criteria.

( idx <- which( a > 100000) )
  if(length(idx) > 0)
    b.sub <- b[idx,]

If you know that your query criteria will return results, this could be distilled to one line of code (note; do not run if you already ran the above example).

b.sub <- b[which( gArea(b, byid=TRUE) > 100000),]

Check results

dim(b)
summary(gArea(b, byid=TRUE))

dim(b.sub)
summary(gArea(b.sub, byid=TRUE))

Now, across your two datasets, replacing polygons in shpA that meet the criteria of pol_buffered.

pol_buffered = gBuffer(shpA, byid=TRUE, width=-30)
  idx <- which(gArea(pol_buffered, byid=TRUE) > 10000)
shpA[idx,] <- pol_buffered[idx,]