3

可能重复：
加快 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) )

}
4

3 回答 3

12

使用过滤器，您甚至可以在没有任何循环的情况下进行计算（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你经常循环使用它。然后是过滤器。您的计算只不过是n在每个 time value 总结过去的值j。因此，我们有类似的东西

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但是，请注意 R 帮助中的以下注释filter：

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

于 2013-01-25T10:55:41.963 回答
3

矢量化计算并不总是意味着使用 *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名称用于两种不同的事物。

于 2013-01-25T10:48:18.077 回答
0

假设您的显式循环是正确的计算，请尝试以下操作：

 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 之间的第 iter-th 索引绝对差。然后将其添加到 signal[j]

于 2013-01-25T10:27:52.417 回答