Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How to exclude points beyond a line in R

I'm working with two dataframes in R: a "red" dataframe and a "black" dataframe. In both there are two columns representing the coordinates.
I used a plot to explain what I want to do.

I would like to select all the points from the "red" dataframe that are beyond the "black" line. E.g. all the points excluded from the area of the polygon delimited by the black points.

enter image description here

like image 399
Giordano Avatar asked Sep 01 '25 11:09

Giordano


2 Answers

My previous answer only showed how not to draw the points outside the polygon. To actually identify the points outside the polygon, you can use the function pip2d from the package ptinpoly. It returns negative values for points outside the polygon.

Example:

library(ptinpoly)
poly.vertices <- data.frame(x=c(20,40,80,50,40,30), y=c(30,20,70,80,50,60))
p <- data.frame(x=runif(100, min=0, max=100), y=runif(100, min=0, max=100))

outside <- (pip2d(as.matrix(poly.vertices), as.matrix(p)) < 0)

plot(p$x, p$y, col=ifelse(outside, "red", "black"))
polygon(poly.vertices$x, poly.vertices$y, border="blue", col=NA)

The same should be achieved with the function PtInPoly from the package DescTools, which returns zero for points outside the polygon. The implementation of ptinpoly, however, has the advantage of implmenting the particularly efficient algorithm described in

J. Liu, Y.Q. Chen, J.M. Maisog, G. Luta: "A new point containment test algorithm based on preprocessing and determining triangles." Computer-Aided Design, Volume 42, Issue 12, December 2010, Pages 1143-1150

Edit: Out of curiosity, I have compared the runtime of ptinpoly::pip2d and DescTools::PtInPoly with microbenchmark and N=50000 points, and pip2d turned out to be considerably faster:

> microbenchmark(outside.pip2d(), outside.PtInPoly())
Unit: milliseconds
              expr       min        lq      mean   median        uq      max
   outside.pip2d()  3.375084  3.421631  4.459051  3.48939  4.251395 65.97793
outside.PtInPoly() 27.537927 27.666688 28.739288 27.97984 28.514595 90.11313

neval 100 100

enter image description here

like image 142
cdalitz Avatar answered Sep 02 '25 23:09

cdalitz


You could use the sf package to define a convex hull and intersect your target points with that polygon.

Define a convex hull based on black:

library(sf)

set.seed(99)
red <- data.frame(x = runif(100,-10,10), y = runif(100,-4,4))
black <- data.frame(x = runif(100,-8,8), y = runif(100,-4,3))

# Convert df to point feature
blk <- st_as_sf(black, coords = c("x", "y"))
# Convert to multipoint
blk_mp <- st_combine(blk)
# Define convex hull
blk_poly <- st_convex_hull(blk_mp)

plot(black)
points(red, col = "red")
plot(blk_poly, add = TRUE)

enter image description here

Intersecting red with the convex hull returns red within that polygon:

rd <- st_as_sf(red, coords = c("x", "y"))
rd_inside <- st_intersection(rd, blk_poly)

plot(black)
points(red)
plot(blk_poly, add = TRUE)
plot(rd_inside, pch = 24, col = "red", bg = "red", add = TRUE)

enter image description here

like image 29
Skaqqs Avatar answered Sep 02 '25 23:09

Skaqqs