2

R で big.matrix オブジェクトの基本的な C++ コードを実装しようとしています。Rcpp パッケージを使用しています。ここでデモを読み、 rcpp-devel リストで見つけた別の単純な関数を適用しました。

#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
4

1 に答える 1

2

問題は転置アルゴリズムにあります。すべての行と列をループしているため、swaps が再スワップされています。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)
##
*/

の出力と の出力を比較すると、2 番目のバージョンが正しく機能することがわかります。ここで、BigTranspose2はと に相当します。t(M)Mmatrixb1b2

> 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

これは正方行列では正しく機能しますが、任意の次元の行列で機能するには、いくつかの変更を加える必要があることに注意してください。

于 2014-10-30T18:25:13.617 に答える