R中多个big.matrix对象的逐元素平均值

Question

R中多个big.matrix对象的逐元素平均值

3

我有17个文件支持的大矩阵对象（维度为10985 x 52598，每个大小为4.3GB），我想计算逐元素平均值。结果可以存储在另一个大矩阵（gcm.res.outputM）中。

biganalytics :: apply（）不起作用，因为MARGIN只能设置为1或2。我尝试使用两个循环，如下所示

gcm.res.outputM <- filebacked.big.matrix(10958, 52598, separated = FALSE, backingfile = "gcm.res.outputM.bin", backingpath = NULL, descriptorfile = "gcm.res.outputM.desc", binarydescriptor = FALSE)

for(i in 1:10958){
   for(j in 1:52598){
    t <- rbind(gcm.res.output1[i,j], gcm.res.output2[i,j],gcm.res.output3[i,j], gcm.res.output4[i,j],
           gcm.res.output5[i,j], gcm.res.output6[i,j],gcm.res.output7[i,j], gcm.res.output8[i,j],
           gcm.res.output9[i,j], gcm.res.output10[i,j],gcm.res.output11[i,j], gcm.res.output12[i,j],
           gcm.res.output13[i,j], gcm.res.output14[i,j],gcm.res.output15[i,j], gcm.res.output16[i,j],
           gcm.res.output17[i,j])
    tM <- apply(t, 2, mean, na.rm = TRUE)
    gcm.res.outputM[i,j] <- tM
    }
}

这个计算大约需要每行i1.5分钟，因此需要运行约11天。

有没有任何关于如何加快这个计算的想法？我正在使用一台64x Windows10机器，具有16GB的RAM。

谢谢！

- Jdh

看一下这个网址 http://winvector.github.io/Accumulation/。使用 data.table 包可能会有所帮助。 - Tung

有很多方法可以实现。您可以使用Rcpp实现，或者在所有矩阵的列块上使用R函数，或逐个添加矩阵。这些解决方案在实现的易用性、速度和内存使用量之间进行权衡。您有多少RAM？ - F. Privé

感谢您的评论。@F Privé，我有16GB的RAM。使用big.matrix对象的优点是我不必直接将整个矩阵加载到RAM中。我可以一次在整个值行上进行上述计算，因此摆脱第二个for循环。如果您有更多建议，我很乐意听取！ - Jdh

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- F. Privé · Accepted Answer

您可以使用以下Rcpp代码：

// [[Rcpp::depends(BH, bigmemory, RcppEigen)]]
#include <bigmemory/MatrixAccessor.hpp>
#include <RcppEigen.h>
using namespace Eigen;
using namespace Rcpp;

// [[Rcpp::export]]
void add_to(XPtr<BigMatrix> xptr_from, XPtr<BigMatrix> xptr_to) {

  Map<MatrixXd> bm_from((double *)xptr_from->matrix(),
                        xptr_from->nrow(), xptr_from->ncol());
  Map<MatrixXd> bm_to((double *)xptr_to->matrix(),
                      xptr_to->nrow(), xptr_to->ncol());

  bm_to += bm_from;
}

// [[Rcpp::export]]
void div_by(XPtr<BigMatrix> xptr, double val) {

  Map<MatrixXd> bm((double *)xptr->matrix(),
                   xptr->nrow(), xptr->ncol());

  bm /= val;
}

如果您有一个大小相同的big.matrix对象列表，则可以执行以下操作：

library(bigmemory)
bm_list <- lapply(1:5, function(i) big.matrix(1000, 500, init = i))
res <- deepcopy(bm_list[[1]])
lapply(bm_list[-1], function(bm) add_to(bm@address, res@address))
res[1:5, 1:5]  # verif
div_by(res@address, length(bm_list))
res[1:5, 1:5]  # verif