I am trying to do a point-in-poly analysis using R. I have 4.5 million points and 5k polygons. This is computationally intensive. I was pointed to parallel processing through another post and tried to run point.in.poly through parallel processing.

Without parallel

table1 = point.in.poly(spdf, all_polygons)

With parallel


ptm <- proc.time()
foreach(i=1:length(split_df), .combine = rbind, .packages='spatialEco','sp')  %dopar%  {
  myvars = c("LONG_X","LAT_Y")
  producers_lat_long = split_df[i][myvars]
  spdf <- SpatialPointsDataFrame(coords = producers_lat_long, data = split_df,
                               proj4string = CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"))
  clip1 = point.in.poly(spdf, all_polygons)
proc.time() - ptm

Without parallel results in about 7 minute processing time for a 100k subset. With parallel takes longer, then freezes and crashed my computer. Note that I am on a four core computer.

How might I speed up this function? Without parallel processing, my rough estimation is that this process could take 6 hours.

  • 1
    I would try sf::st_join(pts, pols). The heavy computations are done in C++ so should be acceptable re running time. – TimSalabim Aug 2 '18 at 14:26
  • How old is your installation of r and r-spatial? As far as i know spatial indixes are pretty new to r and the sf functions. Such an index should speed this up much faster than parallel processing. Using spatial indexes this task should run in the range of seconds. – Matte Aug 2 '18 at 14:40
  • What you seem to have done there is run the whole thing nrow times across your CPUs. To get the win from foreach you have to split your task up inside the loop. For example, loop from 1:10 and then point.in.poly(spdf[a:b,], all_polygons), computing a and b from the loop index i. – Spacedman Aug 2 '18 at 15:30
  • @TimSalabim st_intersects could be quicker since it only returns the indexes and doesn't try and join the data frames? Note to macworthy - this is all using the new sf classes rather than the old SpatialWhateverDataFrame classes. – Spacedman Aug 2 '18 at 15:32
  • 1
    Make a fully reproducible example - use sample data from the R packages so we can all test this. At the moment we have no idea what getValues does, for example, so we can't tell where your problem is. – Spacedman Aug 2 '18 at 19:23

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.