임의의 상태를 RcppEigen을 사용하여 setRandom에 전달하십시오.

Question

임의의 상태를 Eigen의 setRandom에 RcppEigen으로 전달하는 방법이 있습니까? 아니면 runif을 사용해야합니까? 여기

은 예이다 :

// [[Rcpp::depends(RcppEigen)]] 
#include <RcppEigen.h> 
using namespace Rcpp; 
using Eigen::MatrixXd; 
using Eigen::VectorXd; 

// [[Rcpp::export]] 
NumericVector fx() { 
    RNGScope scope; 
    MatrixXd x(3,2); 
    x=x.setRandom(); 
    x.col(1)=as<VectorXd>(runif(3,0,1)); 

return wrap(x); 
} 

테스트 그것을 : 2 열 (즉, runif)을 재현하는 방법

set.seed(42); fx() 
#   [,1]  [,2] 
#[1,] -0.8105760 0.9148060 
#[2,] 0.6498853 0.9370754 
#[3,] 0.6221027 0.2861395 

set.seed(42); fx() 
#   [,1]  [,2] 
#[1,] -0.9449154 0.9148060 
#[2,] 0.8063267 0.9370754 
#[3,] -0.0673205 0.2861395 

참고하지만, 1 열 (즉, setRandom)이 아니다.

Answer 1

Eigen의 RNG는 R의 RNG와 직교합니다. 우리는 RNGScope을 가지고 R을 처리하고 알고있는 것처럼 set.seed()을 허용합니다. 아이겐의 RNG로하는 일은 별개입니다.

특히 Eigen의 RNG를 R에 대한 사용자 제공 RNG로 등록해야합니다. 기본적으로 R은 R 사용자가 R의 RNG를 기대할 때 Eigen을 인식하지 못합니다.

여기에 문제는 Eigen RNG를 초기화하지 못한 바닐라 'Eigen'프로그래밍 문제인 것 같습니다. Rcpp와는 아무런 관련이 없습니다.

편집 : 다음은 수정 된 예입니다. Eigen은 시스템 RNG를 사용하기 때문에 시드를 사용하려면 srand()이 필요합니다. 영업 이익의 의견에 응답

// [[Rcpp::depends(RcppEigen)]] 
#include <RcppEigen.h> 
using namespace Rcpp; 
using Eigen::MatrixXd; 
using Eigen::VectorXd; 

// [[Rcpp::export]] 

    RNGScope scope; 
    MatrixXd x(3,2); 
    srand(seed); 
    x=x.setRandom(); 
    x.col(1)=as<VectorXd>(runif(3,0,1)); 

    return wrap(x); 
} 

편집 : 2 :

R> sourceCpp("/tmp/roland.cpp") 
R> set.seed(42); fx(42) 
      [,1]  [,2] 
[1,] -0.933060 0.914806 
[2,] -0.340072 0.937075 
[3,] 0.381271 0.286140 
R> set.seed(42); fx(42) 
      [,1]  [,2] 
[1,] -0.933060 0.914806 
[2,] -0.340072 0.937075 
[3,] 0.381271 0.286140 
R> 

그것은 코드의 수정 된 버전을 사용

내가 Rcpp::runif()와 사이의 속도 차이의 대부분을 기대하지 않을 것이다

srand() Eigen에서 만든 호출이므로 srand()에 문제가 발생했으며 시스템간에 다르게 동작 할 수 있습니다.

빠른 데모 스크립트 :

이

#include <RcppEigen.h> 

// [[Rcpp::depends(RcppEigen)]] 

// [[Rcpp::export]] 
Rcpp::NumericVector v1(int n) { 
    return Rcpp::runif(n); 
} 

// [[Rcpp::export]] 
Rcpp::NumericVector v2(int n) { 
    Eigen::VectorXd x(n); 
    x = x.setRandom(); 
    return Rcpp::wrap(x); 
} 

/*** R 
library(rbenchmark) 
N <- 1e7 
benchmark(v1(N), v2(N)) 
*/ 

R> sourceCpp("/tmp/roland.cpp") 

R> library(rbenchmark) 

R> N <- 1e7 

R> benchmark(v1(N), v2(N)) 
    test replications elapsed relative user.self sys.self user.child sys.child 
1 v1(N)   100 12.633 1.000 11.356 1.261   0   0 
2 v2(N)   100 17.222 1.363 13.981 3.198   0   0 
R> 

을 생산하고 RcppEigen도 단지 만들기 벡터의 간단 설치 여기에 느린입니다 있습니다. 하지만 여기서 마이크로 초를 말하고 있습니다. 아마도 다른 어떤 병목 현상이있을 가능성이 높은 실제 응용 프로그램에서 걱정하지 않아도됩니다.