在 R 中计算 3D 数组一维乘积的快速方法

Fast way to calculate product over one dimension of a 3D array in R

我有一个三维数组(例如维度 = 4000 x 4000 x 2)。现在我想计算第三个维度的乘积以获得一个二维数组(维度 = 4000 x 4000)作为结果。

我尝试在 apply() 函数中使用 prod() 计算乘积;然而,这非常耗时。因此,我想知道是否有更快更有效的方法来进行此类计算?

apply() 方法:

A <- array(runif(4000*4000*2),dim=c(4000,4000,2))
system.time(apply(A, c(1,2), prod))

这里是一个较小的例子,数组 B:

B <- array(c(1,2,1,2,3,4,3,4),dim=c(2,2,2))

结果B_res:

B_res <-  array(c(3,3,8,8),dim=c(2,2))

更新: 正如@42 所提到的,这可以通过元素明智(手动)乘法来完成,例如:B_res <- B[,,1]*B[,,2]。但是,第三维的大小可能在 2 到 x 之间。因此手动编码 B[,,1]*B[,,2]... *B[,,x] 可能不可行。这里循环计算可能是一种可能的解决方案:

array_prod <- function(C){
  C_res <- C[,,1]
  for(i in 2:dim(C)[3]){
    C_res <- C_res*C[,,i]
  }
  return(C_res)
}

这里是三种方法的比较(应用、手动逐元素和循环乘法):

A <- array(runif(400*400*10),dim=c(400,400,10))
system.time(apply(A, c(1,2), prod)); system.time(A[,,1]*A[,,2]*A[,,3]*A[,,4]*A[,,5]*A[,,6]*A[,,7]*A[,,8]*A[,,9]*A[,,10]); system.time(array_prod(A))
  user  system elapsed 
  0.492   0.021   0.512 
   user  system elapsed 
  0.031   0.000   0.032 
   user  system elapsed 
  0.032   0.001   0.032 

...这表明 apply 函数明显比其他两种方法慢得多,这两种方法基本上快得差不多。

这表明通过将前两个维度留空并使用 * 运算符,可以使用 R 中称为向量化的方法来实现逐元素数组乘法。也可以设置 TRUE 来表示特定维度的所有实例:

A <- array( 1:(4*4*2),dim=c(4,4,2))
apply(A, c(1,2), prod)
#============
     [,1] [,2] [,3] [,4]
[1,]   17  105  225  377
[2,]   36  132  260  420
[3,]   57  161  297  465
[4,]   80  192  336  512
#=============
A[ , , 1]*A[ , , 2]
     [,1] [,2] [,3] [,4]
[1,]   17  105  225  377
[2,]   36  132  260  420
[3,]   57  161  297  465
[4,]   80  192  336  512

这表明性能提高了 100 倍(虽然我厌倦了等待 apply to 运行 的 4000x4000 版本,所以我只显示该示例中矢量化方法的结果:)

> A <- array(runif(400*400*2),dim=c(400,400,2))
> system.time(apply(A, c(1,2), prod)); system.time(A[,,1]*A[,,2])
   user  system elapsed 
  0.448   0.018   0.452   # the apply timings
   user  system elapsed 
  0.005   0.000   0.004   # the vectorised operation

> A <- array(runif(4000*4000*2),dim=c(4000,4000,2))
>  system.time(A[,,1]*A[,,2])
   user  system elapsed 
  0.525   0.096   0.604