1

NARcppEigen で実装されている関数内のベクトルから効率的に値を削除する必要があります。もちろんforループでもできますが、もっと効率の良い方法はないかと。

次に例を示します。

library(RcppEigen)
library(inline)

incl <- '
using  Eigen::Map;
using  Eigen::VectorXd;
typedef  Map<VectorXd>  MapVecd;
'

body <- '
const MapVecd         x(as<MapVecd>(xx)), y(as<MapVecd>(yy));
VectorXd              x1(x), y1(y);
int                   k(0);
for (int i = 0; i < x.rows(); ++i) {
 if (x.coeff(i)==x.coeff(i) && y.coeff(i)==y.coeff(i)) {
  x1(k) = x.coeff(i);
  y1(k) = y.coeff(i);
  k++;
 };
};
x1.conservativeResize(k);
y1.conservativeResize(k);
return Rcpp::List::create(Rcpp::Named("x") = x1,
                          Rcpp::Named("y") = y1);
'

na.omit.cpp <- cxxfunction(signature(xx = "Vector", yy= "Vector"), 
                   body, "RcppEigen", incl)

na.omit.cpp(c(1.5, NaN, 7, NA), c(7.0, 1, NA, 3))
#$x
#[1] 1.5
#
#$y
#[1] 7

私の使用例では、これをループ (Rcpp 関数内) で約 100 万回行う必要があり、ベクトルは非常に長くなる可能性があります (1000 要素と仮定しましょう)。

PS: を使用してすべてのNA/NaN値を検索するルートも調査しましx.array()==x.array()たが、Eigen でサブセット化するために結果を使用する方法を見つけることができませんでした。

4

2 に答える 2

12

forおそらく私は質問を正しく理解していませんが、Rcpp 内では、これをループ よりも効率的に行う方法がわかりません。forループは一般的に R では非効率的です。これは、R でループを反復処理するには多くの重いインタープリター機構が必要になるためです。しかし、C++ レベルに到達すると、これは当てはまりません。ネイティブにベクトル化された R 関数でさえ、最終的にはforC のループで実装されます。したがって、これをより効率的にする唯一の方法は、並列で実行することです。

たとえば、単一のベクトルから値na.omit.cppを省略する単純な関数を次に示します。NA

rcppfun<-"
Rcpp::NumericVector naomit(Rcpp::NumericVector x){
std::vector<double> r(x.size());
int k=0;
  for (int i = 0; i < x.size(); ++i) {
    if (x[i]==x[i]) {
    r[k] = x[i];
    k++;
   }
  }
 r.resize(k);
 return Rcpp::wrap(r);    
}"

na.omit.cpp<-cppFunction(rcppfun)

これは、R の組み込みよりもさらに高速に実行されna.omitます。

> set.seed(123)
> x<-1:10000
> x[sample(10000,1000)]<-NA
> y1<-na.omit(x)
> y2<-na.omit.cpp(x)
> all(y1==y2)
[1] TRUE
> require(microbenchmark)
> microbenchmark(na.omit(x),na.omit.cpp(x))
Unit: microseconds
           expr     min       lq   median      uq      max neval
     na.omit(x) 290.157 363.9935 376.4400 401.750 6547.447   100
 na.omit.cpp(x) 107.524 168.1955 173.6035 210.524  222.564   100
于 2013-10-03T12:06:31.240 に答える