[Rcpp-devel] speed
Two other pieces of information. I profiled the execution of fxR and it does seem to be spending its time in dbinom and dpois.
summaryRprof("Rprof.out")
$by.self
self.time self.pct total.time total.pct
"dbinom" 0.32 39.0 0.34 41.5
"dpois" 0.28 34.1 0.30 36.6
"fxR" 0.14 17.1 0.82 100.0
":" 0.04 4.9 0.04 4.9
"sum" 0.04 4.9 0.04 4.9
"str" 0.00 0.0 0.82 100.0
$by.total
total.time total.pct self.time self.pct
"fxR" 0.82 100.0 0.14 17.1
"str" 0.82 100.0 0.00 0.0
"dbinom" 0.34 41.5 0.32 39.0
"dpois" 0.30 36.6 0.28 34.1
":" 0.04 4.9 0.04 4.9
"sum" 0.04 4.9 0.04 4.9
$sampling.time
[1] 0.82
Like Dirk I was able to get about a 65-70% improvement in speed by
calling the R API functions Rf_dbinom and Rf_dpois and unrolling the
loops a bit (see enclosed). (Also my code seems to pass the test
although that makes me wonder why the fx version also passes the
test.)
I would interpret the information from the profiling above to indicate
that there are not big gains to be realized after this because so much
of the elapsed time will be tied up in the dbinom and dpois
evaluations.
On Wed, Aug 4, 2010 at 3:10 PM, Dirk Eddelbuettel <edd at debian.org> wrote:
Richard, Great question. One obvious consideration for improvement is to not call dbinom and dpois from R, but rather using the C entry points provided by Rmath.h. With that, I get these performances:
require(rbenchmark)
Loading required package: rbenchmark
benchmark(fx(0:100),
+ ? ? ? ? ? fx2(0:100), + ? ? ? ? ? fxR(0:100)) ? ? ? ?test replications elapsed relative user.self 1 ?fx(0:100) ? ? ? ? ?100 ?89.426 3.811850 ? ? 89.42 2 fx2(0:100) ? ? ? ? ?100 ?23.460 1.000000 ? ? 23.46 3 fxR(0:100) ? ? ? ? ?100 ?41.699 1.777451 ? ? 41.69
# deleted empty sys.self, user,child, sys.child columns here
The bad news is that I don't get the same results:
all.equal(fx(0:10), fxR(0:10))
[1] TRUE
all.equal(fx2(0:10), fxR(0:10))
[1] "Mean relative difference: 0.2535309"
But as Doug suggested, there may be a need for reworking things and checks
anyway. Or maybe I just introduced a bug -- dunno. I merely meant to help on
pointing out that Rmath.h is there too, and I did this as a quick and dirty
wrap around the 'atomistic' Rmath functions. Maybe someone would want to
contribute vectorised versions of these ? ?Patches welcome, as they say...
Dirk
inc2 <- '
#include <Rmath.h>
NumericVector wrap_dbinom(IntegerVector c, int k, double f) {
?NumericVector x(c.size());
?for (int i=0; i<c.size(); i++) x[i] = dbinom(c[i], k, f, false);
?return x;
}
NumericVector wrap_dpois(IntegerVector c, int l) {
?NumericVector x(c.size());
?for (int i=0; i<c.size(); i++) x[i] = dpois(c[i], l, false);
?return x;
}
'
fx2 <- cxxfunction(signature(Kr="integer"), '
? ?IntegerVector K(Kr);
? ?NumericMatrix bpsum(K.size(), K.size());
? ?for(int i=0; i<K.size(); i++) {
? ? ? ?for(int j=0; j<K.size(); j++) {
? ? ? ? ? ?IntegerVector Ki = K[i];
? ? ? ? ? ?IntegerVector cmin = seq_len(Ki.size()+1);
? ? ? ? ? ?IntegerVector cmin0 = cmin-1;
? ? ? ? ? ?NumericVector bin = wrap_dbinom(cmin0, K[j], 0.5);
? ? ? ? ? ?NumericVector pois = wrap_dpois(K[j]-cmin0, 1.5);
? ? ? ? ? ?NumericVector bp = bin * pois;
? ? ? ? ? ?bpsum(i, j) = std::accumulate(bp.begin(), bp.end(), 0.0);
? ? ? ? ? ?}
? ? ? ?}
? ?return bpsum;
? ?', plugin="Rcpp", include=inc2)
--
Dirk Eddelbuettel | edd at debian.org | http://dirk.eddelbuettel.com
_______________________________________________ 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 -------------- A non-text attachment was scrubbed... Name: foo.Rout Type: application/octet-stream Size: 2899 bytes Desc: not available URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20100804/50ec7c70/attachment-0001.obj>