logo标签列表

在R中高效地复制矩阵

rmatrixreplication
16
我有一个矩阵并希望高效地将其复制n次（其中n是数据集中的观测数）。例如，如果我有一个矩阵 A
A <- matrix(1:15, nrow=3)
则我想要的输出形式为
rbind(A, A, A, ...) #n次
显然，有许多方法可以构造这样一个大矩阵，例如使用for循环或apply或类似的函数。但是，“矩阵复制函数”的调用发生在我的优化算法的核心位置，在运行程序期间会被调用数万次，因此循环、apply类型的函数以及任何类似的方法都不够高效。 （这样的解决方案基本上意味着对n进行循环执行数万次，这显然是低效的。）我已经尝试过使用普通的rep函数，但没有找到一种安排rep输出的方式使其成为所需格式的矩阵。
解决方案do.call("rbind", replicate(n, A, simplify=F))也太低效了，因为在这种情况下rbind使用得太频繁。（然后，我的程序总运行时间的约30%用于执行rbind。）
有人知道更好的解决方案吗？
rbind只在do.call的方式下使用了一次。可能是复制操作导致了性能问题。 - Matthew Plourde
我用 Rprof 进行了测试，结果显示 rbind 所需的时间是 replicate 的两倍左右。这个结果让我感到惊讶。 - Wolfgang Pößnecker
6个回答
26
两个更多的解决方案：
第一个是对问题中示例的修改。
do.call("rbind", rep(list(A), n))

第二种方法涉及到展开矩阵，复制它并重新组装它。
matrix(rep(t(A),n), ncol=ncol(A), byrow=TRUE)

由于要求效率，需要进行基准测试。

library("rbenchmark")
A <- matrix(1:15, nrow=3)
n <- 10

benchmark(rbind(A, A, A, A, A, A, A, A, A, A),
          do.call("rbind", replicate(n, A, simplify=FALSE)),
          do.call("rbind", rep(list(A), n)),
          apply(A, 2, rep, n),
          matrix(rep(t(A),n), ncol=ncol(A), byrow=TRUE),
          order="relative", replications=100000)

这将会给出：

                                                 test replications elapsed
1                 rbind(A, A, A, A, A, A, A, A, A, A)       100000    0.91
3                   do.call("rbind", rep(list(A), n))       100000    1.42
5  matrix(rep(t(A), n), ncol = ncol(A), byrow = TRUE)       100000    2.20
2 do.call("rbind", replicate(n, A, simplify = FALSE))       100000    3.03
4                                 apply(A, 2, rep, n)       100000    7.75
  relative user.self sys.self user.child sys.child
1    1.000      0.91        0         NA        NA
3    1.560      1.42        0         NA        NA
5    2.418      2.19        0         NA        NA
2    3.330      3.03        0         NA        NA
4    8.516      7.73        0         NA        NA

