I'm figuring out how to do a Intersection (Spatial Join) between point and polygons from shapefiles. My idea is to get the closest points and those points that match completely inside the polygons. In ArcMap there's a function for match option named CLOSEST and they have defined by:

The feature in the join features that is closest to a target feature is matched. It is possible that two or more join features are the same distance away from the target feature. When this situation occurs, one of the join features is randomly selected as the matching feature.

I have a function to intersect points into polygons. It was contributed by Lyndon Estes at the r-sig-geo list and the code works very well when all the polygons have filled all areas. The second case is known as a Spatial join distance and in ArcMap is know as INTERSECT when match_option is CLOSEST, as ArcMap does. So, you can modify the minimal distance between the point and the polygon when the area is not filled by all polygons, but I don't have access to ArcMap.

Here's the data (direct to my public ftp) and the function of the first INTERSECT:

xy.map <- readShapeSpatial("http://www.udec.cl/~jbustosm/points.shp")
manzana.map <- readShapeSpatial("http://www.udec.cl/~jbustosm/manzanas_from.shp" )

IntersectPtWithPoly <- function(x, y) { 
# Extracts values from a SpatialPolygonDataFrame with SpatialPointsDataFrame, and appends table (similar to 
# ArcGIS intersect)
# Args: 
#   x: SpatialPoints*Frame
#   y: SpatialPolygonsDataFrame
# Returns:
# SpatialPointsDataFrame with appended table of polygon attributes

  # Set up overlay with new column of join IDs in x
  z <- overlay(y, x)

  # Bind captured data to points dataframe
  x2 <- cbind(x, z)

  # Make it back into a SpatialPointsDataFrame 
  # Account for different coordinate variable names 
  if(("coords.x1" %in% colnames(x2)) & ("coords.x2" %in% colnames(x2))) {
    coordinates(x2) <- ~coords.x1 + coords.x2  
  } else if(("x" %in% colnames(x2)) & ("x" %in% colnames(x2))) {
    coordinates(x2) <- ~x + y 

  # Reassign its projection if it has one
  if(is.na(CRSargs(x@proj4string)) == "FALSE") {
    x2@proj4string <- x@proj4string  

test<-IntersectPtWithPoly (xy.map,manzana.map)

Sharing some ideas with Lyndon Estes, he told me this:

I think the easiest thing to do would be to put a buffer around each of the points (you could specify 50 m if it is in projected coordinates), converting them to polygons, and then your task becomes an intersection of two different polygon objects.

I haven't done this type of operation in R, but I suspect you could find your answer with the following functions:

library(sp) ?over

library(rgeos) ?gBuffer ?gIntersects

Maybe someone else who has a better idea on polygon to polygon intersects/overlays could suggest the method.

I know that this functions could help to achieve it.



On the other hand, I think perhaps making the points bigger could help but I'm not sure!


1 Answer 1


You could do polygon to polygon overlays using sp with rgeos. You'd need to load rgeos after you load sp.

>over(polygon1, polygon2)
  • 1
    over is provided by sp, not rgeos
    – Mike T
    Jan 12, 2012 at 1:01
  • Got it. Thanks. What I wanted to say was that you need rgeos to complete some of the over methods.
    – R.K.
    Jan 12, 2012 at 1:37
  • Hi @R.K., the output of over(polygon1, polygon2) does not bring the ids of polygon1. How can I get that information into the output ? Oct 29, 2015 at 8:41

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.