[Rcpp-devel] Suggestions for possible optimizations

[Alessandro Mammana]

Dear Rcpp developers,

I am optimizing some R code by replacing some loops with Rcpp code.
After some profiling, it looks like the cpp function below is still
the limiting step. Does anybody see some possible source of
inefficiency or something easy to improve? Or is it already as fast as
it can get? Thanks a lot.

#include <Rcpp.h>

using namespace Rcpp;

// [[Rcpp::export]]
List llik2posteriors(NumericMatrix lliks){
    NumericMatrix posteriors(lliks.nrow(), lliks.ncol());
    double llik = 0;
    double cmax;
    double psum;
    for (int i = 0; i < lliks.ncol(); ++i){
        MatrixColumn<REALSXP> llikrow = lliks.column(i);
        MatrixColumn<REALSXP> postrow = posteriors.column(i);
        cmax = max(llikrow);
        llik += cmax;
        postrow = exp(llikrow - cmax);
        psum = sum(postrow);
        llik += log(psum);
        postrow = postrow/psum;
    }

    return List::create(Named("posteriors")=posteriors, Named("tot_llik")=llik);
}

Alessandro