Как преобразовать список пространственных точек SF в маршрутизируемый граф
У меня есть sf dataframe
объект с серией точек, представляющих форму автобусного маршрута. Я хотел бы превратить этот объект в маршрутизируемый граф, чтобы я мог оценить время, необходимое для перехода от точкиc
к t
.
Вот что я пробовал использовать dodgr
package, но я не уверен, что делаю здесь неправильно:
library(dodgr)
graph <- weight_streetnet(mydata, wt_profile = "motorcar", type_col="highway" , id_col = "id")
Ошибка в check_highway_osmid(x, wt_profile): укажите type_col, который будет использоваться для взвешивания streetnet
Воспроизводимые данные
Данные выглядят как на изображении ниже
mydata <- structure(list(shape_id = c(52421L, 52421L, 52421L, 52421L, 52421L,
52421L, 52421L, 52421L, 52421L, 52421L, 52421L, 52421L, 52421L,
52421L, 52421L, 52421L, 52421L, 52421L, 52421L, 52421L), length = structure(c(0.191422504106197,
0.191422504106197, 0.191422504106197, 0.191422504106197, 0.191422504106197,
0.191422504106197, 0.191422504106197, 0.191422504106197, 0.191422504106197,
0.191422504106197, 0.191422504106197, 0.191422504106197, 0.191422504106197,
0.191422504106197, 0.191422504106197, 0.191422504106197, 0.191422504106197,
0.191422504106197, 0.191422504106197, 0.191422504106197), units = structure(list(
numerator = "km", denominator = character(0)), class = "symbolic_units"), class = "units"),
geometry = structure(list(structure(c(-46.5623281998182,
-23.5213458001468), class = c("XY", "POINT", "sfg")), structure(c(-46.562221,
-23.52129), class = c("XY", "POINT", "sfg")), structure(c(-46.562121,
-23.521235), class = c("XY", "POINT", "sfg")), structure(c(-46.5620233332577,
-23.5211840000609), class = c("XY", "POINT", "sfg")), structure(c(-46.561925666591,
-23.5211330000609), class = c("XY", "POINT", "sfg")), structure(c(-46.561828,
-23.521082), class = c("XY", "POINT", "sfg")), structure(c(-46.5618098335317,
-23.5212126666783), class = c("XY", "POINT", "sfg")), structure(c(-46.5617916670273,
-23.5213433333544), class = c("XY", "POINT", "sfg")), structure(c(-46.5617735004869,
-23.5214740000284), class = c("XY", "POINT", "sfg")), structure(c(-46.5617553339104,
-23.5216046667004), class = c("XY", "POINT", "sfg")), structure(c(-46.5617371672978,
-23.5217353333702), class = c("XY", "POINT", "sfg")), structure(c(-46.5617190006492,
-23.5218660000379), class = c("XY", "POINT", "sfg")), structure(c(-46.5617008339645,
-23.5219966667036), class = c("XY", "POINT", "sfg")), structure(c(-46.5616826672438,
-23.5221273333671), class = c("XY", "POINT", "sfg")), structure(c(-46.5616645004869,
-23.5222580000284), class = c("XY", "POINT", "sfg")), structure(c(-46.5616463336941,
-23.5223886666877), class = c("XY", "POINT", "sfg")), structure(c(-46.5616281668651,
-23.5225193333449), class = c("XY", "POINT", "sfg")), structure(c(-46.56161,
-23.52265), class = c("XY", "POINT", "sfg")), structure(c(-46.5617355000207,
-23.5226427501509), class = c("XY", "POINT", "sfg")), structure(c(-46.5618610000276,
-23.5226355002012), class = c("XY", "POINT", "sfg"))), class = c("sfc_POINT",
"sfc"), precision = 0, bbox = structure(c(xmin = -46.5623281998182,
ymin = -23.52265, xmax = -46.56161, ymax = -23.521082), class = "bbox"), crs = structure(list(
epsg = 4326L, proj4string = "+proj=longlat +datum=WGS84 +no_defs"), class = "crs"), n_empty = 0L),
id = c("a", "b", "c", "d", "e", "f", "g", "h", "i", "j",
"k", "l", "m", "n", "o", "p", "q", "r", "s", "t"), speed_kmh = c(11,
11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11,
11, 11, 11, 11)), sf_column = "geometry", agr = structure(c(shape_id = NA_integer_,
length = NA_integer_, id = NA_integer_, speed_kmh = NA_integer_
), class = "factor", .Label = c("constant", "aggregate", "identity"
)), row.names = c("1.13", "1.14", "1.15", "1.16", "1.17", "1.18",
"1.19", "1.20", "1.21", "1.22", "1.23", "1.24", "1.25", "1.26",
"1.27", "1.28", "1.29", "1.30", "1.31", "1.32"), class = c("sf",
"data.table", "data.frame"))
3 ответа
В weight_streetnet
функция действительно предназначена только для обработки реальных уличных сетей, как правило, созданных osmdata::osmdata_sf/sp/sc()
функции. Тем не менее, его можно настроить для обработки подобных случаев. Главное, что нужно, - это преобразовать точки во что-то, что знает о ребрах между ними, напримерsf::LINESTRING
объект:
x <- sf::st_combine (mydata) %>%
sf::st_cast ("LINESTRING") %>%
sf::st_sf ()
Это дает однорядный объект, который затем может быть преобразован в dodgr
формат, а id
значения сопоставлены с краями
net <- weight_streetnet (x, type_col = "shape_id", id_col = "id", wt_profile = 1)
net$from_id <- mydata$id [as.integer (net$from_id)]
net$to_id <- mydata$id [as.integer (net$to_id)]
В таком случае, dodgr
будут рассчитаны и вставлены расстояния непосредственно из географических координат. Ваши расстояния также могут быть вставлены и использованы для трассировки, заменивd_weighted
ценности:
net$d_weighted <- as.numeric (mydata$length [1])
dodgr_dists (net, from = "c", to = "t") # 236.0481
Если вы действительно хотите, чтобы ваши расстояния представляли абсолютные расстояния, используемые для расчета окончательного результата, просто замените $d
ценности
net$d <- net$d_weighted
dodgr_dists (net, from = "c", to = "t") # 3.254183
Обратите внимание, что для "простых" задач, подобных этой, igraph
обычно будет быстрее, потому что он вычисляет маршруты с использованием единственного набора весов. Единственное реальное преимуществоdodgr
в этом контексте - возможность использовать "двойные веса" - $d_weighted
а также $d
значения - такие, что маршрут рассчитывается согласно $d_weighted
, а конечные расстояния согласно $d
.
Если вы хотите включить это в "аккуратный" рабочий процесс, вы также можете рассмотреть возможность сочетания sf
а также tidygraph
. Последний предлагает удобную структуру для сетей / графиков в видеtbl_graph
класс, подклассы которого igraph
(следовательно, вы можете использовать tbl_graph
объекты внутри всех igraph
функционирует как igraph
объект). Однако вы можете анализировать свои узлы и ребра как тибблы и использовать функции какfilter()
, select()
, mutate()
и так далее. Конечно, эти таблицы могут также содержать столбец списка геометрии, который мы знаем изsf
, добавляя географическую информацию к узлам и ребрам.
Подход еще далек от совершенства, и улучшения были бы очень желательны, но все же он показывает другой способ решения проблемы.
# Load libraries.
library(tidyverse)
library(sf)
library(tidygraph)
library(igraph)
library(units)
Как и в других ответах, нам нужно создать ребра между узлами. Пока я предполагаю, что точки просто связаны в алфавитном порядке. Дляtidygraph
Однако, похоже, нам нужны числовые идентификаторы, а не символы.
# Add a numeric ID column to the nodes.
nodes <- mydata %>%
rename(id_chr = id) %>%
rowid_to_column("id") %>%
select(id, id_chr, everything())
# Define the source node of each edge, and the target node of each edge.
sources <- nodes %>% slice(-n())
targets <- nodes %>% slice(-1)
# Write a function to create lines between data frames of source and target points.
pt2l <- function(x, y) { st_linestring(rbind(st_coordinates(x), st_coordinates(y))) }
# Create the edges.
edges <- tibble(
from = sources %>% pull(id),
to = targets %>% pull(id),
length = sources %>% pull(length),
speed = sources %>% pull(speed_kmh),
geometry = map2(st_geometry(sources), st_geometry(targets), pt2l)
) %>% st_as_sf() %>% st_set_crs(st_crs(nodes))
# Add a time column to the edges.
edges <- edges %>%
mutate(speed = set_units(speed, "km/h")) %>%
mutate(time = length / speed)
# Clean up the nodes data.
nodes <- nodes %>%
select(-length, -speed_kmh)
# Create the tbl_graph object out of the nodes and edges.
# Providing the edges as sf object is problematic for tidygraph, unfortunately.
# Therefore, we have to provide them as a tibble.
graph <- tbl_graph(nodes = nodes, edges = as_tibble(edges), directed = FALSE)
Это дает нам следующее tbl_graph
объект:
# A tbl_graph: 20 nodes and 19 edges
#
# An undirected simple graph with 1 component
#
# Node Data: 20 x 4 (active)
id id_chr shape_id geometry
<int> <chr> <int> <POINT [°]>
1 1 a 52421 (-46.56233 -23.52135)
2 2 b 52421 (-46.56222 -23.52129)
3 3 c 52421 (-46.56212 -23.52124)
4 4 d 52421 (-46.56202 -23.52118)
5 5 e 52421 (-46.56193 -23.52113)
6 6 f 52421 (-46.56183 -23.52108)
# … with 14 more rows
#
# Edge Data: 19 x 6
from to length speed geometry time
<int> <int> [km] [km/h] <LINESTRING [°]> [h]
1 1 2 0.1914225 11 (-46.56233 -23.52135, -46.56222 -23.5… 0.017402…
2 2 3 0.1914225 11 (-46.56222 -23.52129, -46.56212 -23.5… 0.017402…
3 3 4 0.1914225 11 (-46.56212 -23.52124, -46.56202 -23.5… 0.017402…
# … with 16 more rows
Теперь у нас есть все в структуре графа, мы можем выбрать узел, из которого мы хотим выполнить маршрут, и узел, к которому мы хотим выполнить маршрут, и найти кратчайший путь между ними со временем в пути в качестве переменной веса, используя shortest_path
функция от igraph
. Теперь мы просто работаем с маршрутом "один-к-одному" (от "c" до "t"), но он будет таким же для "один ко многим", "многие к одному" или "многие ко многим".
# Select the node from which and to which the shortest path should be found.
from_node <- graph %>%
activate(nodes) %>%
filter(id_chr == "c") %>%
pull(id)
to_node <- graph %>%
activate(nodes) %>%
filter(id_chr == "t") %>%
pull(id)
# Find the shortest path between these nodes
path <- shortest_paths(
graph = graph,
from = from_node,
to = to_node,
output = 'both',
weights = graph %>% activate(edges) %>% pull(time)
)
Результирующий путь - это список с узлами и ребрами, составляющими путь.
$vpath
$vpath[[1]]
+ 18/20 vertices, from e43a089:
[1] 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
$epath
$epath[[1]]
+ 17/19 edges from e43a089:
[1] 3-- 4 4-- 5 5-- 6 6-- 7 7-- 8 8-- 9 9--10 10--11 11--12 12--13
[11] 13--14 14--15 15--16 16--17 17--18 18--19 19--20
Мы можем создать подграф исходного графа, который содержит только узлы и ребра кратчайшего пути.
path_graph <- graph %>%
subgraph.edges(eids = path$epath %>% unlist()) %>%
as_tbl_graph()
# A tbl_graph: 18 nodes and 17 edges
#
# An undirected simple graph with 1 component
#
# Node Data: 18 x 4 (active)
id id_chr shape_id geometry
<int> <chr> <int> <POINT [°]>
1 3 c 52421 (-46.56212 -23.52124)
2 4 d 52421 (-46.56202 -23.52118)
3 5 e 52421 (-46.56193 -23.52113)
4 6 f 52421 (-46.56183 -23.52108)
5 7 g 52421 (-46.56181 -23.52121)
6 8 h 52421 (-46.56179 -23.52134)
# … with 12 more rows
#
# Edge Data: 17 x 6
from to length speed geometry time
<int> <int> [km] [km/h] <LINESTRING [°]> [h]
1 1 2 0.1914225 11 (-46.56212 -23.52124, -46.56202 -23.5… 0.017402…
2 2 3 0.1914225 11 (-46.56202 -23.52118, -46.56193 -23.5… 0.017402…
3 3 4 0.1914225 11 (-46.56193 -23.52113, -46.56183 -23.5… 0.017402…
# … with 14 more rows
Здесь происходит то, что мне не нравится. Tidygraph/igraph, похоже, имеет внутреннюю структуру идентификатора узла, и вы видите, что в подграфеfrom
а также to
столбцы в данных egdes не совпадают с нашими id
столбца в данных узлов больше, а вместо этого просто ссылаются на номера строк данных узлов. Я не знаю, как это исправить.
В любом случае, теперь у нас есть путь от "c" к "t" в качестве подграфа, и мы можем легко его проанализировать. Например, подсчитав общее время прохождения пути (как был вопрос).
path_graph %>%
activate(edges) %>%
as_tibble() %>%
summarise(total_time = sum(time))
# A tibble: 1 x 1
total_time
[h]
1 0.2958348
Но его также легко построить с сохранением географической информации (просто экспортируя узлы и ребра как объекты SF).
ggplot() +
geom_sf(data = graph %>% activate(edges) %>% as_tibble() %>% st_as_sf(), col = 'darkgrey') +
geom_sf(data = graph %>% activate(nodes) %>% as_tibble() %>% st_as_sf(), col = 'darkgrey', size = 0.5) +
geom_sf(data = path_graph %>% activate(edges) %>% as_tibble() %>% st_as_sf(), lwd = 1, col = 'firebrick') +
geom_sf(data = path_graph %>% activate(nodes) %>% filter(id %in% c(from_node, to_node)) %>% as_tibble() %>% st_as_sf(), size = 2)
Об этом подходе tidygraph-sf может появиться сообщение в блоге r-space;)
Я думаю, вы можете решить эту проблему, преобразовав ваши данные в объект igraph и используя функции библиотеки igraph. Вам необходимо установить Edges и Vertex, а также значения веса. В igraph Edge - это ссылка, представляющая соединение между двумя узлами (Source и Target). В этом случае ссылка - это "улица", а точки - это узлы.
library(igraph)
GraphResult <- data.frame(Source = c(NULL),
Target = c(NULL),
weight = c(NULL))
for (i in 1:(dim(mydata)[1] - 1)) {
TempGraphResult <- data.frame(Source = c(0),
Target = c(0),
weight = c(0))
TempGraphResult$Source[1] <- mydata$id[i]
TempGraphResult$Target[1] <- mydata$id[i + 1]
TempGraphResult$weight[1] <- mydata$length[i]
GraphResult <- rbind(GraphResult, TempGraphResult) }
MyIgraph <- graph_from_data_frame(GraphResult)
#In this case works perfectly. But if you have more weight variables and even
#additional variables for the nodes, igraph have functions for constructing the
#igraph object
distances(MyIgraph, "c", "t") #returns 3.254183. Seems correct (0.1914225*17)
SquareMatrix <- distances(MyIgraph)
#*distances() is a function from igraph that performs the routing calculations.
Возможно создание более сложных сетей и расчет маршрутов. Например, вы можете задать направление дорог.
Может, ловкач справится с проблемой, но я не уверен.