计算数据帧R中多个距离的经度纬度

Calculate distance longitude latitude of multiple in dataframe R

我不确定这里是否适合提出我的问题(我是 R 和这个网站的新手)。 我的问题如下:如何计算经度和纬度点之间的距离?

我在此站点上搜索了我的问题的答案,但答案只考虑了 2 分(而我的数据集包含超过 207000 行)。

我有一个名为 'adsb_relevant_columns_correct_timedifference' 的数据框,其中包含以下列:呼号、高度、速度、方向、Date_Time、纬度、经度。

Callsign Altitude  Speed Direction   Date_Time             Latitude     Longitude
A118       18000    110  340         2017-11-06 22:28:09   70.6086      58.2959
A118       18500    120  339         2017-11-06 22:29:09   72.1508      58.7894
B222       18500    150  350         2017-11-08 07:28:09   71.1689      59.1234
D123       19000    150  110         2018-05-29 15:13:27   69.4523      68.1235

我想计算每次测量之间的距离(以米为单位)(每行是一个新测量值)并将其添加为名为 'Distance' 的新列。第一个距离计算应该出现在第二行,因为为了以后的目的。因此,第 'Distance' 列的第一行可以是零或 NA,这无关紧要。

所以,我想知道第一次测量(Lat 70.6086 和 Long 58.2959)和第二次测量(Lat 72.1508 和 58.7894)之间的距离。然后是第二次和第三次测量之间的距离。然后在第三次和第四次之间,以此类推超过207000次测量。

预期的输出应该是这样的:

Callsign Altitude  Speed Direction   Date_Time             Latitude     Longitude  Distance 
A118       18000    110  340         2017-11-06 22:28:09   70.6086      58.2959    NA  
A118       18500    120  339         2017-11-06 22:29:09   72.1508      58.7894    172000
B222       18500    150  350         2017-11-08 07:28:09   71.1689      59.1234    110000
D123       19000    150  110         2018-05-29 15:13:27   69.4523      68.1235    387000

我在 R 中找到了 distm 函数,为此我可以手动完成两次测量而不是完整的数据集。

distm(p1, p2, fun = distHaversine)

我试过下面的代码

adsb_relevant_columns_correct_timedifference <- mutate(adsb_relevant_columns_correct_timedifference, Distance =
distm(c(adsb_relevant_columns_correct_timedifference$Longitude, adsb_relevant_columns_correct_timedifference$Latitude),
      c(lag(adsb_relevant_columns_correct_timedifference$Longitude, adsb_relevant_columns_correct_timedifference$Latitude)), fun = distCosine))

但是,我得到了以下错误

Error in mutate_impl(.data, dots) : Evaluation error: Wrong length for a vector, should be 2.

很抱歉我的解释很长,但我希望我的问题很清楚。有人可以告诉我如何计算几个测量值之间的距离并将其作为新列添加到我的数据框中吗?

使用distm也是一种选择。它产生一个距离矩阵,但是:

library(geosphere)
p <- cbind(df$Longitude, df$Latitude)
distm(head(p, -1), tail(p, -1), fun = distHaversine)
#          [,1]     [,2]     [,3]
# [1,] 172569.2  69279.8 394651.3
# [2,]      0.0 109928.5 454096.2
# [3,] 109928.5      0.0 387356.2

那你可以

diag(distm(head(p, -1), tail(p, -1), fun = distHaversine))
# [1] 172569.2 109928.5 387356.2

您也可以直接使用 distHaversine 函数来避免计算整个矩阵:

distHaversine(head(p, -1), tail(p, -1))
# [1] 172569.2 109928.5 387356.2

所以

df$Distance <- c(NA, distHaversine(head(p, -1), tail(p, -1)))

首先我们加载数据和库

library(geosphere)
df <- data.frame(read.table(text ="
Callsign Altitude  Speed Direction   Date_Time             Latitude     Longitude
A118       18000    110  340         2017-11-06T22:28:09   70.6086      58.2959
A118       18500    120  339         2017-11-06T22:29:09   72.1508      58.7894
B222       18500    150  350         2017-11-08T07:28:09   71.1689      59.1234
D123       19000    150  110         2018-05-29T15:13:27   69.4523      68.1235"
, header=TRUE))

接下来我们创建一个包含两点之间距离的新列,从第二行开始。

df$distance[2:nrow(df)] <- sapply(2:nrow(df), 
    function(x) distm(df[x-1,c('Longitude', 'Latitude')], df[x,c('Longitude', 'Latitude')], fun = distHaversine))

这样做会导致

