Distance from point to nearest polygon in R

Question

I am working on a project at the moment where I have a point function - a point function includes 142 points - and a polygon (about 10). I want to calculate the distance between each point and the closest polygon function in R.

My current approach is tedious and a bit long. I am currently planning to calculate the distance between each point and each individual polygon. For example, I would calculate the distance between 142 points and polygon A, the distance between 142 points and polygon B, the distance between 142 points and polygon C, etc. Here is a sample code for one of these distance calculations:

dist_cen_polya <- dist2Line(centroids_coor, polygonA_shp) 

After performing these calculations, I would write code to select the minimum / closest distance between each individual point and the nearest polygon. The problem is that this procedure is tedious.

Does anyone know a package / code that minimizes runtime / computation? Would I really like to use a package that compares one to point to the closest polygon function, or calculates the distance between a point and all the polygons of interest?

Thanks.

Answer 1

Here I use the gDistance function in the rgeos topology library. I use double brute force loop, but it is surprisingly fast. It takes less than 2 seconds for 142 points and 10 polygons. I am sure there is a more elgetic way to execute a loop.

  require(rgeos) # CREATE SOME DATA USING meuse DATASET data(meuse) coordinates(meuse) <- ~x+y pts <- meuse[sample(1:dim(meuse)[1],142),] data(meuse.grid) coordinates(meuse.grid) = c("x", "y") gridded(meuse.grid) <- TRUE meuse.grid[["idist"]] = 1 - meuse.grid[["dist"]] polys <- as(meuse.grid, "SpatialPolygonsDataFrame") polys <- polys[sample(1:dim(polys)[1],10),] plot(polys) plot(pts,pch=19,cex=1.25,add=TRUE) # LOOP USING gDistance, DISTANCES STORED IN LIST OBJECT Fdist <- list() for(i in 1:dim(pts)[1]) { pDist <- vector() for(j in 1:dim(polys)[1]) { pDist <- append(pDist, gDistance(pts[i,],polys[j,])) } Fdist[[i]] <- pDist } # RETURN POLYGON (NUMBER) WITH THE SMALLEST DISTANCE FOR EACH POINT ( min.dist <- unlist(lapply(Fdist, FUN=function(x) which(x == min(x))[1])) ) # RETURN DISTANCE TO NEAREST POLYGON ( PolyDist <- unlist(lapply(Fdist, FUN=function(x) min(x)[1])) ) # CREATE POLYGON-ID AND MINIMUM DISTANCE COLUMNS IN POINT FEATURE CLASS pts@data <- data.frame( pts@data , PolyID=min.dist, PDist=PolyDist) # PLOT RESULTS require(classInt) ( cuts <- classIntervals( pts@data $PDist, 10, style="quantile") ) plotclr <- colorRampPalette(c("cyan", "yellow", "red"))( 20 ) colcode <- findColours(cuts, plotclr) plot(polys,col="black") plot(pts, col=colcode, pch=19, add=TRUE) 

The vector min.dist represents the line number of the polygon. For example, you can multiply the nearest polygons using this vector as such.

 near.polys <- polys[unique(min.dist),] 

The PolyDist vector contains the actual Cartesian minimum distances in projection feature units.

Answer 2

In the landfill you have quite a few lines. The distance between the polygon is zero if the point lies inside the polygon or lies on an edge.

So you are looking for two cases:

• Check if the point is inside any polygon (remember that it can be more than one polygon
• Get a collection of all the edges and calculate each point distance from the edge. The closest distance gives you the distance from the polygon to which the edge belongs.

So, this is a simple algorithm that, if we look at 10 edges per polygon, takes O (10 * 10) * 142 for all of your points. This makes 100 * 142 = 14,200 operations. => O (m * deltaE) * n (m is the number of polygons, deltaE is the average number of edges per polygon, n is the number of points).

Now we are checking to see if we can speed it up. The first thing that comes to mind is that we can use bounding boxes or bounding circles for each polygon.

Another idea is to prepare the nearest edges for each polygon for a set of angles. For example, if you have 8 angles (every 45 °), you can remove all edges from the list, which will be replaced by another edge (therefore, any point on the far edge will always be a greater distance than any point on the other edge of the same polygon.

This way, as a rule, you can significantly reduce complexity. If you think of a rectangle, you only have one or two edges instead of 4. If you think of a regular 8x polygon, you can end up with one or two and a maximum of 3 edges per polygon.

If you add a normal vector to each edge, you can calculate whether the point can be inside, and you need to perform a scan or any check (or you know its konvex).

There are also mapping indices, such as dividing two-dimensional space into x and y in an equidistant mannor. In this case, you only need to check the polygons lying in 9 sectors.

In the next version, an R-tree can be used so that each bounding box (circle) of each node should be checked for the minimum expected distance and maximum expected distance. Therefore, you do not need to check node polygons that result in a much larger minimum distance than the maximum distances of another node.

Another thing is if you have a tree like yours, for example, with map data. On the street map you always have the world → region → country → county → city → urban sector → ...

In this case, you can find the closest location on the world map containing millions of polygons in a reasonable amount of time, mostly <10ms.

So you have a lot of options here. And pre-processing the list of polygons and extracting the corresponding edges either using binary spatial partitions of the polygon trees, or using the angular approach, or even with something even more attractive. It's up to you to decide. I expect you to end up doing something in the logarithmic range, like O (log (n) * log (deltaE)), becoming O (log (n)) as medium complexity.