[GIS] Spatial distance join between polygon and points using R

attribute-joinsrshapefile

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 ArcGIS 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 kindly 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 ArcGIS is know as INTERSECT when match_option is CLOSEST, as ArcGIS 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 ArcGIS.

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

library(rgeos)
library(sp) 
library(maptools)
library(rgdal)
library(sp)
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  
  }
  return(x2)
}


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.

library(sp)
?over

library(rgeos)
?gBuffer
?gIntersects

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

Best Answer

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

>library(rgeos)
>over(polygon1, polygon2)