[R-sig-ME] Computational speed - MCMCglmm/lmer
Jarrod Hadfield
j.hadfield at ed.ac.uk
Thu Jun 24 10:33:33 CEST 2010
Hi Dave,
Sorry for the late response on this - I've been doing field work. I
think there is little option for speeding up the time per iteration in
MCMCglmm without considerable extra work. Most cpu are spent either
solving the mixed model equations or generating random (normal)
deviates. It may be possible to multi-thread the latter but how easy
this would be and how much of an increase in speed you would get I'm
not sure. For ZIP models at least, the long computing time is a
result of poor mixing. The next version (as you know) will allow
hurdle models to be fitted. These have much better mixing properties
and in terms of model fit are very similar to ZIP models.
Cheers,
Jarrod
On 19 Jun 2010, at 16:42, David Atkins wrote:
>
> Hi all--
>
> I use (g)lmer and MCMCglmm on a weekly basis, and I am wondering
> about options for speeding up their computations. This is primarily
> an issue with MCMCglmm, given the many necessary MCMC iterations to
> get to convergence on some problems. But, even with glmer(), I have
> runs that get into 20-30 minutes.
>
> First, let me be very clear that this is in no way a criticism of
> Doug's and Jarrod's work (package developers for lme4 and MCMCglmm,
> respectively). Their code has probably brought models/data into
> range that would not have been possible.
>
> Second, I have included link to data and script below, along with
> some timings on my computer: Mac Book Pro, 2.5GHz, with 4GB RAM.
> Here's sessionInfo from my runs:
>
> > sessionInfo()
> R version 2.11.1 (2010-05-31)
> i386-apple-darwin9.8.0
>
> locale:
> [1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
>
> attached base packages:
> [1] stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] MCMCglmm_2.05 corpcor_1.5.6 ape_2.5-3
> [4] coda_0.13-5 tensorA_0.35 lme4_0.999375-34
> [7] Matrix_0.999375-41 lattice_0.19-7
>
> loaded via a namespace (and not attached):
> [1] gee_4.13-14 grid_2.11.1 nlme_3.1-96 stats4_2.11.1
> tools_2.11.1
>
> Specific questions:
>
> 1. Would be curious to know timings on other people's set-ups.
> Jarrod and I had an exchange one time where he was gracious enough
> to run a zero-inflated model where I was concerned about
> convergence. He ran a model with 1.3M iterations, which I think
> would take a number of days, if not a week on my computer. This was
> part of what got me thinking about this. Thus, my first interest is
> whether there is an "optimal" hardware/OS configuration, or does it
> matter?
>
> Some things I see in R-help archives:
>
> 2. 32 vs. 64-bit: Seems like this is mostly an issue of data/model
> size and whether you need to access more than 4GB of RAM. AFAICS,
> 64-bit processors are not necessarily faster.
>
> 3. "Optimized" BLAS: There's a bit of discussion about optimized
> BLAS (basis linear algebra... something). However, these
> discussions note that there is no generally superior BLAS. Not sure
> whether specific BLAS might be optimized for GLMM computations.
>
> 4. Parallel computing: With multi-core computers, looks like there
> are some avenues for splitting intensive computations across
> processors.
>
> Finally, I'm asking here b/c I run into these issues with GLMM (and
> zero-inflated mixed models), though most of the discussion I've seen
> thus far about computation speed has been on R-help.
>
> The data below are self-reported drinks (alcohol) from college
> students for up to the last 90 days. Distribution of counts is zero-
> inflated. I run a Poisson GLMM with glmer, over-dispersed Poisson
> GLMM with MCMCglmm, and then zero-inflated OD Poisson with MCMCglmm
> and provide timings for my set-up.
>
> Okay, any and all thoughts welcomed.
>
> thanks, Dave
>
>
> ### thinking through computational speed with lmer and MCMCglmm
> #
> ### read drinking data
> drink.df <- read.table(file = "http://depts.washington.edu/cshrb/newweb/stats%20documents/tlfb.txt
> ",
> header = TRUE, sep = "\t")
> str(drink.df) # 57K rows
>
> ### id is id variable (shocking)
> ### gender and weekday are 0/1 indicator variables
> ### drinks has number of drinks consumed
>
> ### distribution of outcome
> table(drink.df$drinks)
> plot(table(drink.df$drinks), lwd=2) # zero-inflated
>
> ### how many people?
> length(unique(drink.df$id)) # 990
> sort(table(drink.df$id))
>
> ### NOTE: most people have max of 90, which they should
> ### two folks with 180 and 435 (prob data errors)
> ### long negative tail down from 90
> #
> ### NOTE: for this purpose, not worrying about the 180/435
>
> ### speed tests
> #
> ### fit random intercept and random slope for weekday
> ### fixed effects for gender and weekday and interaction
> #
> ### Poisson GLMM with glmer()
> library(lme4)
>
> system.time(
> drk.glmer <- glmer(drinks ~ weekday*gender + (weekday | id),
> data = drink.df, family = poisson,
> verbose = TRUE)
> )
> summary(drk.glmer)
>
> ### timing
> #
> ### user system elapsed
> ### 36.326 9.013 45.316
>
> ### over-dispersed Poisson GLMM with MCMCglmm()
> library(MCMCglmm)
>
> prior <- list(R = list(V = 1, n = 1),
> G = list(G1 = list(V = diag(2), n = 2)))
> system.time(
> drk.mcmc <- MCMCglmm(drinks ~ weekday*gender,
> random = ~ us(1 + weekday):id,
> data = drink.df, family = "poisson",
> prior = prior)
> )
> summary(drk.mcmc) # NOTE: using summary.MCMCglmm in v 2.05 of package
>
> ### timing
> #
> ### user system elapsed
> ### 1034.317 165.128 1203.536
>
> ### zero-inflated, over-dispersed Poisson GLMM with MCMCglmm()
> #
> ### NOTE: haven't run the following yet, other than a quick "toy
> run" to
> ### sure that it is set up correctly.
> ### NOTE: this only has random intercept in each portion of the model
>
> prior2 <- list(R = list(V = diag(2), n = 1, fix = 2),
> G = list(G1 = list(V = 1, n = 1),
> G2 = list(V = 1, n = 1)))
> system.time(
> drk.zimcmc <- MCMCglmm(drinks ~ -1 + trait*weekday*gender,
> random = ~ idh(at.level(trait, 1)):id + idh(at.level(trait,
> 2)):id,
> rcov = ~ idh(trait):units,
> data = drink.df, family = "zipoisson",
> prior = prior2)
> )
> summary(drk.zimcmc)
>
> ### timing
> #
> ### user system elapsed
> ### 2105.366 544.881 2640.030
>
> --
> Dave Atkins, PhD
> Research Associate Professor
> Department of Psychiatry and Behavioral Science
> University of Washington
> datkins at u.washington.edu
>
> Center for the Study of Health and Risk Behaviors (CSHRB)
> 1100 NE 45th Street, Suite 300
> Seattle, WA 98105
> 206-616-3879
> http://depts.washington.edu/cshrb/
> (Mon-Wed)
>
> Center for Healthcare Improvement, for Addictions, Mental Illness,
> Medically Vulnerable Populations (CHAMMP)
> 325 9th Avenue, 2HH-15
> Box 359911
> Seattle, WA 98104?
> 206-897-4210
> http://www.chammp.org
> (Thurs)
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.
More information about the R-sig-mixed-models
mailing list