[R-sig-ME] Numerical integration for cross-classified random effects in lme4
datkins at u.washington.edu
Wed Oct 27 21:39:35 CEST 2010
Just to chime in here, I think the devil is almost always in the details
of the specific data and model. Thus, "Stata can 'do' AGQ for multiple
random-effects" really does not mean too much outside of a specific set
of data and model.
I recently assisted with an Epidemiology colleague who was trying to run
a cross-classified GLMM with logit outcome and ~45K participants. After
24 hours Stata had made it through 4-5 iterations. lmer() fit the model
in 5 minutes.
At the same time, that example does not mean lmer() is universally
superior for all models and datasets. However, if you're talking
cross-classified data... it probably is. ;)
As for MCMC, I would strongly recommend taking a look at MCMCglmm, which
is a fully Bayesian package for generalized linear mixed models and very
For what it's worth...
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
Center for Healthcare Improvement, for Addictions, Mental Illness,
Medically Vulnerable Populations (CHAMMP)
325 9th Avenue, 2HH-15
Seattle, WA 98104
Thanks for the clarification, Doug.
I asked a colleague about this the other day, and she indicated that it
is possible to use numerical integration for multiple random effects in
Stata . . . but she wasn't sure if that would apply to cross-classified
models or perhaps just hierarchically nested random effects.
On a related note, I'd be interested in using MCMC estimation, and I
remember reading a while ago that you were thinking about integrating
that option in lme4. Since then, I've seen some evidence that there's
an mcmcsamp option, but I haven't been able to find any representative
scripts so far. Do you know of any examples that would be good to follow?
--- On Mon, 10/25/10, Douglas Bates <bates at stat.wisc.edu> wrote:
From: Douglas Bates <bates at stat.wisc.edu>
Subject: Re: [R-sig-ME] Numerical integration for cross-classified
random effects in lme4
To: "Jeremy Koster" <helixed2 at yahoo.com>
Cc: r-sig-mixed-models at r-project.org
Date: Monday, October 25, 2010, 2:49 PM
On Sun, Oct 24, 2010 at 6:46 PM, Jeremy Koster <helixed2 at yahoo.com> wrote:
> Hi folks,
> This has likely been asked before, so please feel free to link any
relevant posts . . .
> I am looking at a network of 300+ households. The outcome variable
is dichotomous -- i.e., does Household A share food with Household B?
Because food can flow both ways, one row in the dataset is for Household
A to Household B, and then there's another row for B to A. So I'd like
to add a random effect for the dyad (DyadID). In addition, some
households just generally give more, and others generally receive more,
so I'd like to add random effects for donating households (GNO) and
receiving households (RNO).
> My impression was that it's possible to fit cross-classified random
effects using numerical integration in lme4. However, when I try to fit
an empty model using this code:
> model.empty <- glmer (Sharing ~ 1 + (1|DyadID) + (1|GNO) +
(1|RNO),family = binomial, nAGQ=10)
> Then I get the following error message:
> Error in validObject(.Object) :
> invalid class "mer" object: AGQ method requires a single grouping
As the message indicates, adaptive Gauss-Hermite quadrature is only
available for models with random effects defined with respect to a
single grouping factor.
When there is only one grouping factor the observations can be split
according to the levels of the grouping factor and the integral
defining the likelihood of the parameters can be expressed as the
product of a number of low-dimensional integrals. You need low
dimensional, preferably one-dimensional, integrals before you can hope
to apply AGQ. For high-dimensional integrals that number of
evaluations of the conditional mean that would be required for a
single evaluation of the likelihood of the parameters would be
The Laplace approximation, which does require optimization of the
unscaled conditional density, but only requires evaluation of the
conditional mean at that point, is feasible for models with crossed
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models