Как преобразовать список пространственных точек 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.

Возможно создание более сложных сетей и расчет маршрутов. Например, вы можете задать направление дорог.

Может, ловкач справится с проблемой, но я не уверен.

Другие вопросы по тегам