[R-sig-ME] Numerical integration for cross-classified random effects in lme4

[David Atkins]

Jeremy--

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 
good.

For what it's worth...

cheers, Dave
-- 
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
http://www.chammp.org
(Thurs)

Jeremy wrote:

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?

Thanks again.

--- 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 
factor

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
prohibitive.

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
random effects.




	[[alternative HTML version deleted]]