加快所有行对的点积并添加截距项

Speed up taking dot product of all pairs of rows and adding intercept term

我正在对 83 个数据点应用 100 个不同的模型模拟,并检查每个数据点的每个模型模拟的估计范围。

每个计算本身都是 342 个变量乘以 342 个系数的乘积,再加上一个截距。我在下面编写的代码可以正确计算这些值,但是速度慢得令人发指。有什么办法可以提高处理速度吗?

spec = read.csv(spectra)
coef = read.csv(coefficents)

shell     = matrix(data=NA,ncol=101,nrow=nrow(spec))
shell     = as.data.frame(shell)
heading = paste("Sim_",seq(1,100,1),sep="")
names(shell)[1] = "Filename"
names(shell)[2:101] = heading

shell[1] = spec[1]

for (i in 1:nrow(spec))
{
  for (j in 1:100)
  {
    shell[i,j+1] = sum(spec[i,2:341]*coef[j,3:342]) + coef[j,2]
  }
}

实际上,您正在执行 spec 的矩阵乘法与 coef 的转置,然后向每一列添加一个常数。您应该通过使用内置矩阵乘法函数 %*% 和向量化操作来缩放列来获得加速:

out <- cbind(spec[,1], t(t(spec[,2:341] %*% t(coef[1:100,3:342])) + coef[1:100,2]))

此 1-liner 的输出与原始 post 中代码的输出相同(稍作修改以获取和输出矩阵而不设置任何名称):

OP <- function(spec, coef) {
  shell = matrix(data=NA,ncol=101,nrow=nrow(spec))
  shell[,1] <- spec[,1]
  for (i in 1:nrow(spec)) {
    for (j in 1:100) {
      shell[i,j+1] = sum(spec[i,2:341]*coef[j,3:342]) + coef[j,2]
    }
  }
  shell
}
all.equal(out, OP(spec, coef))
# [1] TRUE

就运行时间而言,矢量化操作产生了显着的好处 (38x),即使对于这个相对较小的示例(spec 中有 1000 行)也是如此:

system.time(cbind(spec[,1], t(t(spec[,2:341] %*% t(coef[1:100,3:342])) + coef[1:100,2])))
#    user  system elapsed 
#   0.028   0.001   0.030 
system.time(OP(spec, coef))
#    user  system elapsed 
#   0.927   0.224   1.161 

数据:

set.seed(144)
spec <- matrix(rnorm(1000*341), nrow=1000)
coef <- matrix(rnorm(100*342), nrow=100)