Быстрое подмножество матрицы с RcppGSL

Question

После просмотра этого поста я попытался установить матрицу с Rcpp,

С RcppArmadillo:

// [[Rcpp::depends(RcppArmadillo)]]
#include "RcppArmadillo.h"
// [[Rcpp::export]]
arma::mat submatrix(const arma::mat& m1in, int fromin, int toin){
  arma::mat s1 = m1in.cols(fromin-1,toin-1);
  return(s1);
}

затем submatrix(M, 1, 900) немного быстрее чем M[,1:900],

С RcppGSL:

#include <RcppGSL.h>
#include <gsl/gsl_matrix.h>
// [[Rcpp::export]]
gsl_matrix_const_view submatrix(const RcppGSL::Matrix & X, int k1, int k2, int n1, int n2) {
  return gsl_matrix_const_submatrix(X, k1, k2, n1, n2);
}

Вот submatrix(M, 0, 0, 1000, 900) медленнее, чем M[,1:900]:

> microbenchmark(M[,1:900], submatrix(M, 0, 0, 1000, 900))
Unit: milliseconds
                          expr       min       lq     mean   median       uq      max neval
                    M[, 1:900]  8.035749 10.20265 13.25657 11.75554 14.27586 117.2533   100
 submatrix(M, 0, 0, 1000, 900) 16.597605 19.55858 23.04454 21.52959 23.98431 141.6158   100

Есть ли более быстрый способ подмножества матрицы с RcppGSL?

Answer 1

Я думаю, что причина в том, что ваша матрица не передана по ссылке (возможно, потому что матрицы R и матрицы GSL не совместимы).

Чтобы доказать мою точку зрения, проверьте это:

// [[Rcpp::depends(RcppGSL)]]
#include <RcppGSL.h>
#include <gsl/gsl_matrix.h>

// [[Rcpp::export]]
gsl_matrix_const_view submatrix(RcppGSL::Matrix & X, int k1, int k2, int n1, int n2) {
  X(0, 0) = 1;
  return gsl_matrix_const_submatrix(X, k1, k2, n1, n2);
}

/*** R
M <- matrix(0, 1000, 1000)
test <- submatrix(M, 0, 0, 1000, 900)
M[1, 1]
*/

Если я прав, у вас будет одна и та же проблема каждый раз, когда вы используете RcppGSL. Может быть, существует представление карты матрицы (как в Eigen), чтобы использовать вместо & (Я не знаю GSL).