df
#  Callsign Altitude Speed Direction           Date_Time Latitude Longitude distance
#1     A118    18000   110       340 2017-11-06T22:28:09  70.6086   58.2959       NA
#2     A118    18500   120       339 2017-11-06T22:29:09  72.1508   58.7894 172569.2
#3     B222    18500   150       350 2017-11-08T07:28:09  71.1689   59.1234 109928.5
#4     D123    19000   150       110 2018-05-29T15:13:27  69.4523   68.1235 387356.2

您可以使用 sf 包来计算与 st_distance() 的距离。为了让一行到下一行之间的距离,你可以做类似的事情。 (我在这里假设 data 是你的 data.frame。

library('sf')    

# some points
poi <- data.frame(id = 1:4, 
                  long = c(103.864325, 103.874916, 103.989693, 103.789615), 
                  lat = c(1.281949, 1.305004, 1.359878, 1.404454),
                  name = c('Gardens by the Bay', 'National Stadium', 'Airport', 'Zoo'))

# coerce to sf object
poi <- st_as_sf(x = poi, coords = c('long', 'lat'))

# duplicate object with shifted rows
poi2 <- poi[c(1, 1:(nrow(poi) - 1)), ]

# compute pairwise distance
distance <- st_distance(x = poi, y = poi2)

# extract only the diagonal corresponding to the distance between row r and row r+1
distance <- diag(distance)

# add result to data
poi$distance <- distance

# first distance to NA
poi$distance[1] <- NA

您可以使用 distHaversine 函数代替 distm。在您的 mutate 调用中,您不应重复数据框并使用 $ 运算符,mutate 已经知道在哪里查找列。发生错误是因为您需要使用 cbind 而不是 c,因为 c 创建一个长向量,只是将列堆叠在一起,而 cbind 创建一个包含两列的数据框(在这种情况下你想要什么)。

library(geosphere)
library(dplyr)

mutate(mydata, 
       Distance = distHaversine(cbind(Longitude, Latitude),
                                cbind(lag(Longitude), lag(Latitude))))

#   Callsign Altitude Speed Direction           Date_Time Latitude Longitude Distance
# 1     A118    18000   110       340 2017-11-06T22:28:09  70.6086   58.2959       NA
# 2     A118    18500   120       339 2017-11-06T22:29:09  72.1508   58.7894 172569.2
# 3     B222    18500   150       350 2017-11-08T07:28:09  71.1689   59.1234 109928.5
# 4     D123    19000   150       110 2018-05-29T15:13:27  69.4523   68.1235 387356.2

使用 distCosine 有点棘手,因为如果缺少输入的纬度或经度之一,它不会 return NA。因此我稍微修改了函数并解决了问题:

modified_distCosine <- function(Longitude1, Latitude1, Longitude2, Latitude2) {
  if (any(is.na(c(Longitude1, Latitude1, Longitude2, Latitude2)))) {
    NA
  } else {
    distCosine(c(Longitude1, Latitude1), c(Longitude2, Latitude2))
  }
}

mutate(mydata, 
       Distance = mapply(modified_distCosine, 
                         Longitude, Latitude, lag(Longitude), lag(Latitude)))

#   Callsign Altitude Speed Direction           Date_Time Latitude Longitude Distance
# 1     A118    18000   110       340 2017-11-06T22:28:09  70.6086   58.2959       NA
# 2     A118    18500   120       339 2017-11-06T22:29:09  72.1508   58.7894 172569.2
# 3     B222    18500   150       350 2017-11-08T07:28:09  71.1689   59.1234 109928.5
# 4     D123    19000   150       110 2018-05-29T15:13:27  69.4523   68.1235 387356.2

这里我使用 mapply 来应用参数为 Longitude, Latitude, lag(Longitude), lag(Latitude) 的修改后的函数。
我很确定必须有一种更优雅的方式,但至少这是有效的。

数据

mydata <- structure(list(Callsign = c("A118", "A118", "B222", "D123"), 
                         Altitude = c(18000L, 18500L, 18500L, 19000L), 
                         Speed = c(110L, 120L, 150L, 150L), 
                         Direction = c(340L, 339L, 350L, 110L), 
                         Date_Time = c("2017-11-06T22:28:09", "2017-11-06T22:29:09", "2017-11-08T07:28:09", "2018-05-29T15:13:27"), 
                         Latitude = c(70.6086, 72.1508, 71.1689, 69.4523), 
                         Longitude = c(58.2959, 58.7894, 59.1234, 68.1235)), 
                    .Names = c("Callsign", "Altitude", "Speed", "Direction", "Date_Time", "Latitude", "Longitude"), 
                    class = "data.frame", row.names = c(NA, -4L))