0

I am trying to join spatial point data to spatial polygon data.

A while back, this very helpful thread was posted: How to pick up the information for the nearest associated polygon to points using R?

and cross-linked: https://gis.stackexchange.com/questions/18726/spatial-distance-join-between-polygon-and-points-using-r

I'm trying to use the code posted by the OP:

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

ky.map <- readShapeSpatial("~/KYcounties.shp") 

points.xy <- readShapeSpatial('~/points.shp') # Sample data of the points file 
                                              # is below.

The function I am using is

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)
}

and is called as

test <- IntersectPtWithPoly(points.xy, ky.map)

Which works pretty well, except I get:

Warning message: 'overlay' is deprecated. Use 'over' instead.

As well as a bunch of NAs for the joined polygon fields.

Fair enough:

z <- over(y, x) 

But that gets me:

Error in data.frame(..., check.names = FALSE) : arguments imply differing
number of rows: 96, 121

Below is a sample of the point data I am using (dropped into QGIS & exported as a shapefile):

structure(list(ID = c(5L, 11L, 33L, 41L, 73L, 119L, 175L, 206L, 
216L, 229L), Latitude = c(38.239829, 38.9231075, 38.04169866, 
38.13953849, 38.14800534, 37.85059283, 38.17820492, 38.00833019, 
38.22189947, 36.87255699), Longitude = c(-85.746211, -84.4108925, 
-84.49738797, -85.6776972, -85.02502406, -84.57452618, -85.56667237, 
-84.42084641, -85.48653626, -84.22021486)), .Names = c("ID", 
"Latitude", "Longitude"), row.names = c(NA, 10L), class = "data.frame")

I've shared the shapefile on Google Drive: https://docs.google.com/file/d/0B1bNxpx4XS8dSFJlOTlQSUw0LTg/edit?usp=sharing

Thoughts?

Community
  • 1
  • 1
NiuBiBang
  • 628
  • 1
  • 15
  • 30
  • Please add example data to your question, without it it is impossible to reproduce your problem. – Paul Hiemstra Aug 13 '13 at 18:51
  • Might be possible with `spatstat` package `dist` or `crossdist` functions. – Carl Witthoft Aug 13 '13 at 19:44
  • @CarlWitthoft it's perfectly possible with `sp` classes and `rgeos`, but we need more info from the OP to know what is going on. – Simon O'Hanlon Aug 13 '13 at 20:18
  • You haven't added the download link for the shapefile. Also, to notify someone you have left a comment you have to include their handle in the comment ( e.g. @SimonO101 ) – Simon O'Hanlon Aug 14 '13 at 09:32
  • @SimonO101 I've added a google drive link for the shapefile. Thanks for the tutorial on @ handles -- this is my first time using the internet. Ever. – NiuBiBang Aug 14 '13 at 13:46
  • @DaNiu **1)** Your now deleted comment in reply to someone else did not include a handle so they were not notified that you had updated your question. I was trying to help you out. **2)** For a good example of a completely self-contained reproducible example that doesn't require potential answerers to do a lot of messing around see [**my question here**](http://stackoverflow.com/q/11179666/1478381) **3)** Your sarcasm has completely evaporated my willingness to help you. Ever. – Simon O'Hanlon Aug 14 '13 at 13:57
  • @SimonO101 1) i didn't leave any comment. all i did was update the post. 2) thnx for the link, I will try to reproduce similarly shortly. 3) my sarcasm is an acquired taste. – NiuBiBang Aug 14 '13 at 14:19
  • Possible duplicate of [arguments imply differing number of rows: 8, 20](http://stackoverflow.com/questions/26147558/arguments-imply-differing-number-of-rows-8-20) – rafa.pereira Oct 26 '15 at 23:04

1 Answers1

0

There is a simple solution using the spatialEco library.

library(spatialEco)

# intersect points in polygon
  pts <- point.in.poly(pts, ply)

# check plot
  plot(ply)
  plot(a, add=T)

# convert to data frame, keeping your data
  pts <- as.data.frame(pts)
rafa.pereira
  • 13,251
  • 6
  • 71
  • 109