You are calling back to R using Function. This is expensive.
What is ker ? Can you implement it in C++. This is a wild guess, but
that is I think where the bottleneck is.
Romain
Le 04/12/12 20:14, Honglang Wang a écrit :
Dear All,
I have tried out the first example by using RcppArmadillo, but I am not
sure whether the code is efficient or not. And I did the comparison of
the computation time.
1) R code using for loop in R: 87.22s
2) R code using apply: 77.86s
3) RcppArmadillo by using for loop in C++: 53.102s
4) RcppArmadillo together with apply in R: 47.310s
It is kind of not so big increase. I am wondering whether I used an
inefficient way for the C++ coding:
// [[Rcpp::depends(RcppArmadillo)]]
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::export]]
List betahat(Function ker, double t0, NumericMatrix Xr, NumericMatrix
yr, NumericVector tr, double h, int m) {
int n = Xr.nrow(), p = Xr.ncol();
arma::mat X(Xr.begin(), n, p, false);
arma::mat y(yr.begin(), n, 1, false);
arma::colvec t(tr.begin(), tr.size(), false);
arma::mat T = X;
T.each_col() %= (t-t0)/h;
arma::mat D = arma::join_rows(X,T);
arma::vec kk =as<arma::vec>(ker(tr-t0,h));
arma::mat W = (arma::diagmat(kk))/m;
arma::mat Inv = arma::trans(D)*W*D;
arma::vec betahat = arma::inv(Inv)*arma::trans(D)*W*y;
arma::colvec betahat0(betahat.begin(),betahat.size()/2,false);
return List::create(Named("betahat") = betahat0);
}
Anyone could help me with how to increase the efficiency of the coding?
Thanks.
Best wishes!
Honglang Wang
Office C402 Wells Hall
Department of Statistics and Probability
Michigan State University
1579 I Spartan Village, East Lansing, MI 48823
[email protected] <mailto:[email protected]>
--
Romain Francois
Professional R Enthusiast
+33(0) 6 28 91 30 30
R Graph Gallery: http://gallery.r-enthusiasts.com
`- http://bit.ly/SweN1Z : SuperStorm Sandy
blog: http://romainfrancois.blog.free.fr
|- http://bit.ly/RE6sYH : OOP with Rcpp modules
`- http://bit.ly/Thw7IK : Rcpp modules more flexible
_______________________________________________
Rcpp-devel mailing list
[email protected]
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
