如何加速这种双 for 循环？

Posted 2023-03-12

【中文标题】如何加速这种双 for 循环？

【英文标题】：How to speed this kind of double for-loop?

【发布时间】：2011-05-14 00:33:57

【问题描述】：

我正在用 R 编写一个expectation-maximization algorithm。为了加快计算速度，我想对这个瓶颈进行矢量化处理。我知道 N 大约是 k 的一百倍。

MyLoglik = 0
for (i in c(1:N))

 for (j in c(1:k))
 
  MyLoglik = MyLoglik + MyTau[i,j]*log(MyP[j]*MyF(MyD[i,], MyMu[j,], MyS[[j]]))
 

还有这个矩阵列表：

MyDf.list <- vector("list", k)
for(i in 1:k)

 MyDf.list[[i]] <- matrix(0,d,d)
 for (j in c(1:N))
 
  MyDf.list[[i]] = MyDf.list[[i]] + MyTau[j,i]*as.numeric((MyD[j,]-MyMu[i,])) %*% t(as.numeric(MyD[j,]-MyMu[i,]))  
 
 MyDf.list[[i]] = MyDf.list[[i]] / MyM[i]

---

我使用以下方法加快了速度：

MyLoglik = 0
for (j in c(1:k))

 MyR= apply(MyD, 1, function(x) log(MyP[j]*MyF(x, MyMu[j,], MyS[[j]])))
 MyLoglik = MyLoglik + sum(MyTau[,j]*MyR)

和：

d = dim(MyD)[2]
MyDf.list <- vector("list", k)
for(i in 1:k)

 MyDf.list[[i]] <- matrix(0,d,d)
 MyR= apply(MyD, 1, function(x) as.numeric((x-MyMu[i,])) %*% t(as.numeric(x-MyMu[i,])))
 MyDf.list[[i]] = matrix(rowSums(t(MyTau[,i]*t(MyR))) / MyM[i],d,d)

【参考方案1】：

对于第一个，我假设 MyF 是您创建的函数？如果您可以确保它将您的矩阵和列表作为输入并输出一个矩阵，您可以执行以下操作：

MyLoglik = sum(MyTau%*%log(MyP)) + sum(MyTau*log(MyF(MyD, MyMu, MyS)))

对于第二个，我认为由于您将其作为列表进行，因此矢量化将更加困难。也许你可以有一个 3 维数组而不是矩阵列表？这样 MyDf.array[i,j,k] 的维度为 N, d, d（或 d, d, N）。

【讨论】：

第一个建议不错！关于第二个，我认为 list 是像在 Matlab 中一样获取数组的唯一方法。

检查 ?array - 它可以处理多个维度。

【参考方案2】：

我什至不愿意过早地提出这个建议，但这是在 R 中构建 C 扩展可能有意义的事情。对于具有定义（已知）大小的矩阵（您在此处拥有！），C 扩展并不那么难构建，我保证！这里最讨厌的一点可能是传入'myF'

我的 R 知识已经过时了，但是 for 循环（尤其是像这个！）曾经是残酷的。

也许计时并找出慢的部分会有所帮助吗？是我的F吗？如果你把它改成一个身份呢？

【讨论】：

感谢您的建议。 myF 是对 Multivariate normal distribution 的密度函数（***上的 pdf 缩写）的调用。这是一个快速的单线。这实际上是 N 上的循环，它占用了大部分时间。 N 为 500，而 k 为 4。

我还建议粘贴整个示例并提交到 R 列表。我不再做很多 R，但 E-M 是一个很好理解的领域！ （尤其是像 mvnorm 这样的简单 e-m 函数！）我猜这是 Solved (tm)。

【参考方案3】：

如果事情是对称的，您可以减少在内部循环中完成的工作：A[i,j] = A[j,i]

【讨论】：

感谢您的建议。不幸的是，这里没有对称性。