gBuffer(rGEOS) causing R to crash

Question

I have downloaded a US federal lands data layer (fedlanp010g.shp.tar.gz) from http://nationalmap.gov/small_scale/atlasftp.html?openChapters=chpbound#chpbound

In R, I am trying to create a buffer around each polygon, However gBuffer crashes when the code below is run. I added a pause (sys.sleep(1)) to add a delay, and I received the error message "Value of SET_STRING_ELT() must be "CHARSXP" not a 'NULL'" at row 14.

Here is the code I have been using (I am processing with a loop):

library(rgeos)
library(raster)

fed = readOGR("downloaded file from above here")
i = 1
temp = fed[i,]
buff = gBuffer(temp,byid=TRUE,width=25)
result = erase(buff,temp)

for(i in 2:nrow(fed)){
    temp = fed[i,]
    buff = gBuffer(temp,byid=TRUE,width=25)
    temp_buff = erase(buff,temp)
    result = rbind(result,temp_buff)
  } 

Answer 1

Is R actually crashing or is the code just failing? These are very different outcomes. I believe that your problem was that your data is in a geographic and not planar projection. I reprojected the data, from your link, to the USGS definition of Albers and gBuffer ran with no problems. You will need rgdal along with sp for the projection transformation.

library(rgeos)
library(sp)
library(rgdal)

setwd("D:/TMP/test")
fed <- readOGR(getwd(), "fedlanp010g")

fed <- spTransform(fed, CRS( " +proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23 +lon_0=-96 +x_0=0 +y_0=0+datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0" ) )  

fed.buf <- gBuffer(fed, byid = TRUE, width = 25)

You may be getting into trouble with how the for loop is defined and the use of erase. Since you did not provide a reproducible example it is difficult to say what the exact issue is as it could be package, code or data related.

Try this code, where each buffered polygon is stored in a list and then combined using do.call with rbind. However, in this case, I am not sure what happens with intersecting geometry.

First, add required libraries and create some example polygon data.

library(rgeos)
library(sp)

p1=Polygons(list(Polygon(cbind(c(2,4,4,1,2),c(2,3,5,4,2)))), "1")
p2=Polygons(list(Polygon(cbind(c(5,4,2,5),c(2,3,2,2)))), "2")
p3=Polygons(list(Polygon(cbind(c(4,4,5,10,4),c(5,3,2,5,5)))), "3") 
fed = SpatialPolygonsDataFrame(SpatialPolygons(list(p1,p2,p3)), data.frame(ID=c("1","2","3")))

Here is a for loop that buffers each polygon, adds to a list object and then uses do.call to rbind the buffered polygons together.

result <- list()
for(i in 1:length(fed)){result[[i]] <- gBuffer(fed[i,], byid=TRUE, width = -0.05)}
result <- do.call("rbind", result)

plot(result) 

I am wondering why you are applying a for loop. You get exactly the same results just using gBuffer with the entire spatial object and the byid=TRUE argument.

fed.buf <- gBuffer(fed, byid=TRUE, width = -0.05)
plot(fed.buf)