[R-sig-ME] lme varFix under ML fit does not match coefficients standard error

Dimitris Rizopoulos d@r|zopou|o@ @end|ng |rom er@@mu@mc@n|
Tue Feb 20 21:40:32 CET 2024


In principle, the standard errors should be different with REML and ML. When I try different datasets, I see differences. For example,

fm_reml <- lme(distance ~ age, data = Orthodont, random = ~ age | Subject)
fm_ml <- lme(distance ~ age, data = Orthodont, random = ~ age | Subject,
             method = "ML")

coef(summary(fm_reml))
coef(summary(fm_ml))

Best,
Dimitris


-----Original Message-----
From: R-sig-mixed-models <r-sig-mixed-models-bounces using r-project.org> On Behalf Of Vaida, Florin via R-sig-mixed-models
Sent: Tuesday, February 20, 2024 9:25 PM
To: Ben Bolker <bbolker using gmail.com>; r-sig-mixed-models using r-project.org
Subject: Re: [R-sig-ME] lme varFix under ML fit does not match coefficients standard error



Waarschuwing: Deze e-mail is afkomstig van buiten de organisatie. Klik niet op links en open geen bijlagen, tenzij u de afzender herkent en weet dat de inhoud veilig is.
Caution: This email originated from outside of the organization. Do not click links or open attachments unless you recognize the sender and know the content is safe.



Thanks Ben, interesting!
Any particular reason for this default choice?
It sounds like separating parameter estimation (ML) from SE of the parameters (REML), but presumably when someone chooses ML for estimation they assume this goes to how the SE's are computed too.

Florin

On 2/18/24, 11:52 AM, "R-sig-mixed-models on behalf of Ben Bolker" <r-sig-mixed-models-bounces using r-project.org on behalf of bbolker using gmail.com> wrote:

       From ?summary.lme:

    adjustSigma: an optional logical value.  If ‘TRUE’ and the estimation
               method used to obtain ‘object’ was maximum likelihood, the
               residual standard error is multiplied by sqrt(nobs/(nobs -
               npar)), converting it to a REML-like estimate.  This argument
               is only used when a single fitted object is passed to the
               function.  Default is ‘TRUE’.




    On 2024-02-18 2:34 p.m., Vaida, Florin via R-sig-mixed-models wrote:
    > Hello all,
    >
    > This is probably known, but it’s news to me: the standard errors printed for the lme model fit run under method=”ML” are in fact those computed under method=”REML”.
    > Is this the expected behavior?  And if so, are there any reasons for this choice?
    > Reproducible example below.
    >
    > Thanks,
    > Florin
    >
    > library(nlme)
    > fit.reml =  lme(log(conc) ~ Time, random=~1|Subject, data=Glucose, na.action=na.omit, method="REML")
    > (se.reml = summary(fit.reml)$tTable[,2])
    > ## (Intercept)        Time
    > ## 0.019457141 0.005829144
    > (se.reml = sqrt(diag(summary(fit.reml)$varFix)))
    > ## (Intercept)        Time
    > ## 0.019457141 0.005829144
    > fit.ml =  lme(log(conc) ~ Time, random=~1|Subject, data=Glucose, na.action=na.omit, method="ML")
    > (se.ml = summary(fit.ml)$tTable[,2]) # they match the REML SE’s
    > ## (Intercept)        Time
    > ## 0.019457141 0.005829144
    > (se.ml = sqrt(diag(summary(fit.ml)$varFix))) # they do not match the tTable SE’s
    > ## (Intercept)        Time
    > ## 0.019405324 0.005813621
    >
    >   [[alternative HTML version deleted]]
    >
    > _______________________________________________
    > R-sig-mixed-models using r-project.org mailing list
    > https://urldefense.com/v3/__https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models__;!!LLK065n_VXAQ!maZY0riqf0o-YADm5wrg-eNFM4LHerS92ofHlC_Pv5pnBNSjJYPepDr4S5Y16jYfN9zHl3B6SPmJz71d$

    _______________________________________________
    R-sig-mixed-models using r-project.org mailing list
    https://urldefense.com/v3/__https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models__;!!LLK065n_VXAQ!maZY0riqf0o-YADm5wrg-eNFM4LHerS92ofHlC_Pv5pnBNSjJYPepDr4S5Y16jYfN9zHl3B6SPmJz71d$


_______________________________________________
R-sig-mixed-models using r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


More information about the R-sig-mixed-models mailing list