R代码运行太慢，如何加速和重写这段代码答案

【问题标题】：R code runs too slow, how to speed and rewrite this codeR代码运行太慢，如何加速和重写这段代码
【发布时间】：2016-04-15 02:19:23
【问题描述】：

stu.csv 包含 850,000 行和 3 列。第 2 列是 ID 的经度，第 3 列是 ID 的纬度。 stu.csv 文件中的数据是这样的：

   ID    longitude    latitude  
  156   41.88367183 12.48777756
  187   41.92854333 12.46903667
  297   41.89106861 12.49270456
  89    41.79317669 12.43212196
  79    41.90027472 12.46274618
  ...       ...         ...

伪代码如下。它旨在计算地球表面两个ID之间的距离和经纬度，并输出任意两个ID的累积和：

  dlon = lon2 - lon1
  dlat = lat2 - lat1
  a = (sin(dlat/2))^2 + cos(lat1) * cos(lat2) * (sin(dlon/2))^2
  c = 2 * atan2( sqrt(a), sqrt(1-a) )
  distance = 6371000 * c (where 6371000 is the radius of the Earth)

这段代码如下，但是运行速度太慢了。如何加速和重写代码？谢谢。

    stu<-read.table("stu.csv",header=T,sep=",");

    ## stu has 850,000 rows and 3 columns.

    m<-nrow(stu);

    distance<-0;

    for (i in 1:(m-1))
    {
      for (j in (i+1))
      {     
        dlon = stu[j,2] - stu[i,2];
        dlat = stu[j,3] - stu[i,3];
        a = (sin(dlat/2))^2 + cos(stu[i,3]) * cos(stu[j,3]) * (sin(dlon/2))^2;
        c = 2 * atan2( sqrt(a), sqrt(1-a) );
        distance <-distance+6371000 * c;
       }
    }

    distance

【问题讨论】：

当 N=850000 时，O(N^2) 算法总是需要很长时间......
@chinsoon12 感谢您的回复。我描述了问题
除了用纯 C 编码并且仍然需要等待之外，您可能需要降低维数
只是为了澄清-您想要df中每个点与每个其他点之间的距离总和（这就是您的代码现在试图做的）？我在想你可能真的想要从一行到下一行的累积距离？
结果必须有多准确？如果所有点的经纬度相对（主观判断）彼此接近，您可能想提出一个新函数来近似距离计算，但取消昂贵的三角计算（即sin，cos）。用一个合理的常数或一些可以向量化的线性计算代替它们。

标签： r performance

【解决方案1】：

对于你的情况，如果是累积距离，我们可以向量化：

x <- read.table(text = "ID    longitude    latitude  
156   41.88367183 12.48777756
187   41.92854333 12.46903667
297   41.89106861 12.49270456
89    41.79317669 12.43212196
79    41.90027472 12.46274618", header= TRUE)


x$laglon <- dplyr::lead(x$longitude, 1)
x$laglat <- dplyr::lead(x$latitude, 1)


distfunc <- function(long, lat, newlong, newlat){
  dlon = newlong - long
  dlat = newlat - lat
  a = (sin(dlat/2))^2 + cos(lat) * cos(newlat) * (sin(dlon/2))^2
  c = 2 * atan2( sqrt(a), sqrt(1-a) )
  6371000 * c 
}

distfunc(x$longitude, x$latitude, x$laglon, x$laglat)
308784.6 281639.6 730496.0 705004.2       NA

取输出的cumsum得到总距离。

在一百万行上，我的系统大约需要 0.4 秒

【讨论】：

如何解决这个错误？ stackoverflow.com/questions/36744127/…

【解决方案2】：

您正在寻找的内容称为“矢量化循环”。见this related question。

基本思想是，在一个循环中，CPU 必须在处理第二个单元之前完成第一个单元的处理，除非它可以保证第一个单元的处理方式不会影响第二个单元的状态细胞。但如果它是向量计算，则存在这种保证，并且它可以同时处理尽可能多的元素，从而提高速度。（还有其他原因可以更好地工作，但这是基本动机。）

请参阅 R 中的 this introduction 到 apply，了解如何在没有循环的情况下重写代码。（大部分计算您应该能够保持原样。）

【讨论】：