2

I'm trying to get the area of all the polygons in a SpatialPolygonsDataFrame: parcProj. 

library("rgeos")
gArea(parcProj[2,])

works fine but:

apply(parcProj,1,gArea)

or 

sapply(parcProj@polygons,gArea)

don't.

Any standard way to do an apply besides the obvious:

for(i in 1:length(parcProj)) areas <- c(areas,gArea(parcProj[i,]))
Improve this question

2 Answers 2

Reset to default
3

apply is not really the right tool, since this is a data frame (of sorts) not an array.

How about 

gArea(parcProj, byid = TRUE)
Improve this answer
0

Not that @mdsumner's answer is not perfectly correct, and honestly easier to understand, but why not just pull the areas from the "polygons" and "area" slot(s)? 

sapply(slot(parcProj, "polygons"), function(i) slot(i, "area"))

I do notice a failure on both slot values and values returned by gArea in relation to holes (No Data area) verses real polygons. As illustrated by this code, gArea returns the correct area for a polygon with a hole whereas the area slot is not correct (as expected). However, when the "hole" is an actual polygon neither gArea or the area slot return the correct values. This hints at a bug where gArea is not returning the correct area with nested polygons. 

library(rgeos)
library(sp)

Sr1 <- Polygon(cbind(c(2,4,4,1,2),c(2,3,5,4,2)))
Sr2 <- Polygon(cbind(c(5,4,2,5),c(2,3,2,2)))
Sr3 <- Polygon(cbind(c(4,4,5,10,4),c(5,3,2,5,5)))
Sr4 <- Polygon(cbind(c(5,6,6,5,5),c(4,4,3,3,4)), hole = TRUE)

Srs1 <- Polygons(list(Sr1), "s1")
Srs2 <- Polygons(list(Sr2), "s2")
Srs3 <- Polygons(list(Sr3, Sr4), "s3hole")

#### with hole
poly.hole <- SpatialPolygons(list(Srs1,Srs2,Srs3), 1:3)

sapply(slot(poly.hole, "polygons"), function(i) slot(i, "area"))
gArea(poly.hole, byid = TRUE)

#### without hole
Srs3 <- Polygons(list(Sr3), "s3")
poly.nohole <- SpatialPolygons(list(Srs1,Srs2,Srs3), 1:3)

sapply(slot(poly.nohole, "polygons"), function(i) slot(i, "area"))
gArea(poly.nohole, byid = TRUE)

#### with real polygon in place of hole
Sr4 <- Polygon(cbind(c(5,6,6,5,5),c(4,4,3,3,4)))
Srs4 <- Polygons(list(Sr4), "s4")
poly <- SpatialPolygons(list(Srs1,Srs2,Srs3,Srs4), 1:4)

sapply(slot(poly, "polygons"), function(i) slot(i, "area"))
gArea(poly, byid = TRUE)

# Plot data
par(mfrow=c(2,2))
  plot(poly.hole, col=1:3, pbg="white");title("hole")
  plot(poly.nohole, col=1:3);title("no hole")
  plot(poly, col=1:4);title("real poly")
Improve this answer
4
  • My understanding is that the area slot is not accounting for holes in polygons: it is attached to a given polygon as set of points (including a polygon drawing a hole)
    – cmbarbu
    Commented May 27, 2015 at 0:24
  • You are correct. However, I did find that there is a possibility that, with nested polygons, gArea is returning the incorrect area as well. I should note that I rarely see true no data areas (holes) in polygon data but polygon nesting is very common.
    – Jeffrey Evans
    Commented May 27, 2015 at 1:33
  • 1
    This is not a bug: in your "real polygon in place of hole" you are not defining the extra polygon as a hole. It is desirable that the polygons can be superimposed if not specified that one is the hole of an other.
    – cmbarbu
    Commented May 27, 2015 at 8:41
  • @cmbarbu, this is interesting, I was wondering if this was the case.
    – Jeffrey Evans
    Commented May 27, 2015 at 13:34

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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