Как я могу создать объект linnet, начиная с объекта sf с геометрическим столбцом LINESTRING?

В настоящий момент я работаю над проектом с событиями точечных паттернов в линейной сети (автомобильные аварии) и читаю главу 17 spatstat книги: «Пространственные точечные паттерны: методология и приложения с R».

Авторы книги объясняют, что они определили новый класс объектов под названием lpp для анализа точечных паттернов в линейной сети. Каркас каждого lpp объекта - это объект linnet, и есть несколько функций для создания объекта linnet. Для моего приложения соответствующими функциями являются linnet и as.linnet. Функция linnet создает линейный сетевой объект из пространственного положения каждой вершины и информации о том, какие вершины соединены ребром, в то время как функция as.linnet может применяться к psp объекту, который преобразуется в linnet объекты, определяя связность сети, используя заданный порог расстояния.

Причина, по которой я задаю этот вопрос, заключается в том, что я не знаю, как эффективно создать объект linnet, начиная с объекта sf с геометрией LINESTRING. Насколько мне известно, можно преобразовать объект sf в объект sp (т.е. объект SpatialLines), затем я могу преобразовать объект sp в объект psp (с помощью функции as.psp), а затем я могу преобразовать объект psp в объект объект linnet с помощью функции as.psp.linnet (которая определена в пакете maptools). Основная проблема с этим подходом (как заявили авторы пакета в своей книге) заключается в том, что предполагаемая сеть неверна каждый раз, когда в моих сетевых данных возникает эстакада или подземный переход, поскольку соответствующая сеть будет создавать искусственные пересечения в новой сети. Более того, как написали авторы в своей книге, код становится экспоненциально медленнее.

Следующий код представляет собой упрощенную версию того, что я делал до сих пор, но я думаю, что должен быть более простой и лучший способ создать объект linnet из объекта sf. Я бы использовал функцию linnet, но проблема в том, что я не знаю, как создать (разреженную) матрицу смежности для соответствующих вершин сети или матрицу связей между краями сети.

# packages
library(sf)
#> Linking to GEOS 3.6.1, GDAL 2.2.3, PROJ 4.9.3
library(spatstat)
#> Loading required package: spatstat.data
#> Loading required package: nlme
#> Loading required package: rpart
#> 
#> spatstat 1.61-0       (nickname: 'Puppy zoomies') 
#> For an introduction to spatstat, type 'beginner'
#> 
#> Note: spatstat version 1.61-0 is out of date by more than 11 weeks; a newer version should be available.
library(maptools)
#> Loading required package: sp
#> Checking rgeos availability: TRUE
library(osmdata)
#> Data (c) OpenStreetMap contributors, ODbL 1.0. http://www.openstreetmap.org/copyright


# download data
iow_polygon <- getbb("Isle of Wight, South East, England", format_out = "sf_polygon", featuretype = "state") %>% 
  st_transform(crs = 27700)
iow_highways <- st_read("https://download.geofabrik.de/europe/great-britain/england/isle-of-wight-latest.osm.pbf", layer = "lines") %>% 
  st_transform(crs = 27700)
#> Reading layer `lines' from data source `https://download.geofabrik.de/europe/great-britain/england/isle-of-wight-latest.osm.pbf' using driver `OSM'
#> Simple feature collection with 44800 features and 9 fields
#> geometry type:  LINESTRING
#> dimension:      XY
#> bbox:           xmin: -5.716262 ymin: 43.35489 xmax: 1.92832 ymax: 51.16517
#> epsg (SRID):    4326
#> proj4string:    +proj=longlat +datum=WGS84 +no_defs

# subset the data otherwise the code takes ages
iow_highways <- iow_highways[iow_polygon, ] %>% 
  subset(grepl(pattern = c("primary|secondary|tertiary"), x = highway))

# transform as sp 
iow_highways_sp <- as(iow_highways %>% st_geometry(), "Spatial")

# transform as psp
iow_highways_psp <- as.psp(iow_highways_sp)

# transform as linnet
iow_highways_linnet <- as.linnet.psp(iow_highways_psp, sparse = TRUE)

Я могу извлечь координаты каждой вершины сети

stplanr::line2points(iow_highways)
#> Simple feature collection with 2814 features and 1 field
#> geometry type:  POINT
#> dimension:      XY
#> bbox:           xmin: 430780.7 ymin: 75702.05 xmax: 464851.7 ymax: 96103.72
#> epsg (SRID):    27700
#> proj4string:    +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +towgs84=446.448,-125.157,542.06,0.15,0.247,0.842,-20.489 +units=m +no_defs
#> First 10 features:
#>    id                  geometry
#> 1   1 POINT (464851.7 87789.73)
#> 2   1 POINT (464435.4 88250.85)
#> 3   2 POINT (464390.9 87412.27)
#> 4   2 POINT (464851.7 87789.73)
#> 5   3 POINT (462574.6 88987.62)
#> 6   3 POINT (462334.6 88709.92)
#> 7   4 POINT (464066.9 87576.84)
#> 8   4 POINT (464390.9 87412.27)
#> 9   5   POINT (464420 88227.79)
#> 10  5 POINT (464398.7 88225.33)  

но тогда я не знаю, как построить матрицу смежности.

Создано 2 декабря 2019 г. пакетом REPEX (v0.3.0)


person agila    schedule 02.12.2019    source источник


Ответы (2)


Я не понимаю, почему вы используете формат psp на пути к linnet. Попробуйте заменить последние две строки вашего первого фрагмента кода на:

iow_highways_linnet <- as.linnet.SpatialLines(iow_highways_sp)

Это преобразует SpatialLines напрямую в linnet и объединяет линии, которые имеют общую вершину. Я не думаю, что подземный переход будет соединен с путепроводом, если обе линии не будут иметь вершину в точке пересечения. См. Пример ниже:

l1 <- sf::st_linestring(matrix(c(-1,1,-1,1,1,1,-1,-1), ncol = 2))
l2 <- sf::st_linestring(matrix(c(-1,-1,1,1,2,1,-1,-2), ncol = 2))
l_sf <- sf::st_sf(id = 1:2, geom = sf::st_sfc(l1,l2))
l_sp <- sf::as_Spatial(l_sf)
l <- maptools::as.linnet.SpatialLines(l_sp)
plot(l)

person Ege Rubak    schedule 04.12.2019
comment
Извините, я не видел вашего ответа. Я никогда не использовал функцию as.linnet.SpatialLines, потому что не знал, что она существует. В любом случае я только что протестировал ваш подход на примере IOW, и он работает без сбоев. Спасибо! - person agila; 06.12.2019

Просто подтверждаю, что пакет spatstat не предоставляет функций для обработки других форматов; мы ожидаем, что maptools или другие пакеты предоставят код преобразования формата; это еще не доступно для sf объектных форматов, предположительно потому, что sf относительно новый.

Ключевой вопрос заключается в том, содержит ли объект sf с геометрией LINESTRING достаточно информации для определения возможности подключения к сети. Если это так, то я предлагаю вам создать матрицу из двух столбцов, в которой перечислены все пары вершин, соединенных ребрами, и вызвать spatstat::linnet. Если нет, то имеющихся данных недостаточно ...

Наконец, обратите внимание, что текущая разрабатываемая версия spatstat (1.61-0.061), доступная на GitHub, намного быстрее чем текущая версия (1.61-0) для многих операций в линейных сетях. Скоро он будет опубликован.

person Adrian Baddeley    schedule 03.12.2019