Rcpp Rcpp에서 임의의 다 변수 법선 벡터를 생성하는 방법은 무엇입니까?

Question

큰 무작위 다변량 (6 차원 이상) 정상 샘플을 생성하고 싶습니다. R에서는 많은 패키지가 rmnorm, rmvn과 같이 이것을 할 수 있습니다 ... 그러나 문제는 속도입니다! 그래서 Rcpp를 통해 C 코드를 작성하려고했습니다. 온라인으로 튜토리얼을 살펴 보았지만 STL 라이브러리가 아닌 다 변수 배포를위한 "설탕"이 없다.

도움을 주시면 감사하겠습니다.

감사합니다.

Answer 1

다 변수 (cholesky, svd 등)를 근사화하고 Eigen (RccpEigen) 또는 Armadillo (RcppArmadillo 사용)를 사용하여 프로그래밍하는 좋은 알고리즘을 찾지 못하면 Rcpp가 도움이되는지 확신 할 수 없습니다.

가 여기에 콜레의 분해를 사용하여 하나의 방법이다 (Rcpp) 아르마 이제

#include <RcppArmadillo.h> 

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

// [[Rcpp::export]] 

using namespace arma; 
using namespace Rcpp; 

mat mvrnormArma(int n, mat sigma) { 
    int ncols = sigma.n_cols; 
    mat Y = randn(n, ncols); 
    return Y * chol(sigma); 
} 

everythings이

sigma <- matrix(c(1, 0.9, -0.3, 0.9, 1, -0.4, -0.3, -0.4, 1), ncol = 3) 
cor(mvrnormR(100, sigma)) 
cor(MASS::mvrnorm(100, mu = rep(0, 3), sigma)) 
cor(mvrnormArma(100, sigma)) 

를 작업하는 경우 또한 확인할 수 있습니다 순수 R

mvrnormR <- function(n, sigma) { 
    ncols <- ncol(sigma) 
    matrix(rnorm(n * ncols), ncol = ncols) %*% chol(sigma) 
} 

의 원래 구현

이제 벤치마킹하자

require(bencharmk) 
benchmark(mvrnormR(1e4, sigma), 
      MASS::mvrnorm(1e4, mu = rep(0, 3), sigma), 
      mvrnormArma(1e4, sigma), 
      columns=c('test', 'replications', 'relative', 'elapsed')) 


## 2 MASS::mvrnorm(10000, mu = rep(0, 3), sigma)   100 
## 3     mvrnormArma(10000, sigma)   100 
## 1      mvrnormR(10000, sigma)   100 
## relative elapsed 
## 2 3.135 2.295 
## 3 1.000 0.732 
## 1 1.807 1.323 

이 예제에서는 단위 분산 및 널 평균을 사용하는 정규 분포를 사용했지만 사용자 정의 평균 및 분산을 사용하여 쉽게 가우스 분포로 일반화 할 수있었습니다.

희망이 있습니다.