R Rcpp big.matrix加入

Question

我想在R中为big.matrix对象实现一些基本的C++代码。我正在使用Rcpp包，已经阅读演示here，甚至应用了另一个简单函数，我在rcpp-devel list ：

#include "bigmemory/BigMatrix.h" 
#include "bigmemory/MatrixAccessor.hpp" 
#include <Rcpp.h> 

using namespace Rcpp; 

// [[Rcpp::export]] 
void fun(SEXP A) { 
    Rcpp::XPtr<BigMatrix> bigMat(A); 
    MatrixAccessor<int> Am(*bigMat); 

    int nrows = bigMat->nrow(); 
    int ncolumns = bigMat->ncol(); 
    for (int j = 0; j < ncolumns; j++){ 
     for (int i = 1; i < nrows; i++){ 
       Am[j][i] = Am[j][i] + Am[j][i-1]; 
      } 
    } 
    return; 
} 

// [[Rcpp::export]] 
void BigTranspose(SEXP A) 
{ 
    Rcpp::XPtr<BigMatrix> pMat(A); 
    MatrixAccessor<int> mat(*pMat); 

    int r = pMat->nrow(); 
    int c = pMat->ncol(); 

    for(int i=0; i<r; ++i) 
     for(int j=0; j<c; ++j) 
     std::swap(mat[j][i], mat[i][j]); 

    return; 
} 

这fun功能工作完全正常，修改big.matrix对象。

a <- matrix(seq(25), 5,5) 
> a 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 6 11 16 21 
[2,] 2 7 12 17 22 
[3,] 3 8 13 18 23 
[4,] 4 9 14 19 24 
[5,] 5 10 15 20 25 
> fun(b@address) 
> head(b) 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 6 11 16 21 
[2,] 3 13 23 33 43 
[3,] 6 21 36 51 66 
[4,] 10 30 50 70 90 
[5,] 15 40 65 90 115 

但是，当我尝试一个简单的矩阵矩阵转置函数矩阵不被修改。为什么fun函数可以工作，但不是我的'BigTranspose`？

a <- matrix(seq(25), 5,5) 
> a 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 6 11 16 21 
[2,] 2 7 12 17 22 
[3,] 3 8 13 18 23 
[4,] 4 9 14 19 24 
[5,] 5 10 15 20 25 
b <- as.big.matrix(a) 
BigTranspose(b@address) 
> head(b) 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 6 11 16 21 
[2,] 2 7 12 17 22 
[3,] 3 8 13 18 23 
[4,] 4 9 14 19 24 
[5,] 5 10 15 20 25 

Answer 1

的问题是在你的转置算法;您正在循环所有行和列，因此您的swap正在重新交换。您可以通过j的低指数只是设置的j = i+1代替j = 0解决这个问题：

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

// [[Rcpp::export]] 
void BigTranspose(SEXP A) 
{ 
    Rcpp::XPtr<BigMatrix> pMat(A); 
    MatrixAccessor<int> mat(*pMat); 

    int r = pMat->nrow(); 
    int c = pMat->ncol(); 

    for(int i=0; i<r; i++){ 
     for(int j=0; j<c; j++){ 
     std::swap(mat[j][i], mat[i][j]); 
     } 
    } 
    return; 
} 

// [[Rcpp::export]] 
void BigTranspose2(SEXP A) 
{ 
    Rcpp::XPtr<BigMatrix> pMat(A); 
    MatrixAccessor<int> mat(*pMat); 

    int r = pMat->nrow(); 
    int c = pMat->ncol(); 

    for(int i=0; i<r; i++){ 
     for(int j=(i+1); j<c; j++){ 
     std::swap(mat[j][i], mat[i][j]); 
     } 
    } 
    return; 
} 

/*** R 
## 
b1 <- as.big.matrix(a) 
head(b1) 
## 
BigTranspose(b1@address) 
head(b1) 
## 
## 
b2 <- as.big.matrix(a) 
head(b2) 
## 
BigTranspose2(b2@address) 
head(b2) 
## 
## 
M <- matrix(1:25,ncol=5) 
t(M) 
## 
*/ 

你可以看到，第二个版本由BigTranspose2输出与的t(M)，其中M是matrix比较正常工作的b1和b2相当于：

> Rcpp::sourceCpp('bigMatTranspose.cpp') 

> b1 <- as.big.matrix(a) 

> head(b1) 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 6 11 16 21 
[2,] 2 7 12 17 22 
[3,] 3 8 13 18 23 
[4,] 4 9 14 19 24 
[5,] 5 10 15 20 25 

> ## 
> BigTranspose(b1@address) 

> head(b1) 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 6 11 16 21 
[2,] 2 7 12 17 22 
[3,] 3 8 13 18 23 
[4,] 4 9 14 19 24 
[5,] 5 10 15 20 25 

> ## 
> ## 
> b2 <- as.big.matrix(a) 

> head(b2) 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 6 11 16 21 
[2,] 2 7 12 17 22 
[3,] 3 8 13 18 23 
[4,] 4 9 14 19 24 
[5,] 5 10 15 20 25 

> ## 
> BigTranspose2(b2@address) 

> head(b2) 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 2 3 4 5 
[2,] 6 7 8 9 10 
[3,] 11 12 13 14 15 
[4,] 16 17 18 19 20 
[5,] 21 22 23 24 25 

> ## 
> ## 
> M <- matrix(1:25,ncol=5) 

> t(M) 
    [,1] [,2] [,3] [,4] [,5] 
[1,] 1 2 3 4 5 
[2,] 6 7 8 9 10 
[3,] 11 12 13 14 15 
[4,] 16 17 18 19 20 
[5,] 21 22 23 24 25 

只要注意这将正常工作的方阵，但你将不得不作出一些莫它可以与任意维度的矩阵一起工作。