0

I have a large set of points (latitude and longitude) constrained by a complex polygon. However some of the points are not within the bounds of the polygon, and I would like to subset these points from the original data frame (not the ppp object, described below).

#Simple example of a polygon and points.
ex.poly <- data.frame(x=c(0,5,5,2.5,0), y=c(0,0,5,2.5,5))
df <- data.frame(x=c(0.5, 2.5, 4.5, 2.5), y=c(4,1,4, 4))

bound <- owin(poly=data.frame(x=ex.poly$x, y=ex.poly$y))

test.ppp <- ppp(x=df$x, y=df$y, window=bound)

#plotting example, to show the one out of the bound owin object
plot(bound)
points(df$x, df$y)

The error message 1 point was rejected as lying outside the specified window comes up, as expected. How would I subset the original dataframe df to find out which point(s) is rejected?

user3389288
  • 992
  • 9
  • 30
  • I came across [this](http://robinlovelace.net/r/2014/07/29/clipping-with-r.html) a while ago. This link may be useful. – jazzurro Oct 12 '14 at 03:43
  • You should mention that you're using the package `spatstat` . IIRC that package has tools to tell you whether a point is in the interior or not, so take a look at the help pages. – Carl Witthoft Oct 12 '14 at 12:29
  • The spatstat command to use is inside.owin –  Oct 13 '14 at 06:19

2 Answers2

4

Generate the boundary and list of points

ex.poly <- data.frame(x=c(0,5,5,2.5,0), y=c(0,0,5,2.5,5))
df <- data.frame(x=c(0.5, 2.5, 4.5, 2.5), y=c(4,1,4, 4))

Identify and plot points inside and outside the polygon using package spatstat

library(spatstat)
bound <- owin(poly=data.frame(x=ex.poly$x, y=ex.poly$y))
isin<-inside.owin(x=df$x,y=df$y,w=bound)
point_in <- df[isin,]
point_out <- df[!isin,]
plot(bound)
points(df)
points(point_out,col="red",cex = 3 )
points(point_in,col="green",cex = 3 )

spatstat in out of polygon

Alternatively, using package splancs

library(splancs)
pts<-as.matrix(df)
poly<-as.matrix(ex.poly)
point_in<-pip(pts,poly,out=FALSE)
point_out<-pip(pts,poly,out=TRUE)
plot(pts, ylim=c(0,5), xlim=c(0,5))
polygon (poly)
points(point_out,col="red",cex = 3 )
points(point_in,col="green",cex = 3 )

PointsInOutOfPolygon

0

If your main goal is to "find out which point(s) is rejected", you can simply access the rejects attribute of the test.ppp you already created:

exterior_points <- attr(test.ppp, "rejects")

exterior_points is now a separate ppp that contains only the point(s) that are outside the window you specified when creating test.ppp.

Bryan
  • 1,771
  • 4
  • 17
  • 30