[Rcpp-devel] RcppArmadillo -- crude benchmark of various mcmcvs rcppbugs
Whit Armstrong
armstrong.whit at gmail.com
Fri Apr 27 17:34:13 CEST 2012
The windows build was just accepted this morning by CRAN.
It should be appearing on the mirrors shortly.
It appears the mac builds are being held back because gcc 4.2 has some
template issues for modern c++ code. I'll try to find some #ifdef's
to work around this issue.
-Whit
On Fri, Apr 27, 2012 at 11:26 AM, Smith, Dale <Dale.Smith at fiserv.com> wrote:
> As the MCMCpack and rcppbugs packages were built on Linux, I can't install them on my Windows box. I'll have to get the sources and try to make them when I have time.
>
> Thanks,
> Dale Smith, Ph.D.
> Senior Financial Quantitative Analyst
> Risk & Compliance
> Fiserv.
> 107 Technology Park
> Norcross, GA 30092
> Office: 678-375-5315
> Mobile: 678-982-6599
> Mail: dale.smith at fiserv.com
> www.fiserv.com
>
>
> -----Original Message-----
> From: rcpp-devel-bounces at r-forge.wu-wien.ac.at [mailto:rcpp-devel-bounces at r-forge.wu-wien.ac.at] On Behalf Of Dirk Eddelbuettel
> Sent: Thursday, April 26, 2012 11:28 PM
> To: Whit Armstrong
> Cc: rcpp-devel
> Subject: Re: [Rcpp-devel] RcppArmadillo -- crude benchmark of various mcmcvs rcppbugs
>
>
> On 26 April 2012 at 15:58, Whit Armstrong wrote:
> | Slightly off topic, but Dirk kindly allowed me to make a quick
> | benchmark post on rcppbugs.
>
> "Allowed" ? Come on, benchmarks are always ok, as are posts on Rcpp. This is after all the Rcpp list...
> |
> | Just a short post to the list showing some performance numbers (for
> | running a simple linear regression) of rcppbugs (with an RcppArmadillo
> | core) vs MCMCpack (Scythe core) vs JAGS (?? core).
> |
> | I was fairly optimistic about the speed of rcppbugs, yet this exceeds
> | my expectations by so much that I fear there may be a bug in my
> | benchmark code.
> |
> | I'm gearing up for my presentation at R/Finance, so I would be
> | grateful to anyone who wants to kick the tires a bit. The core has
> | been tested (we've been using the pure c++ implementation for over a
> | year), but the R binding (via rcppbugs) is completely new.
> |
> | Thanks in advance. Code below (you'll need to install MCMCpack,
> | rjags, and rcppbugs (now on CRAN) to run it).
> |
> | The quick summary is a 5x speedup over MCMCpack and a 77x speedup (I'm
> | not sure I believe this yet) over JAGS.
>
> A bit more modest here (two-year old i7 running Ubuntu amd64)
>
> R> print(round(all.times.ratio,2))
> time ratio
> lm 0.09 0.08
> rcppbugs 1.18 1.00
> mcmcpack 4.05 3.45
> jags 57.06 48.56
>
> But still pretty good compared to MCMCpack. That is possibly legit -- modern
> C++ versus a C++ design from a decade ago. Jags is a little harder to
> believe. Maybe time to write an email to Martyn?
>
> Thanks for sharing this. Looking forward to your talk next month.
>
> Dirk
>
> | Comments welcome.
> |
> | -Whit
> |
> |
> | library(rcppbugs,quietly=TRUE,verbose=FALSE)
> | library(rjags,quietly=TRUE,verbose=FALSE)
> | library(MCMCpack,quietly=TRUE,verbose=FALSE)
> |
> | NR <- 1e2L
> | NC <- 2L
> |
> | y <- rnorm(NR,1) + 10
> | X <- matrix(nr=NR,nc=NC)
> | X[,1] <- 1
> | X[,2] <- y + rnorm(NR)/2 - 10
> |
> |
> | lm.time <- system.time(lm.res <- lm.fit(X,y))
> | ##print(coef(lm.res))
> |
> | b <- mcmc.normal(rnorm(NC),mu=0,tau=0.0001)
> | tau.y <- mcmc.uniform(sd(as.vector(y)),lower=0,upper=100)
> | ##y.hat <- deterministic(function(X,b) { X %*% b }, X, b) y.hat <-
> | linear(X,b) y.lik <- mcmc.normal(y,mu=y.hat,tau=tau.y,observed=TRUE)
> | m <- create.model(b, tau.y, y.hat, y.lik)
> |
> |
> | iterations <- 1e5L
> | burn <- iterations
> | adapt <- 1e3L
> | thin <- 10L
> |
> | cat("running rcppbugs model...\n")
> | rcppbugs.time <- system.time(ans <- run.model(m,
> | iterations=iterations, burn=burn, adapt=adapt, thin=thin))
> | rcppbugs.coefs <- apply(ans[["b"]],2,mean)
> | ##print(rcppbugs.coefs)
> | print(rcppbugs.time)
> |
> | cat("running MCMCpack model...\n")
> | mcmcpack.time <- system.time(mcmcpack.out <- MCMCregress(y ~
> | X.2,data=data.frame(y=y,X=X),burnin = burn, mcmc = iterations,
> | thin=thin))
> | print(mcmcpack.time)
> |
> | cat("running jags model...\n")
> |
> | bug.model <- '
> | model {
> | for (i in 1:NR){
> | y[i] ~ dnorm(y.hat[i], tau.y)
> | y.hat[i] <- inprod(b, X[i,])
> | }
> | for (j in 1:NC){
> | b[j] ~ dnorm(0, .0001)
> | }
> | tau.y ~ dunif(0, 100)
> | }
> | '
> | bug.file <- tempfile(fileext=".bug")
> | cat(bug.model,file=bug.file)
> |
> | jags.time1 <- system.time(
> | jags <- jags.model(bug.file,
> | data = list('X' = X,
> | 'y' = y,
> | 'NR' = NR,
> | 'NC' = NC),
> | n.chains = 1,
> | n.adapt = adapt,
> | quiet=TRUE
> | ))
> | jags.time2 <- system.time(update(jags, n.iter=burn,
> | progress.bar="none"))
> | jags.time3 <- system.time(jags.trace <-
> | jags.samples(jags,c('b','tau.y'),n.iter=iterations,thin=thin,progress.
> | bar="none"))
> |
> | jags.time <- jags.time1 + jags.time2 + jags.time3
> | print(jags.time)
> | jags.coefs <- apply(jags.trace$b,1,mean)
> |
> | ##print(jags.trace)
> |
> | cat("rcppbugs speedup:\n")
> | all.times <-
> | rbind(lm.time["elapsed"],rcppbugs.time["elapsed"],mcmcpack.time["elaps
> | ed"],jags.time["elapsed"]) all.times.ratio <-
> | cbind(all.times,all.times/rcppbugs.time["elapsed"])
> | colnames(all.times.ratio) <- c("time","ratio")
> | rownames(all.times.ratio) <- c("lm","rcppbugs","mcmcpack","jags")
> | print(round(all.times.ratio,2))
> |
> | cat("coef comparison:\n")
> | coef.compare <-
> | rbind(coef(lm.res),rcppbugs.coefs,summary(mcmcpack.out)$statistics[,"M
> | ean"][1:2],jags.coefs)
> | rownames(coef.compare) <- c("lm","rcppbugs","mcmcpack","jags")
> | print(coef.compare)
> | _______________________________________________
> | 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-deve
> | l
> --
> R/Finance 2012 Conference on May 11 and 12, 2012 at UIC in Chicago, IL See agenda, registration details and more at http://www.RinFinance.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
More information about the Rcpp-devel
mailing list