Расстояния точек между рядами с sf

У меня есть несколько траекторий, сохраненных в простой функции (sf) типа POINT, Я хотел бы рассчитать евклидово расстояние между последующими местоположениями (то есть рядами). До сих пор я "вручную" рассчитывал расстояния, используя формулу Пифагора для расчета евклидовых расстояний в двумерном пространстве. Мне было интересно, смогу ли я сделать то же самое, используя функцию sf::st_distance(), Вот быстрый пример:

library(sf)
library(dplyr)

set.seed(1)

df <- data.frame(
  gr = c(rep("a",5),rep("b",5)),
  x  = rnorm(10),
  y = rnorm(10)
 )

df <- st_as_sf(df,coords = c("x","y"),remove = F)


df %>%
  group_by(gr) %>%
  mutate(
    dist = sqrt((lead(x)-x)^2+(lead(y)-y)^2)
  )
#> Simple feature collection with 10 features and 4 fields
#> geometry type:  POINT
#> dimension:      XY
#> bbox:           xmin: -0.8356286 ymin: -2.2147 xmax: 1.595281 ymax: 1.511781
#> epsg (SRID):    NA
#> proj4string:    NA
#> # A tibble: 10 x 5
#> # Groups:   gr [2]
#>    gr         x       y   dist                 geometry
#>    <fct>  <dbl>   <dbl>  <dbl>                  <POINT>
#>  1 a     -0.626  1.51    1.38     (-0.6264538 1.511781)
#>  2 a      0.184  0.390   1.44     (0.1836433 0.3898432)
#>  3 a     -0.836 -0.621   2.91   (-0.8356286 -0.6212406)
#>  4 a      1.60  -2.21    3.57        (1.595281 -2.2147)
#>  5 a      0.330  1.12   NA         (0.3295078 1.124931)
#>  6 b     -0.820 -0.0449  1.31  (-0.8204684 -0.04493361)
#>  7 b      0.487 -0.0162  0.992  (0.4874291 -0.01619026)
#>  8 b      0.738  0.944   0.204    (0.7383247 0.9438362)
#>  9 b      0.576  0.821   0.910    (0.5757814 0.8212212)
#> 10 b     -0.305  0.594  NA       (-0.3053884 0.5939013)

Я хотел бы рассчитать dist с sf::st_distance(), Как бы я пошел по этому поводу?

Источник

2 ответа

Решение

Первое, что нужно знать о sf является то, что столбец геометрии (один из класса sfc) хранится в виде списка-столбца внутри фрейма данных. Ключ, который обычно делает что-либо со списком-столбцом, - это либо использовать purrr::map и друзья или использовать функцию, которая принимает списочные столбцы в качестве аргументов. В случае st_distance его аргументы могут быть объектом sf (датафрейм), sfc (столбец геометрии), или даже sfg (один ряд geom), поэтому нет необходимости map и друзья. Решение должно выглядеть примерно так:

df %>%
  group_by(gr) %>%
  mutate(
    dist = st_distance(geometry)
  )

Тем не менее, это не работает. После некоторых исследований мы обнаруживаем две проблемы. Первый, st_distance возвращает матрицу расстояний, а не одно значение. Чтобы решить эту проблему, мы используем by_element = T аргумент st_distance,

Далее мы не можем просто сделать dist = st_distance(geometry, lead(geometry), by_element = T) так как lead работает только для векторных столбцов, но не для столбцов списка.

Чтобы решить эту вторую проблему, мы сами создаем ведущую колонку, используя geometry[row_number() + 1],

Вот полное решение:

library(sf)
library(dplyr)

df %>%
  group_by(gr) %>%
  mutate(
    lead = geometry[row_number() + 1],
    dist = st_distance(geometry, lead, by_element = T),
  )
#> Simple feature collection with 10 features and 4 fields
#> Active geometry column: geometry
#> geometry type:  POINT
#> dimension:      XY
#> bbox:           xmin: -0.8356286 ymin: -2.2147 xmax: 1.595281 ymax: 1.511781
#> epsg (SRID):    4326
#> proj4string:    +proj=longlat +datum=WGS84 +no_defs
#> # A tibble: 10 x 6
#> # Groups:   gr [2]
#>    gr         x       y  dist                       geometry
#>    <fct>  <dbl>   <dbl> <dbl>         <sf_geometry [degree]>
#>  1 a     -0.626  1.51   1.38     POINT (-0.6264538 1.511781)
#>  2 a      0.184  0.390  1.44     POINT (0.1836433 0.3898432)
#>  3 a     -0.836 -0.621  2.91   POINT (-0.8356286 -0.6212406)
#>  4 a      1.60  -2.21   3.57        POINT (1.595281 -2.2147)
#>  5 a      0.330  1.12   0         POINT (0.3295078 1.124931)
#>  6 b     -0.820 -0.0449 1.31  POINT (-0.8204684 -0.04493361)
#>  7 b      0.487 -0.0162 0.992  POINT (0.4874291 -0.01619026)
#>  8 b      0.738  0.944  0.204    POINT (0.7383247 0.9438362)
#>  9 b      0.576  0.821  0.910    POINT (0.5757814 0.8212212)
#> 10 b     -0.305  0.594  0       POINT (-0.3053884 0.5939013)
#> # ... with 1 more variable: lead <sf_geometry [degree]>
Источник

Вот базаRрешение. Хитрость заключается в том, чтобы

  1. split()кадр данных по столбцу группы
  2. По каждой группе
    • рассчитать расстояние между последующими местоположениями, используяhead()иtail()
    • добавитьNAк результату (поскольку последняя позиция не имеет последующего местоположения)
  3. Свяжите получившийсяsfобъекты собрать вместе, используяrbindвdo.call().
      split(df, df$gr) |> 
  lapply(\(x){
    x$dist <- c(st_distance(head(x,-1), tail(x,-1),by_element = TRUE), NA)
    x
  }) |> 
  do.call(rbind, args = _)
Источник
Другие вопросы по тегам