在R中有效地计算三维数组的行和。

Consider the array a:

考虑数组:

> a <- array(c(1:9, 1:9), c(3,3,2))
> a
, , 1

     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6    9

, , 2

     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6    9

How do we efficiently compute the row sums of the matrices indexed by the third dimension, such that the result is:

如何有效地计算由第三维索引的矩阵的行和，使结果为:

     [,1] [,2]
[1,]   12   12
[2,]   15   15
[3,]   18   18

? ?

The column sums are easy via the 'dims' argument of colSums():

通过colsum()的“dims”论点，列和很容易:

> colSums(a, dims = 1)

but I cannot find a way to use rowSums() on the array to achieve the desired result, as it has a different interpretation of 'dims' to that of colSums().

但是我无法找到在数组上使用rowsum()来实现所需结果的方法，因为它对dims的解释与colsum()不同。

It is simple to compute the desired row sums using:

用以下方法计算所需的行和是很简单的:

> apply(a, 3, rowSums)
     [,1] [,2]
[1,]   12   12
[2,]   15   15
[3,]   18   18

but that is just hiding the loop. Are there other efficient, truly vectorised, ways of computing the required row sums?

但这只是隐藏的循环。有没有其他有效的、真正的矢量化的方法来计算所需的行和?

4 个解决方案

#1

@Fojtasek's answer mentioned splitting up the array reminded me of the aperm() function which allows one to permute the dimensions of an array. As colSums() works, we can swap the first two dimensions using aperm() and run colSums() on the output.

@Fojtasek的回答提到了分割数组，这让我想起了aperm()函数，它允许一个数组的维数。当colsum()工作时，我们可以使用aperm()交换前两个维度，并在输出上运行colsum()。

> colSums(aperm(a, c(2,1,3)))
     [,1] [,2]
[1,]   12   12
[2,]   15   15
[3,]   18   18

Some comparison timings of this and the other suggested R-based answers:

一些比较时间的这个和其他建议的基于r的答案:

> b <- array(c(1:250000, 1:250000),c(5000,5000,2))
> system.time(rs1 <- apply(b, 3, rowSums))
   user  system elapsed 
  1.831   0.394   2.232 
> system.time(rs2 <- rowSums3d(b))
   user  system elapsed 
  1.134   0.183   1.320 
> system.time(rs3 <- sapply(1:dim(b)[3], function(i) rowSums(b[,,i])))
   user  system elapsed 
  1.556   0.073   1.636
> system.time(rs4 <- colSums(aperm(b, c(2,1,3))))
   user  system elapsed 
  0.860   0.103   0.966

So on my system the aperm() solution appears marginally faster:

所以在我的系统中，aperm()解出现得稍微快一些:

> sessionInfo()
R version 2.12.1 Patched (2011-02-06 r54249)
Platform: x86_64-unknown-linux-gnu (64-bit)

However, rowSums3d() doesn't give the same answers as the other solutions:

然而，rowSums3d()并没有给出与其他解决方案相同的答案:

> all.equal(rs1, rs2)
[1] "Mean relative difference: 0.01999992"
> all.equal(rs1, rs3)
[1] TRUE
> all.equal(rs1, rs4)
[1] TRUE

#2

You could chop up the array into two dimensions, compute row sums on that, and then put the output back together the way you want it. Like so:

你可以把数组分割成二维，计算行和，然后按照你想要的方式把输出放回一起。像这样:

rowSums3d <- function(a){
    m <- matrix(a,ncol=ncol(a))
    rs <- rowSums(m)
    matrix(rs,ncol=2)
}

> a <- array(c(1:250000, 1:250000),c(5000,5000,2))
> system.time(rowSums3d(a))
   user  system elapsed 
   1.73    0.17    1.96 
> system.time(apply(a, 3, rowSums))
   user  system elapsed 
   3.09    0.46    3.74

#3

I don't know about the most efficient way of doing this, but sapply seems to do well

我不知道最有效的方法，但sapply似乎做得很好

a <- array(c(1:9, 1:9), c(3,3,2))
x1 <- sapply(1:dim(a)[3], function(i) rowSums(a[,,i]))
x1
     [,1] [,2]
[1,]   12   12
[2,]   15   15
[3,]   18   18

x2 <- apply(a, 3, rowSums)
all.equal(x1, x2)
[1] TRUE

Which gives a speed improvement as follows:

提高速度如下:

> a <- array(c(1:250000, 1:250000),c(5000,5000,2))

> summary(replicate(10, system.time(rowSums3d(a))[3]))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  2.784   2.799   2.810   2.814   2.821   2.862 

> summary(replicate(10, system.time(apply(a, 3, rowSums))[3]))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  2.730   2.755   2.766   2.776   2.788   2.839 

> summary(replicate(10, system.time( sapply(1:dim(a)[3], function(i) rowSums(a[,,i])) )[3]))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  1.840   1.852   1.867   1.872   1.893   1.914

Timings were done on:

时间是:

# Ubuntu 10.10
# Kernal Linux 2.6.35-27-generic
> sessionInfo()
R version 2.12.1 (2010-12-16)
Platform: x86_64-pc-linux-gnu (64-bit)

#4

If you have a multi-core system you could write a simple C function and make use of the Open MP parallel threading library. I've done something similar for a problem of mine and I get an 8 fold increase on an 8 core system. The code will still work on a single-processor system and even compile on a system without OpenMP, perhaps with a smattering of #ifdef _OPENMP here and there.

如果您有一个多核系统，您可以编写一个简单的C函数并使用Open MP并行线程库。我在我的问题上做过类似的事情，在8个核心系统上我得到了8倍的增长。该代码仍然可以在单处理器系统上运行，甚至可以在没有OpenMP的系统上编译，可能只需要少量的#ifdef _OPENMP就可以了。

Of course its only worth doing if you know that's what's taking most of the time. Do profile your code before optimising.

当然，如果你知道这是大部分时间所需要做的事情，那就值得去做。在优化之前先对代码进行概要分析。

#1

> colSums(aperm(a, c(2,1,3)))
     [,1] [,2]
[1,]   12   12
[2,]   15   15
[3,]   18   18