st_network_paths | R Documentation |
Find shortest paths between nodes in a spatial network
st_network_paths(
x,
from,
to = node_ids(x),
weights = edge_length(),
all = FALSE,
k = 1,
direction = "out",
router = getOption("sfn_default_router", "igraph"),
use_names = FALSE,
return_cost = TRUE,
return_geometry = TRUE,
...
)
x |
An object of class |
from |
The node where the paths should start. Evaluated by
|
to |
The nodes where the paths should end. Evaluated by
|
weights |
The edge weights to be used in the shortest path calculation.
Evaluated by |
all |
Should all shortest paths be returned for each pair of nodes? If
set to |
k |
The number of paths to find. Setting this to any integer higher
than 1 returns not only the shortest path, but also the next k - 1 loopless
shortest paths, which may be longer than the shortest path. Currently, this
is only supported for one-to-one routing, meaning that both the from and to
argument should be of length 1. This argument is ignored if |
direction |
The direction of travel. Defaults to |
router |
The routing backend to use for the shortest path computation.
Currently supported options are |
use_names |
If a column named |
return_cost |
Should the total cost of each path be computed? Defaults
to |
return_geometry |
Should a linestring geometry be constructed for each
path? Defaults to |
... |
Additional arguments passed on to the underlying function of the chosen routing backend. See Details. |
The sfnetworks package does not implement its own routing algorithms to find shortest paths. Instead, it relies on "routing backends", i.e. other R packages that have implemented such algorithms. Currently two different routing backends are supported.
The default is igraph
. This package supports
one-to-many shortest path calculation with the
shortest_paths
function. Note that multiple from nodes
are not supported. If multiple from nodes are given, only the first one is
taken. The igraph router also supports the computation of all shortest path
(see the all
argument) through the
all_shortest_paths
function and of k shortest paths
(see the k
argument) through the
k_shortest_paths
function. In the latter case, only
one-to-one routing is supported, meaning that also only one to node should
be provided. The igraph router does not support dual-weighted routing.
The second supported routing backend is dodgr
. This
package supports many-to-many shortest path calculation with the
dodgr_paths
function. It also supports dual-weighted
routing. The computation of all shortest paths and k shortest paths is
currently not supported by the dodgr router. The dodgr package is a
conditional dependency of sfnetworks. Using the dodgr router requires the
dodgr package to be installed.
The default router can be changed by setting the sfn_default_router
option.
An object of class sf
with one row per requested
path. If return_geometry = FALSE
or edges are spatially implicit, a
tbl_df
is returned instead. If a requested path could
not be found, it is included in the output as an empty path.
Depending on the argument settings, the output may include the following columns:
from
: The index of the node at the start of the path.
to
: The index of the node at the end of the path.
node_path
: A vector containing the indices of all nodes on
the path, in order of visit.
edge_path
: A vector containing the indices of all edges on
the path, in order of visit.
path_found
: A boolean describing if the requested path exists.
cost
: The total cost of the path, obtained by summing the
weights of all visited edges. Included if return_cost = TRUE
.
geometry
: The geometry of the path, obtained by merging the
geometries of all visited edges. Included if return_geometry = TRUE
and the network has spatially explicit edges.
st_network_cost
, st_network_travel
library(sf, quietly = TRUE)
library(tidygraph, quietly = TRUE)
oldpar = par(no.readonly = TRUE)
par(mar = c(1,1,1,1))
net = as_sfnetwork(roxel, directed = FALSE) |>
st_transform(3035)
# Compute the shortest path between two nodes.
# Note that geographic edge length is used as edge weights by default.
paths = st_network_paths(net, from = 495, to = 121)
paths
plot(net, col = "grey")
plot(st_geometry(net)[paths$from], pch = 20, cex = 2, add = TRUE)
plot(st_geometry(paths), col = "orange", lwd = 3, add = TRUE)
# Compute the shortest paths from one to multiple nodes.
# This will return a tibble with one row per path.
paths = st_network_paths(net, from = 495, to = c(121, 131, 141))
paths
plot(net, col = "grey")
plot(st_geometry(net)[paths$from], pch = 20, cex = 2, add = TRUE)
plot(st_geometry(paths), col = "orange", lwd = 3, add = TRUE)
# Compute the shortest path between two spatial point features.
# These are snapped to their nearest node before finding the path.
p1 = st_geometry(net, "nodes")[495] + st_sfc(st_point(c(50, -50)))
st_crs(p1) = st_crs(net)
p2 = st_geometry(net, "nodes")[121] + st_sfc(st_point(c(-10, 100)))
st_crs(p2) = st_crs(net)
paths = st_network_paths(net, from = p1, to = p2)
paths
plot(net, col = "grey")
plot(c(p1, p2), pch = 20, cex = 2, add = TRUE)
plot(st_geometry(net)[paths$from], pch = 4, cex = 2, add = TRUE)
plot(st_geometry(paths), col = "orange", lwd = 3, add = TRUE)
# Use a node type query function to specify destinations.
st_network_paths(net, 1, node_is_adjacent(1))
# Use a spatial edge measure to specify edge weights.
# By default edge_length() is used.
st_network_paths(net, p1, p2, weights = edge_displacement())
# Use a column in the edges table to specify edge weights.
# This uses tidy evaluation.
net |>
activate("edges") |>
mutate(foo = runif(n(), min = 0, max = 1)) |>
st_network_paths(p1, p2, weights = foo)
# Compute the shortest paths without edge weights.
# This is the path with the fewest number of edges, ignoring space.
st_network_paths(net, p1, p2, weights = NA)
# Use the dodgr router for many-to-many routing.
paths = st_network_paths(net,
from = c(1, 2),
to = c(10, 11),
router = "dodgr"
)
# Use the dodgr router for dual-weighted routing.
paths = st_network_paths(net,
from = c(1, 2),
to = c(10, 11),
weights = dual_weights(edge_segment_count(), edge_length()),
router = "dodgr"
)
par(oldpar)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.