Skip to content

Determining in which polygon points lie

3 messages · Dieter Vanderelst, Lionel, Roger Bivand

#
Dear List,

I'm reading in a shapefile of a country with different subregions in the 
shapefile (the shapefile was downloaded from the GADM database of Global 
Administrative Areas). I also have a number of points. Therefore, my 
code looks like this (using the maptools package):

country<-readShapePoly(file)
points<-readShapePoints(file)

Now, I want to know in which of the N regions of the country each of the 
spatial points lies? Is there an easy way to do this? I know there is a 
function points.in.poly in sp. But this functions seems incompatible 
with the objects created by the maptools functions. Also, this function 
requires looping over the regions in a set.

Regards,
Dieter

---
Behavioural, Acoustic and Sensory Ecology Lab
University Bristol, UK

http://bitsofbats.weebly.com/
#
Dear Dieter,

You can download gadm files as R SpatialPolygonsDataFrame. Then using 
load you can load them in your session.
You should keep only the column of interest for you (the N-region 
identifiers) in the SpatialPolygons, you can then use the 'over' 
function in the sp package to get for each polygons (region) the points 
ID (usually their row number, but if you have a SpatialPointsDataFrame 
the first column values).

Sincerely yours,
Lionel
On 14/06/2013 17:51, Dieter Vanderelst wrote:
#
On Fri, 14 Jun 2013, Dieter Vanderelst wrote:

            
This is a FAQ. See:

library(sp)
vignette("over")
?over

Example:

library(maptools)
nc1 <- readShapePoly(system.file("shapes/sids.shp",
  package="maptools")[1],
  proj4string=CRS("+proj=longlat +ellps=clrk66"))
dim(nc1)
set.seed(1)
pts <- spsample(nc1, n=200, type="random")
str(pts %over% as(nc1, "SpatialPolygons"))

Hope this clarifies

Roger