2

Given GIS data of polygons and two points, is there a way to find the minimum number of polygons you need to traverse to walk between the two points? Even an approximate approach would be fine for my purposes if the problem cannot always be solved exactly.

The actual application I am using this for uses tax parcels and pipelines, with the goal of finding minimum number of parcels I need to traverse to lay a pipeline between two points -- so thousands of polygons and thousands of pairs of points. But to provide an example here is some code in R using the sf package. I load North Carolina counties and put two points on it. The straight line path between the two points goes through 5 counties. But if you look at the map you can see there's a path that goes through only 4 counties. I want an approach to calculate that 4 in this case, and more generally.

# install.packages(c('sf','dplyr')) # if needed library(sf) library(dplyr) # use NC counties polygons polys <- read_sf(system.file("shape/nc.shp", package = "sf")) # make two points, goal is to find the minimum # of counties to get between them points <- data.frame(name=c('start','end'),LAT=c(34.623858,35.257017),LONG=c(-78.851954,-78.237195)) points <- st_as_sf(points, coords = c("LONG", "LAT"), crs = 4326) plot(polys$geometry) plot(points$geometry, add=TRUE,col='blue') plot(st_cast(summarize(points,name=first(name)),'LINESTRING'),add=TRUE,col='blue') 
3
  • 2
    It sounds like a minimum spanning tree community.esri.com/thread/… is what you're after, the link is an Esri thread but methods are discussed and evaluated, with a bit of googling you might be able to find the solution you're after. Commented Dec 4, 2019 at 22:12
  • 1
    Yes, minimum spanning tree. QGIS has a pluging for that plugins.qgis.org/plugins/minimum_spanning_tree You would need to convert your parcel centers to points and connect the points with lines which would be easy to do in ArcGIS by creating a TIN and exporting the edges to lines. Commented Dec 4, 2019 at 22:37
  • It looks like you are using R so I've added a tag for that to assist your question reaching the filters of experts in that. Commented Dec 4, 2019 at 22:53

1 Answer 1

2

Helpful comments to my question led me to study minimum spanning trees, which aren't exactly the solution -- spanning trees connect all nodes in the minimum cost way, while I want to connect two given nodes in the minimum cost way. But the insight to consider the data as a graph carries over and leads to a simple solution. Once the GIS data is converted to an adjacency graph there are ready algorithms to find the shortest path between two polygons along the graph. (For example Dijkstra's algorithm, although for my particular problem which has no cost weights it appears that R's igraph package uses a faster algorithm.)

Steps for my solution are:

  1. take the polygon dataset and convert it to a graph using adjacency as the definition of edge, and polygon as the definition of vertex.

  2. Given the two points (start and end) find the vertices that contain them.

  3. Then apply a shortest path algorithm for the two containing vertices.

R code and packages implementing this for the North Carolina example below.

# install.packages(c('sf','dplyr','spdep','igraph')) # if needed library(sf) library(dplyr) library(spdep) library(igraph) # use NC counties polygons for an example polys <- read_sf(system.file("shape/nc.shp", package = "sf")) # pick two points, goal is to find the minimum # of counties to get between them # for these two points, the shortest path between the points (straight line) actually crosses 5 counties, but there are paths of length 4 counties points <- data.frame(name=c('start','end'),LAT=c(34.623858,35.257017),LONG=c(-78.851954,-78.237195)) points <- st_as_sf(points, coords = c("LONG", "LAT"), crs = 4267) # show the problem: plot(polys$geometry) plot(points$geometry, add=TRUE,col='blue') plot(st_cast(summarize(points,name=first(name)),'LINESTRING'),add=TRUE,col='blue') # find the indices of the polygons of the two points indices <- as.integer(st_intersects(points,polys)) indices # points are in county index 94 for start, and 62 for end # convert to graph of polygons (polygons = vertices; edges = adjacency) nc.nb <- spdep::poly2nb(polys) g <- igraph::graph_from_adj_list(nc.nb) plot(polys$geometry) plot(polys[94,]$geometry,add=TRUE,col='red') # starting county plot(polys[62,]$geometry,add=TRUE,col='red') # ending county path <- igraph::shortest_paths(g, from = 94, to = 62) # find shortest path pathint <- as.integer(path$vpath[[1]]) plot(polys[pathint,]$geometry,add=TRUE,col='red') # here's a shortest path plotted in red length(pathint) # gives 4 -- the right answer 

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.