所以最快的方法是使用原始的rbind调用，但这假定n是固定的并且提前已知。如果n不是固定的，则最快的方法是do.call("rbind", rep(list(A), n)。这些是针对3x5矩阵和10次复制的情况。不同大小的矩阵可能会给出不同的排序。
编辑：
对于n = 600，结果的顺序不同（省略显式的rbind版本）：
A <- matrix(1:15, nrow=3)
n <- 600

benchmark(do.call("rbind", replicate(n, A, simplify=FALSE)),
          do.call("rbind", rep(list(A), n)),
          apply(A, 2, rep, n),
          matrix(rep(t(A),n), ncol=ncol(A), byrow=TRUE),
          order="relative", replications=10000)

提供
                                                 test replications elapsed
4  matrix(rep(t(A), n), ncol = ncol(A), byrow = TRUE)        10000    1.74
3                                 apply(A, 2, rep, n)        10000    2.57
2                   do.call("rbind", rep(list(A), n))        10000    2.79
1 do.call("rbind", replicate(n, A, simplify = FALSE))        10000    6.68
  relative user.self sys.self user.child sys.child
4    1.000      1.75        0         NA        NA
3    1.477      2.54        0         NA        NA
2    1.603      2.79        0         NA        NA
1    3.839      6.65        0         NA        NA

如果您使用明确的rbind版本，则比do.call("rbind", rep(list(A), n))版本稍快一些，但差别不大，并且比apply或matrix版本慢。因此，在这种情况下，对任意n的泛化不需要损失速度。

1
哇，非常感谢！然而，我的小基准测试表明，对于更大的n，比如n = 600，矩阵版本比“rbind”调用更有效率。在这种情况下，“matrix(rep(t(...”调用的版本最有效率。 - Wolfgang Pößnecker
10

可能这更有效率：

apply(A, 2, rep, n)
就像我之前说的那样，这种方法不能得出正确的结果。你可以自己试一下：A <- matrix(1:15, nrow=3); n <- 2; rbind(A,A); matrix(rep(A, n), ncol = ncol(A), byrow = TRUE)结果不同... 编辑：为什么我不能在评论中创建换行？ - Wolfgang Pößnecker
@WolfgangPößnecker 对不起，是我的错误。请查看我回答的更新。 - Sven Hohenstein
谢谢，这已经是相当大的改进了。不过我还在寻找更快的解决方案。;) - Wolfgang Pößnecker
很棒的答案，非常简单。 - David Veitch
3

还有另一种方法：

rep(1, n) %x% A
如果您想对复制进行“each”版本，其中 A 的每一行重复 n 次，则代码简单为：A %x% rep(1, n)。 - stofer
1
你可以使用索引。
A[rep(seq(nrow(A)), n), ]
1

我和原帖作者的目的相同，最终更新了@Brian Diggs的比较，包括所有其他发布的答案。希望我做得正确。

#install.packages("rbenchmark")
library("rbenchmark")
A <- matrix(1:15, nrow=3)
n <- 600

benchmark(do.call("rbind", replicate(n, A, simplify=FALSE)),
          do.call("rbind", rep(list(A), n)),
          apply(A, 2, rep, n),
          matrix(rep(t(A),n), ncol=ncol(A), byrow=TRUE),
          A[rep(seq(nrow(A)), n), ],
          rep(1, n) %x% A,
          apply(A, 2, rep, n),
          matrix(rep(as.integer(t(A)),n),nrow=nrow(A)*n,byrow=TRUE),
     order="relative", replications=10000)

#                                                                test replications elapsed relative user.self sys.self user.child sys.child
#5                                          A[rep(seq(nrow(A)), n), ]        10000    0.32    1.000      0.33     0.00         NA        NA
#8 matrix(rep(as.integer(t(A)), n), nrow = nrow(A) * n, byrow = TRUE)        10000    0.36    1.125      0.35     0.02         NA        NA
#4                 matrix(rep(t(A), n), ncol = ncol(A), byrow = TRUE)        10000    0.38    1.188      0.37     0.00         NA        NA
#3                                                apply(A, 2, rep, n)        10000    0.59    1.844      0.56     0.03         NA        NA
#7                                                apply(A, 2, rep, n)        10000    0.61    1.906      0.58     0.03         NA        NA
#6                                                    rep(1, n) %x% A        10000    1.44    4.500      1.42     0.02         NA        NA
#2                                  do.call("rbind", rep(list(A), n))        10000    1.67    5.219      1.67     0.00         NA        NA
#1                do.call("rbind", replicate(n, A, simplify = FALSE))        10000    5.03   15.719      5.02     0.01         NA        NA
0
将其转换为数组，复制内容并创建一个新矩阵，更新行数如何？
A <- matrix(...)
n = 2 # just a test

a = as.integer(A)
multi.a = rep(a,n)
multi.A = matrix(multi.a,nrow=nrow(A)*n,byrow=T)
1
我刚刚尝试了一下，它和上面的答案有相同的问题：不幸的是，你的建议并没有产生正确的结果。 - Wolfgang Pößnecker
1
使用 as.integer(t(A))。然后它就可以工作了：matrix(rep(as.integer(t(A)),n),nrow=nrow(A)*n,byrow=TRUE) - Mark Miller
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接