Intersecting points and overlapping polygons

I’ve been doing some spatial stuff of late and the next little step will involve intersecting points with possibly many overlapping polygons. The sp package has a function called over which returns the polygons that points intersects with. The catch though, is that it only returns the last (highest numerical value) polygon a point overlaps with. So it’s not so useful if you have many overlapping polygons. A little playing, and I’ve overcome that problem…

Here’s a toy example.

Create a couple of polygons and put them into a SpatialPolygons object.

library(sp)

p1 <- matrix(c(1,1,
 2,1,
 4,2,
 3,2), ncol = 2, byrow = TRUE)
p2 <- matrix(c(2.2,1,
 3,1,
 3,2,
 3,3,
 2.8,3), ncol = 2, byrow = TRUE)
p1s <- Polygons(list(Polygon(p1)), 3)
p2s <- Polygons(list(Polygon(p2)), 4)
sps <- SpatialPolygons(list(p1s, p2s))

Define a few points and put them in a SpatialPoints object

point <- matrix(c(2.5, 1.5,
 3.2, 1.75,
 2,3,
 1.5, 1.25, 
 2.75, 2.5), ncol = 2, byrow = TRUE)
points <- SpatialPoints(point)

We can plot them…(not shown)

plot(sps, border = c("black", "blue"))
plot(points, add = TRUE)

As here we have the issue:

over(points, sps)
 1  2  3  4  5 
 2  1 NA  1  2 

only returning a single “hit” per point (point 1 overlaps with both polygons 1 and 2).

To get around this, we can rip the individual polygons out of the SpatialPolygons object and put them in a list, converting the individual polygons into SpatialPolygons along the way:

sps2 <- lapply(sps@polygons, function(x) SpatialPolygons(list(x)))

Then lapply-ing over sps2 we can see which polygons each point intersects…

lapply(sps2, function(x) over(points, x))
[[1]]
 1  2  3  4  5 
 1  1 NA  1 NA 

[[2]]
 1  2  3  4  5 
 1 NA NA NA  1

And now we see that point one overlaps with both polygons.

Code >>>here<<<