I am looking for a fast way to run a st_is_within_distance (or equivalent) in R, to find all the lines (y) within a given distance of each point (x). (This would be an equivalent of ArcGIS's "neartable", for which we do not seem to have an easy equivalent in R's GIS packages, although there are packages for points-to-points, see e.g. Is there an R equivalent of NEAR in ESRI ArcGIS?) The R help for ?st_is_within_distance contains this note:

For most predicates, a spatial index is built on argument x; see http://r-spatial.org/r/2017/06/22/spatial-index.html. Specifically, st_intersects, st_disjoint, st_touches st_crosses, st_within, st_contains, st_contains_properly, st_overlaps, st_equals, st_covers and st_covered_by all build spatial indexes for more efficient geometry calculations. st_relate, st_equals_exact, and st_is_within_distance do not.

So my interpretation is that st_is_within_distance is going to be no faster than calculating all possible pairs of distances and taking the ones within the given cutoff. From my testing, I find this is true.

So I had thought, as a way to improve on the speed of st_is_within_distance, to use an approach based on a spatial buffer: buffer the lines out to the distance I want (15km in my case), then use st_within to see if the points fall into the buffers (and pick which buffers). Conceptually this is similar to what st_is_within_distance is doing, although this gives slightly different answers, probably because the buffer is only an approximation to a "true buffer". Worrying about this turns out to be not important though, because, much to my surprise -- and despite using a spatial index -- this approach is SLOWER than calculating all possible distances! (about 1.7 times slower in my example)

Questions: (1) Is there a way (in theory and/or in R implementation) to use a spatial index to obtain a speed improvement for finding all lines within a given distance of a point? Or is calculating all possible distances the only sure way to accomplish this? (Note that in sf, st_nearest_feature uses a spatial index and is very fast, but that only obtains the closest line to the point, not all lines within a given radius of the point. This suggests to me that a spatial index could perhaps speed up the st_is_within_distance call.)

(2) Why, despite using a spatial index, is the st_within approach slower than calculating all possible distances? I am looking for any intuition to help understand why I was wrong.

  • 1
    How bout buffering the points and then use st_intersects? I am assuming you have more points than lines.
    – TimSalabim
    Sep 24 '20 at 17:09
  • @TimSalabim wow, your approach is way faster (14 seconds from buffer points and st_intersects vs 200 seconds from st_is_within_distance). yes, I'm working with variable number of points (20000 now, up to 2 million later), fixed number of lines (7000 lines). is the speed diff because st_within was preparing the polygons rather than the points when I was buffering the lines? Sep 24 '20 at 20:59
  • 1
    I believe so yes.
    – TimSalabim
    Sep 25 '20 at 5:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.