[Rcpp-commits] r4552 - in pkg/RcppArmadillo: . inst/examples

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Sep 30 03:05:26 CEST 2013


Author: edd
Date: 2013-09-30 03:05:25 +0200 (Mon, 30 Sep 2013)
New Revision: 4552

Modified:
   pkg/RcppArmadillo/ChangeLog
   pkg/RcppArmadillo/inst/examples/fastLm.r
   pkg/RcppArmadillo/inst/examples/varSimulation.r
Log:
use 'const ref' in two example files


Modified: pkg/RcppArmadillo/ChangeLog
===================================================================
--- pkg/RcppArmadillo/ChangeLog	2013-09-30 01:05:03 UTC (rev 4551)
+++ pkg/RcppArmadillo/ChangeLog	2013-09-30 01:05:25 UTC (rev 4552)
@@ -1,3 +1,9 @@
+2013-09-29  Dirk Eddelbuettel  <edd at debian.org>
+
+	* inst/examples/fastLm.r: Added a 'const ref' example
+	* inst/examples/varSimulation.r: Use 'const ref', skip one
+	superfluous transpose
+
 2013-09-28  Dirk Eddelbuettel  <edd at debian.org>
 
 	* DESCRIPTION: Release 0.3.920.1

Modified: pkg/RcppArmadillo/inst/examples/fastLm.r
===================================================================
--- pkg/RcppArmadillo/inst/examples/fastLm.r	2013-09-30 01:05:03 UTC (rev 4551)
+++ pkg/RcppArmadillo/inst/examples/fastLm.r	2013-09-30 01:05:25 UTC (rev 4552)
@@ -67,7 +67,28 @@
 '
 cppFunction(code=src, depends="RcppArmadillo")
 
+src <- '
+Rcpp::List fLmConstRef(const arma::mat & X, const arma::colvec & y) {
+    int df = X.n_rows - X.n_cols;
 
+    // fit model y ~ X, extract residuals
+    arma::colvec coef = arma::solve(X, y);
+    arma::colvec res  = y - X*coef;
+
+    double s2 = std::inner_product(res.begin(), res.end(),
+                                   res.begin(), 0.0)/df;
+    // std.errors of coefficients
+    arma::colvec sderr = arma::sqrt(s2 *
+       arma::diagvec(arma::pinv(arma::trans(X)*X)));
+
+    return Rcpp::List::create(Rcpp::Named("coefficients")=coef,
+                              Rcpp::Named("stderr")      =sderr,
+                              Rcpp::Named("df")          =df);
+}
+'
+cppFunction(code=src, depends="RcppArmadillo")
+
+
 fastLmPureDotCall <- function(X, y) {
     .Call("fastLm", X, y, PACKAGE = "RcppArmadillo")
 }
@@ -79,6 +100,7 @@
 
 res <- benchmark(fLmOneCast(X, y),             	# inline'd above
                  fLmTwoCasts(X, y),            	# inline'd above
+                 fLmConstRef(X, y),            	# inline'd above
                  fastLmPure(X, y),              # similar, but with 2 error checks
                  fastLmPureDotCall(X, y),       # now without the 2 error checks
                  fastLm(frm, data=trees),       # using model matrix
@@ -89,4 +111,4 @@
                  order="relative",
                  replications=5000)
 
-print(res)
+print(res[,1:4])

Modified: pkg/RcppArmadillo/inst/examples/varSimulation.r
===================================================================
--- pkg/RcppArmadillo/inst/examples/varSimulation.r	2013-09-30 01:05:03 UTC (rev 4551)
+++ pkg/RcppArmadillo/inst/examples/varSimulation.r	2013-09-30 01:05:25 UTC (rev 4552)
@@ -52,13 +52,13 @@
 
 ## C++ variant: code passed as a text variable ...
 code <- '
-arma::mat rcppSim(arma::mat coeff, arma::mat errors) {
+arma::mat rcppSim(const arma::mat& coeff, const arma::mat& errors) {
     int m = errors.n_rows;
     int n = errors.n_cols;
     arma::mat simdata(m,n);
     simdata.row(0) = arma::zeros<arma::mat>(1,n);
     for (int row=1; row<m; row++) {
-        simdata.row(row) = simdata.row(row-1)*trans(coeff)+errors.row(row);
+        simdata.row(row) = simdata.row(row-1) * coeff + errors.row(row);
     }
     return simdata;
 }



More information about the Rcpp-commits mailing list