2

我有一个大矩阵,大约有 6000 万行和 150 个列(总共大约 90 亿个元素)。我已将此数据存储在一个big.matrix对象中(来自 package bigmemory)。现在,我希望计算每一行的总和,这是一个问题,因为它big.matrix是面向列的,所以据我所知,所有汇总函数都是面向列的(例如colsumcolmax等)并且没有可用的函数默认用于计算行总和。我当然可以apply(x, 1, sum),但这需要很长时间。我还可以一一遍历列并使用矢量化添加来添加它们:

mysum <- rep(0, nrow(x))
for (i in seq(ncol(x))) 
  mysum <- mysum + x[,i]

但这仍然需要 20 多分钟,而且显然不是最理想的,因为它每次通过循环都会创建一个新的 6000 万元素向量。似乎必须有一些更快的方法来做到这一点。

编辑

我通过一次处理大约一百万行的块,并在这些块上调用 rowSums,然后将结果连接起来,将这个时间缩短到了 10 分钟。不过,我仍然很想知道是否有优化的方法来做到这一点。

4

1 回答 1

2

我已经编写了一些 C++ 代码来执行此操作,改编自bigmemory Rcpp 库

rowSums.cpp

// [[Rcpp::depends(BH)]]
#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::depends(BH, bigmemory)]]
#include <bigmemory/MatrixAccessor.hpp>

#include <numeric>

// Logic for BigRowSums.
template <typename T>
NumericVector BigRowSums(XPtr<BigMatrix> pMat, MatrixAccessor<T> mat) {
    NumericVector rowSums(pMat->nrow(), 0.0);
    NumericVector value(1);
    for (int jj = 0; jj < pMat->ncol(); jj++) {
      for (int ii = 0; ii < pMat->nrow(); ii++) {
        value = mat[jj][ii];
        if (all(!is_na(value))) {
          rowSums[ii] += value[0];
        }   
      }   
    }   
    return rowSums;
}

// Dispatch function for BigRowSums
//
// [[Rcpp::export]]
NumericVector BigRowSums(SEXP pBigMat) {
    XPtr<BigMatrix> xpMat(pBigMat);

    switch(xpMat->matrix_type()) {
      case 1:
        return BigRowSums(xpMat, MatrixAccessor<char>(*xpMat));
      case 2:
        return BigRowSums(xpMat, MatrixAccessor<short>(*xpMat));
      case 4:
        return BigRowSums(xpMat, MatrixAccessor<int>(*xpMat));
      case 6:
        return BigRowSums(xpMat, MatrixAccessor<float>(*xpMat));
      case 8:
        return BigRowSums(xpMat, MatrixAccessor<double>(*xpMat));
      default:
        throw Rcpp::exception("unknown type detected for big.matrix object!");
    }   
}

在 R 中:

library(bigmemory)
library(Rcpp)
sourceCpp("rowSums.cpp")

m <- as.big.matrix(matrix(1:9, 3))
BigRowSums(m@address)
[1] 12 15 18
于 2014-07-11T04:19:13.747 回答