1 2library(Rcpp) 3 4cppFunction(' 5 NumericVector convolveCpp(NumericVector a, NumericVector b) { 6 7 int na = a.size(), nb = b.size(); 8 int nab = na + nb - 1; 9 NumericVector xab(nab); 10 11 for (int i = 0; i < na; i++) 12 for (int j = 0; j < nb; j++) 13 xab[i + j] += a[i] * b[j]; 14 15 return xab; 16 } 17') 18 19convolveCpp(c(1,2,3), matrix(3,3)) 20 21 22cppFunction(depends='RcppArmadillo', code=' 23 List fastLm(NumericVector yr, NumericMatrix Xr) { 24 25 int n = Xr.nrow(), k = Xr.ncol(); 26 27 arma::mat X(Xr.begin(), n, k, false); // reuses memory and avoids copy 28 arma::colvec y(yr.begin(), yr.size(), false); 29 30 arma::colvec coef = arma::solve(X, y); // fit model y ~ X 31 arma::colvec resid = y - X*coef; // residuals 32 33 double sig2 = arma::as_scalar( arma::trans(resid)*resid/(n-k) ); 34 // std.error of estimate 35 arma::colvec stderrest = arma::sqrt( 36 sig2 * arma::diagvec( arma::inv(arma::trans(X)*X)) ); 37 38 return List::create(Named("coefficients") = coef, 39 Named("stderr") = stderrest 40 ); 41} 42') 43 44