【问题标题】:Calculating distances between all possible combinations of lats and longs计算所有可能的经纬度组合之间的距离
【发布时间】:2017-06-14 10:09:43
【问题描述】:

我无法让下面的代码工作。我正在尝试计算我的数据集中所有可能的经纬度组合之间的距离。

我将使用的示例输入数据:

p <- data.frame(lat=runif(6,-90,90), lon=runif(6,-180,180) );

我无法让下面的代码工作。距离函数不起作用,所以我尝试了distm,但这也给了我一条错误消息。错误信息列在代码下方。

d <- setNames(do.call(rbind.data.frame, 
                      combn(1:nrow(p), 2, simplify = FALSE)), 
              c('p1','p2'));
d$dist <- sapply(1:nrow(d), function(r){
    distance(p$lat[d$p1[r]], p$lat[d$p2[r]], p$lon[d$p1[r]], p$lon[d$p2[r]])
})

d$dist <- sapply(1:nrow(d), function(r){
    distm(p$lat[d$p1[r]], p$lat[d$p2[r]], p$lon[d$p1[r]], p$lon[d$p2[r]])
})

#> Error in distm(p$lat[d$p1[r]], p$lat[d$p2[r]], p$lon[d$p1[r]], p$lon[d$p2[r]]) : 
#>  unused argument (p$lon[d$p2[r]])

【问题讨论】:

  • 也许是sapply(combn(nrow(p), 2, simplify = FALSE), function(i){geosphere::distHaversine(p[i[1], 2:1], p[i[2], 2:1])})
  • 这似乎奏效了,但价值确实很大。单位是米吗? @alistaire
  • distHaversine 返回的值与 r(地球半径)参数的单位相同,默认为米。它们应该是很大的数字,因为您的样本数据存在很大差异。

标签: r latitude-longitude


【解决方案1】:

geosphere::distHaversine(以及大多数其他距离函数)是矢量化的,因此您可以一次在所有对上调用它。将所有内容放入一个不错的 data.frame 中,

p <- data.frame(lat = runif(6, -90, 90), 
                lon = runif(6, -180, 180))

# get row indices of pairs
row_pairs <- combn(nrow(p), 2)

# make data.frame of pairs
df_dist <- cbind(x = p[row_pairs[1,],], 
                 y = p[row_pairs[2,],])
# add distance column by calling distHaversine (vectorized) on each pair
df_dist$dist <- geosphere::distHaversine(df_dist[2:1], df_dist[4:3])

df_dist
#>          x.lat      x.lon      y.lat      y.lon     dist
#> 1   -10.281070 -156.30519  -7.027720 -104.76897  5677699
#> 1.1 -10.281070 -156.30519 -51.142344 -100.99517  6750255
#> 1.2 -10.281070 -156.30519  -3.979805 -141.43436  1785251
#> 1.3 -10.281070 -156.30519 -21.239130  -65.97719  9639637
#> 1.4 -10.281070 -156.30519  66.292704 -154.52851  8525401
#> 2    -7.027720 -104.76897 -51.142344 -100.99517  4923176
#> 2.1  -7.027720 -104.76897  -3.979805 -141.43436  4075742
#> 2.2  -7.027720 -104.76897 -21.239130  -65.97719  4459657
#> 2.3  -7.027720 -104.76897  66.292704 -154.52851  9085777
#> 3   -51.142344 -100.99517  -3.979805 -141.43436  6452943
#> 3.1 -51.142344 -100.99517 -21.239130  -65.97719  4502520
#> 3.2 -51.142344 -100.99517  66.292704 -154.52851 13833468
#> 4    -3.979805 -141.43436 -21.239130  -65.97719  8350236
#> 4.1  -3.979805 -141.43436  66.292704 -154.52851  7893225
#> 5   -21.239130  -65.97719  66.292704 -154.52851 12111227

或者,您可以使用geosphere::distm,它为您提供一个距离矩阵,其中包含不同格式的相同数据:

geosphere::distm(p[, 2:1])
#>         [,1]    [,2]     [,3]    [,4]     [,5]     [,6]
#> [1,]       0 5677699  6750255 1785251  9639637  8525401
#> [2,] 5677699       0  4923176 4075742  4459657  9085777
#> [3,] 6750255 4923176        0 6452943  4502520 13833468
#> [4,] 1785251 4075742  6452943       0  8350236  7893225
#> [5,] 9639637 4459657  4502520 8350236        0 12111227
#> [6,] 8525401 9085777 13833468 7893225 12111227        0

?distHaversine 所述,距离以米为单位。随心所欲地转换。另请注意,geosphere 的函数采用 lon/lat,而不是 lat/lon,因此需要反转列才能工作。

【讨论】:

  • 效果很好。我从 csv 中的实际数据框开始,因此我不得不将“p”重命名为该数据框。下一个问题:我只想关注
  • 例如,我有 x.well.label 和 y.well.label。但我只想保持彼此之间
  • 所以对于每个 x.well.label 我只想将 y.well.label 保持在 1320' x.Well.Label x.lat x.lon y.Well.Label y.lat y 。巴特利特 V H 26 32.52477 -92.75865 14037.173260
  • 只是事后的子集,例如p2 &lt;- p[p$distance &lt;= 1320, ]
  • 这似乎不起作用。它给了我一个 0 行的数据表...不知道我做错了什么@alistaire
猜你喜欢
  • 2012-10-13
  • 2012-01-26
  • 1970-01-01
  • 2021-09-05
  • 2018-06-11
  • 1970-01-01
  • 2011-04-11
  • 1970-01-01
  • 2016-06-27
相关资源
最近更新 更多