[Rcpp-devel] How to increase the coding efficiency
Yan Zhou
zhouyan at me.com
Tue Dec 4 21:39:23 CET 2012
You cannot expect C++ to magically make your code faster. If speed is of concern, you need profiling to find the bottleneck instead of blind guessing. I am not quite sure how to profile an Rcpp R program. Maybe someone else can help.
But my intuition is that your program simply need that much time. The reason that R is not much slower is perhaps the linear algebra operations took much of the time, and that is as efficient in R as in C++. If you built R with an optimized Blas, the difference may be even smaller.
The bottom line is that C++ is no magic. The reason that it can be an order of magnitude faster than R sometime is that the close to metal natural enable one to take out some middleware. The computation that requires a bulk of time still need the same time whatever language you use.
And even there is something to improve in coding efficiency, that is not something you can learn through just a few emails. This kind if thing take time, practice and experience.
Yan
On Dec 4, 2012, at 8:27 PM, Honglang Wang <wanghonglang2008 at gmail.com> wrote:
> Since the "ker" function defined in R is very simple and vectorized, I thought it was fine to call this function in C++. And I tried to code this function inside C++, and I found it became ever slower: from 47.310s to 49.536s.
>
> // [[Rcpp::depends(RcppArmadillo)]]
>
> #include <RcppArmadillo.h>
>
> using namespace Rcpp;
>
> // [[Rcpp::export]]
> List betahat(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 tp = t-t0;
> arma::uvec q1 = arma::find(arma::abs(tp)<h);
> arma::uvec q2 = arma::find(arma::abs(tp)>=h);
> arma::vec kk = t;
> kk.elem(q1) = ((1-arma::pow(tp.elem(q1)/h,2))/h)*0.75;
> kk.elem(q2).zeros();
> 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);
> }
>
>
> Best wishes!
>
> Honglang Wang
>
> Office C402 Wells Hall
> Department of Statistics and Probability
> Michigan State University
> 1579 I Spartan Village, East Lansing, MI 48823
> wangho16 at msu.edu
>
>
>
> On Tue, Dec 4, 2012 at 2:19 PM, Romain Francois <romain at r-enthusiasts.com> wrote:
>> 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
>>> wangho16 at msu.edu <mailto:wangho16 at msu.edu>
>>
>>
>> --
>> 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
>> Rcpp-devel at lists.r-forge.r-project.org
>> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
>
> _______________________________________________
> 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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20121204/1e4402da/attachment-0001.html>
More information about the Rcpp-devel
mailing list