Re: [R] generalized linear mixed model by ML

From: Douglas Bates <dmbates_at_gmail.com>
Date: Fri 16 Dec 2005 - 04:51:49 EST

On 12/15/05, Ben Bolker <bolker@ufl.edu> wrote:
> Abderrahim Oulhaj <abderrahim.oulhaj <at> pharmacology.oxford.ac.uk> writes:
>
> >
> > Dear All,
> >
> > I wonder if there is a way to fit a generalized linear mixed models (for
> repeated binomial data) via a direct
> > Maximum Likelihood Approach. The "glmm" in the "repeated" package (Lindsey),
> the "glmmPQL" in the
> > "MASS" package (Ripley) and "glmmGIBBS" (Myle and Calyton) are not using the
> full maximum likelihood as I
> > understand. The "glmmML" of Brostrom uses the "full maximum likelihood" by
> approximating the integral
> > via Gauss- Hermite quadrature. However, glmmML is only valid for the random
> intercept model and the
> > binomial family must be represented only as binary data. Does the lmer do the
> work?

The expression "full maximum likelihood" is confusing. In the multilevel modeling literature it is sometimes used in contrast to restricted maximum likelihood or REML estimation but in that context it represents what statisticians would simply term "maximum likelihood".

I think you are referring to "exact maximum likelihood" where the criterion being optimized is the exact log-likelihood associated with the parameters. I don't believe that can be done in the the evaluation of the (marginal) log-likelihood associated with the parameters to be estimated involves an intractable integral. Various methods for estimating parameters in this model involve approximations to this integral. Even the Bayesian methods such as MCMC methods are using an approximation to this integral (in that case a stochastic approximation).

So, as far as I know, there are no exact maximum likelihood estimation methods for this model except in some special circumstances.
>
> Hmmm. I will be interested to hear what others have to say on
> this topic.
>
> * lmer() in the lme4 package (new version of nlme) can in
> fact do GLMMs with a choice of different
> integration methods (PQL is the default but not the only choice).
>
> * GLMMGibbs [sic] was actually
> using a full likelihood rather than an approximation, but was
> a Bayesian rather than a ML approach [GLMMGibbs is now in the
> "Devel" section of CRAN, apparently because of various unresolved
> compilation/installation problems.]
>
> If you want to fit temporal correlation, as well as
> individual random effects, you may be out of luck: a GEE
> model will probably be your best bet in that case. When I asked
> about the possibility of incorporating temporal/spatial correlation
> structures like those in nlme into lme4, Doug Bates said that he
> wanted to work first on getting the basic framework of the package
> really solid [can't blame him at all, and of course honor and
> glory to him for putting so much work into these tools in the first place]
>

Wow - honor and glory - thanks.

The lmer function in the currently released Matrix package (0.99-3) does provide PQL, Laplace and Adaptive Gauss-Hermite Quadrature (AGQ) methods for estimation of the parameters in a generalized linear mixed model. PQL is the fastest but least accurate method. The Laplacian approximation is actually a special case of the AGQ method where only one quadrature point per random effect is used. The AGQ method is the most accurate and in many ways is the "gold standard". It has a tuning parameter which is the number of quadrature points per random effect. Usually this is an odd number (because you always evaluate at the conditional modes of the random effects which, for odd numbers, is one of the quadrature points). The lmer function allows up to 11 quadrature points per random effect. It would be easier to add the capability for more but adding a whole lot more quadrature points doesn't buy you that much. It is definitely a case of diminishing returns.

The released version of the Matrix package uses a partitioned representation of the mixed model structures and sometimes encounters errors with complex models based on multiple grouping factors that are non-nested. A couple of months ago I started developing an alternative version (for those keeping count, yet this is the sixth implementation - the one I said I wouldn't do) based on what is called a supernodal Cholesky factorization for which Tim Davis recently released code. That is doing fine on the linear mixed models and PQL for the generalized linear mixed model. I'm adding Laplace now and feel I am close to getting that working (although, as Dave Balsinger said, "Programmers spend 95% of their time being `95% done'"). AGQ is more difficult especially for multiple grouping factors. At present AGQ is only available for a univariate random effect associated with a single grouping factor. Once I get that implemented I will release a new Matrix package based on the supernodal factorization.

> good luck,
> Ben Bolker
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Fri Dec 16 05:24:14 2005

This archive was generated by hypermail 2.1.8 : Fri 16 Dec 2005 - 09:45:39 EST