1

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
  ...
}
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,] 

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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