对多个数组应用 R 函数，返回相同大小的数组答案

【问题标题】：Apply an R function over multiple arrays, returning an array of the same size对多个数组应用 R 函数，返回相同大小的数组
【发布时间】：2015-06-20 23:24:08
【问题描述】：

我有两个 2x2 矩阵数组，我想在每对 2x2 矩阵上应用一个函数。这是一个最小的示例，将 A 中的每个矩阵乘以 B 中的相应矩阵：

A <- array(1:20, c(5,2,2))
B <- array(1:20, c(5,2,2))
n <- nrow(A)
# Desired output: array with dimension 5x2x2 that contains 
# the product of each pair of 2x2 matrices in A and B.
C <- aperm(sapply(1:n, function(i) A[i,,]%*%B[i,,], simplify="array"), c(3,1,2))

这需要两个数组，每个数组有 5 个 2x2 矩阵，然后将每对 2x2 矩阵相乘，在 C 中得到所需的结果。

我当前的代码是这个丑陋的最后一行，使用 sapply 循环遍历第一个数组维度并从 A 和 B 中分别拉出每个 2x2 矩阵。然后我需要用 aperm() 置换数组维度以便拥有与原始数组相同的顺序（sapply(...,simplify="array") 使用第三维而不是第一个维度索引每个 2x2 矩阵）。

有没有更好的方法来做到这一点？我讨厌那里的那个丑陋的函数（i），这实际上只是一种伪造 for 循环的方法。而 aperm() 调用使这变得不那么可读。我现在的工作正常；我只是在寻找更像是惯用 R 的东西。

mapply() 将采用多个列表或向量，但它似乎不适用于数组。 plyr 的 aaply() 也很接近，但它不需要多个输入。我最接近的是使用 abind() 和 aaply() 将 A 和 B 打包到一个数组中一次使用 2 个矩阵，但这并不完全有效（它只获取前两个条目；我的索引已关闭）：

aaply(.data=abind(A,B,along=0), 1, function(ab) ab[1,,]%*%ab[2,,])

无论如何，这并不完全干净或清晰！

我试图将其作为一个最小的示例，但我的实际用例需要更复杂的矩阵对函数（而且我也希望将其扩展到两个以上的数组），所以我寻找可以概括和扩展的东西。

【问题讨论】：

标签： arrays r plyr

【解决方案1】：

D <- aaply(abind(A, B, along = 4), 1, function(x) x[,,1] %*% x[,,2])

这是使用abind 和aaply 的有效解决方案。

【讨论】：

【解决方案2】：

有时 for 循环是最容易遵循的。它还可以概括和扩展：

n <- nrow(A)
C <- A
for(i in 1:n) C[i,,] <- A[i,,] %*% B[i,,]

【讨论】：

是的，我应该提到这是另一个明显的解决方案。我希望有一种功能性的方法，但 R 似乎用数组做功能性事情的基础设施比用列表做的要少得多。我想避免循环以使这个比例更好地适应更复杂的 A 和 B 表达式，所以我没有到处都有所有索引。
在不损失可读性的情况下，您可以在 RcppArmadillo 中编写此循环，它与 R/matlab 的矩阵代数和索引非常相似。
循环应该和sapply一样扩展。我见过循环是最快解决方案的情况。在函数中隐藏索引很容易：combine <- function(A, B, FUN = '%*%') { ... } 如果您担心速度，那么将 5 维放在末尾而不是开头。此外，您可以使用矩阵列表来代替数组。

【解决方案3】：

R 的列表基础结构（看起来）比数组要好得多，所以我也可以通过将数组转换为矩阵列表来处理它，如下所示：

A <- alply(A, 1, function(a) matrix(a, ncol=2, nrow=2))
B <- alply(A, 1, function(a) matrix(a, ncol=2, nrow=2))
mapply(function(a,b) a%*%b, A, B, SIMPLIFY=FALSE)

我认为这比我上面所说的更简单，但我仍然希望听到更好的想法。

【讨论】：