1

我有一个带有两个邮政编码和相应的纬度和经度的大型数据集(2.6M 行),我正在尝试计算它们之间的距离。我主要使用该包geosphere来计算邮政编码之间的文森蒂椭圆体距离,但我的数据集需要大量时间。什么是实现这一点的快速方法?

我试过的

library(tidyverse)
library(geosphere)

zipdata <- select(fulldata,originlat,originlong,destlat,destlong)

## Very basic approach
for(i in seq_len(nrow(zipdata))){
  zipdata$dist1[i] <- distm(c(zipdata$originlat[i],zipdata$originlong[i]),
       c(zipdata$destlat[i],zipdata$destlong[i]),
       fun=distVincentyEllipsoid)
}

## Tidyverse approach 
zipdata <- zipdata%>%
 mutate(dist2 = distm(cbind(originlat,originlong), cbind(destlat,destlong), 
   fun = distHaversine))

这两种方法都非常缓慢。我知道 210 万行永远不会是一个“快速”的计算,但我认为它可以做得更快。我在较小的测试数据上尝试了以下方法,但没有任何运气,

library(doParallel)
cores <- 15
cl <- makeCluster(cores)
registerDoParallel(cl)

test <- select(head(fulldata,n=1000),originlat,originlong,destlat,destlong)

foreach(i = seq_len(nrow(test))) %dopar% {
  library(geosphere)
  zipdata$dist1[i] <- distm(c(zipdata$originlat[i],zipdata$originlong[i]),
       c(zipdata$destlat[i],zipdata$destlong[i]),
       fun=distVincentyEllipsoid) 
}
stopCluster(cl)

任何人都可以帮助我以正确的方式使用doParallelwithgeosphere或更好的方式来处理这个问题吗?

编辑:来自(一些)回复的基准

## benchmark
library(microbenchmark)
zipsamp <- sample_n(zip,size=1000000)
microbenchmark(
  dave = {
    # Dave2e
    zipsamp$dist1 <- distHaversine(cbind(zipsamp$patlong,zipsamp$patlat),
                                   cbind(zipsamp$faclong,zipsamp$faclat))
  },
  geohav = {
    zipsamp$dist2 <- geodist(cbind(long=zipsamp$patlong,lat=zipsamp$patlat),
                             cbind(long=zipsamp$faclong,lat=zipsamp$faclat),
                             paired = T,measure = "haversine")
  },
  geovin = {
    zipsamp$dist3 <- geodist(cbind(long=zipsamp$patlong,lat=zipsamp$patlat),
                             cbind(long=zipsamp$faclong,lat=zipsamp$faclat),
                             paired = T,measure = "vincenty")
  },
  geocheap = {
    zipsamp$dist4 <- geodist(cbind(long=zipsamp$patlong,lat=zipsamp$patlat),
                             cbind(long=zipsamp$faclong,lat=zipsamp$faclat),
                             paired = T,measure = "cheap")
  }
,unit = "s",times = 100)

# Unit: seconds
# expr        min         lq       mean     median         uq        max neval  cld
# dave 0.28289613 0.32010753 0.36724810 0.32407858 0.32991396 2.52930556   100    d
# geohav 0.15820531 0.17053853 0.18271300 0.17307864 0.17531687 1.14478521   100  b  
# geovin 0.23401878 0.24261274 0.26612401 0.24572869 0.24800670 1.26936889   100   c 
# geocheap 0.01910599 0.03094614 0.03142404 0.03126502 0.03203542 0.03607961   100 a  

一个简单的all.equal测试表明,对于我的数据集,haversine 方法等于 vincenty 方法,但与geodist包中的“便宜”方法具有“平均相对差异:0.01002573”。

4

3 回答 3

3

R 是一种向量化语言,因此该函数将对向量中的所有元素进行操作。由于您正在计算每一行的原始和目标之间的距离,因此循环是不必要的。矢量化方法大约是循环性能的 1000 倍。
同样distVincentyEllipsoid直接使用(或 distHaveersine 等)并绕过该distm功能也应该可以提高性能。

没有任何示例数据,此代码段未经测试。

library(geosphere)

zipdata <- select(fulldata,originlat,originlong,destlat,destlong)

## Very basic approach
zipdata$dist1 <- distVincentyEllipsoid(c(zipdata$originlong, zipdata$originlat), 
       c(zipdata$destlong, zipdata$destlat))

注意:要使大多数地圈功能正常工作,正确的顺序是:先经度,后纬度。

上面列出的 tidyverse 方法缓慢的原因是该distm函数正在计算每个起点和终点之间的距离,这将导致 200 万乘 200 万元素矩阵。

于 2019-08-21T02:56:53.493 回答
1

我使用@SymbolixAU 的建议来使用该geodist包对我的数据集执行 2.1M 距离计算。我发现它比每次测试的包都要快得多geosphere(我在我的主要问题中添加了其中一个)。measure=cheap选项中的选项使用geodist便宜的标尺方法,该方法在 100 公里以下的距离内具有低错误率。有关更多信息,请参阅地质学家插图。鉴于我的一些距离高于 100 公里,我决定使用文森蒂椭球测量。

于 2019-08-22T12:24:02.207 回答
1

如果您要使用 geosphere,我会使用像 distHaversine 这样的快速近似方法,或者仍然快速且非常精确的 distGeo 方法。(distVincenty* 这些主要是为了好奇而实现的)。

于 2019-08-25T20:01:03.197 回答