R中循环的替代方案[重复]答案

【问题标题】：Alternative to loops in R [duplicate]R中循环的替代方案[重复]
【发布时间】：2013-01-09 06:59:41
【问题描述】：

可能重复：
Speed up the loop operation in R

我有几个关于循环的问题。我知道 R 在矢量化计算中工作得更快，我想更改下面的代码以利用这一点。在论坛上查看其他一些答案，sapply 函数似乎能够替换内部 for 循环，但我正在生成一个零向量，因此出现错误。道仍然是 1000，我认为这是造成问题的原因。

我主要关心的是速度，因为我需要围绕整个算法创建一个循环，并绘制不同的 V 和 n 大小以进行进一步分析。

感谢您的帮助

替代循环

tao = 1000
L = (tao - 1)   
n = 10      
V = 5               
I = 10000                       
V_s = matrix(rnorm(I), I, 1)
V_b = matrix(rnorm(I), I, 1)

signal <- matrix(0, L, 1)  

for( j in (n:L)){

    sapply(((j-n+1):j),function (tao) signal[j] = signal[j] + abs(V_s[tao] - V_b[tao]))

    signal[j] = (signal[j] / (n * V) )

}

原始循环

tao = 1000
L = (tao - 1)   
n = 10      
V = 5               
I = 10000                       
V_s = matrix(rnorm(I), I, 1)
V_b = matrix(rnorm(I), I, 1)

signal <- matrix(0, L, 1)  

for( j in (n:L)){

    for( tao in ((j-n+1):j))    {

        signal[j] = (signal[j] + abs(V_s[tao] - V_b[tao]))

    }
        signal[j] = (signal[j] / (n * V) )

}

【问题讨论】：

您使用矩阵而不是向量是否有原因，而 V_s、V_b 和信号似乎只有一列？
如果您将browser() 放入for(tao) 循环中，您将能够检查函数的内部工作原理并了解发生了什么。
@juba ，我在后面的方法中使用了矩阵格式，但本质上是的，它是一个向量。
关于加速 R 代码的一般建议：stackoverflow.com/a/8474941/636656

标签： r for-loop

【解决方案1】：

使用过滤器，即使没有任何循环，您也可以进行计算（sapply 只不过是一个隐藏循环）。

absdif <- abs(V_s - V_b)
signal <- filter(absdif[1:L], rep(1/(n*V), n), sides=1)
signal[is.na(signal)] <- 0

但是，当您不习惯过滤器时，了解第二行中发生的事情并非易事。让我们仔细看看：

首先我们计算V_s 和V_b 的绝对差值，您循环使用它们。然后是过滤器。您的计算只不过是在每个时间值j 处总结n 过去的值。因此，我们有类似的东西

signal[j] <- sum(absdif[j-n+1:j])

这正是卷积过滤器所做的——将一些值相加——以一般形式通过乘以一些权重。对于所有值，我们选择1/(n*V) 作为权重，这对应于您在外循环中执行的规范化。最后一个参数sides=1 只是告诉过滤器只取过去的值（sides=2 表示sum(absdif[(j-n/2):(j+n/2)])）。

最后一行只是填充了开头的 NA 值（过滤器没有足够的数据来计算总和 - 这等于跳过第一个 n 值）。

最后，一些时机：

您的全循环解决方案：

   User      System       total 
  0.037       0.000       0.037

juba的解决方法：

   User      System       total 
  0.007       0.000       0.008

使用过滤器的解决方案：

   User      System       total 
  0.000       0.000       0.001

请注意，过滤器的概念已经得到很好的研究，并且可以非常快地完成。

编辑：如?filter 中所述，R 不使用标准filter 命令的快速傅立叶变换。通常，FFT 是实现卷积的最有效方式。但是，即使这样也可以通过用

替换 filter 命令来完成

signal <- convolve(absdif[1:L], rep(1/(n*V), n), type='filter')

请注意，现在第一个 n 条目已被删除，而不是设置为 NA。然而，结果是一样的。计时这个时间没有用 - 总时间低于system.time的三位数输出...但是，请注意filter的R帮助中的以下注释：

convolve(, type="filter") 使用 FFT 进行计算，因此对于单变量序列上的长过滤器可能更快，但它不返回时间序列（因此时间对齐不清楚），也不返回处理缺失值。例如，对于长度为 100 的一系列长度为 1000 的过滤器，过滤器更快

【讨论】：

不错的解决方案，不知道filter()的这种用法。
+1，虽然过滤器也可能是变相的循环...
@PaulHiemstra 也许，但它似乎非常有效:)
@Paul 好吧，矢量化始终只是一个隐藏循环，但其间有不同的优化步骤。过滤器应该与abs 运算符处于同一优化级别。但是，考虑到这一点，我们仍然可以做得更好。请稍后查看我的编辑。
感谢您的回答，我运行了您的建议并收到了与我预期不同的输出。我应该用我的两个 for 循环替换你的解决方案？

【解决方案2】：

向量化计算并不总是意味着使用 *apply 函数。

例如，您可以通过将第二个 for 循环替换为矢量索引来简化和加速操作：

for(j in (n:L)){
  sel <- (j-n+1):j
  signal[j] <- sum(abs(V_s[sel] - V_b[sel])) / (n*V)
}

对于这个解决方案，我系统上的执行时间是：

utilisateur     système      écoulé 
      0.008       0.004       0.009

而对于您的 for 循环它是：

utilisateur     système      écoulé 
       0.06        0.00        0.06

顺便说一句，您不应该将tao 名称用于两个不同的事物。

【讨论】：

没错。矢量化几乎总能提高速度。 *apply 函数很少这样做。
这大大缩短了处理时间，感谢您的帮助和意见。

【解决方案3】：

假设您的显式循环是正确的计算，试试这个：

 signal[j]<- signal[j] + 
              sapply((j-n+1):j, 
                   FUN = function(iter){ 
                           abs(V_s[iter] - V_b[iter])
                   }, V_s = V_s, V_b = V_b)

请注意，sapply 返回 V_s 和 V_b 之间的第 i 个索引绝对差。然后将其添加到 signal[j]

【讨论】：