[R] how to obtain EM-estimates of cov(b) and var(e) from lme

Douglas Bates bates at stat.wisc.edu
Sat Dec 8 23:21:02 CET 2001


Han Lai <hlai at whsun1.wh.whoi.edu> writes:

> I have a simple random-coefficients model for m subjects:
> 
>     y = b0 + b1 x + r0 + r1  x + e
> 
> where b0 and b1 are fixed parameters, r0 and r1 are random,
> e ~ N(0,s2 I) and R' = [r0, r1] ~ N(0,T).
> 
> I try to obtain the EM-estimates of s2 and the elements of T by
> 
>     lme(y~x,data=mydata,random= list(group=~x),
>           control=lmeControl(maxIter = 0, niterEM=100,msVerbose = TRUE))
> 
> Does this statement do the job?

The expression "EM-estimates of cov(b) and var(e)" doesn't make sense.
The EM algorithm is one iterative technique to determine the maximum
likelihood (ML) or restricted maximum likelihood (REML) estimates in a
mixed-effects model.  The estimates are the ML or the REML estimates
regardless of whether they were determined by EM iterations or by
other methods.

The EM algorithm tends to be robust to starting values but converges
very slowly near the optimimum.  Convergence is so slow that it is
difficult to determine if the EM algorithm has converged which is why,
in lme, we use EM iterations followed by a general optimization
procedure.  We always check for convergence in the general
optimization procedure.  Setting maxIter = 0 may result in convergence
failure, in which case the fitted model object is not constructed and
you won't be able to apply methods to it.

If for some reason you want to check that the EM iterations converge
you could use what you have above but remove maxIter = 0.

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._



More information about the R-help mailing list