[Rcpp-devel] Best way to create random normal matrix

Dirk Eddelbuettel edd at debian.org
Wed Dec 5 01:03:12 CET 2012


On 4 December 2012 at 15:40, Julian Smith wrote:
| I've been toying around with RcppEigen and I was curious what was the best way
| to create a random matrix of standard normals of size n by n. I have tried in

I'd do it the same way it is done in R -- look eg at MASS::mvrnorm. I assume
you also have a corr.matrix of size N * N (you didn't say...) so draw a
matrix of size N * N of rnorm() and multiply by a suitably decomposition of
the correlation matrix.  

These algorithms are described in a number of places (incl Wikipedia), and
MASS::mvrnorm() can serve as a reference.

And from your error below, something may not be right with your setup of
inline or your compiler -- which we cannot debug from here without extra
information.  

So you may have to step back and make sure you can actually build simpler
programs.

Dirk


| Armadillo but ran into an error. Would prefer to do in RcppEigen.
| 
| (Note: I'm not 100% sure why I get the error with Armadillo either)
| 
| src <- '
| 
| Rcpp::NumericMatrix Xr(Xs);
| 
| int q = Rcpp::as<int>(ys);
| 
| int n = Xr.nrow(), k = Xr.ncol();
| 
| arma::mat X(Xr.begin(), n, k, false);
| 
| arma::mat G, Y, B;
| 
| G = arma::randn(n,q);
| 
| Y = X*G;
| 
| arma::mat Q, R;
| 
| arma::qr(Q,R,Y);
| 
| Q = Q(arma::span::all,arma::span(0,q-1));
| 
| B = arma::trans(Q)*X;
| 
| arma::vec d;
| 
| arma::mat u, v;
| 
| arma::svd(u, d, v, B );
| 
| arma::mat U = Q*u;
| 
| return Rcpp::List::create(Rcpp::Named("d")=d,Rcpp::Named("U")=U,Rcpp::Named
| ("V")=v);'
| 
| rsvd <- cxxfunction(signature(Xs="numeric", ys="integer"), body=src, plugin=
| "RcppArmadillo")
| 
| 
| appsvd <- function(X, q) {
| 
|     .Call("rsvd", X, q, package = "RcppArmadillo")
| 
| }
| 
| 
| [1] Error in .Call("rsvd", X, q, package = "RcppArmadillo") : C symbol name
| "rsvd" not in load table
| 
| 
| 
| ----------------------------------------------------------------------
| _______________________________________________
| Rcpp-devel mailing list
| Rcpp-devel at lists.r-forge.r-project.org
| https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
-- 
Dirk Eddelbuettel | edd at debian.org | http://dirk.eddelbuettel.com  


More information about the Rcpp-devel mailing list