3

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 

library(foreach)
library(doParallel)


registerDoParallel(makeCluster(2))
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)
  getValues(clip1)
}
proc.time() - ptm
endCluster()

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.

Improve this question
10
  • 1
    I would try sf::st_join(pts, pols). The heavy computations are done in C++ so should be acceptable re running time.
    – TimSalabim
    Commented Aug 2, 2018 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
    Commented Aug 2, 2018 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
    Commented Aug 2, 2018 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
    Commented Aug 2, 2018 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
    Commented Aug 2, 2018 at 19:23

0

Reset to default

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.