[R-sig-ME] heterogeneous variance models

Hufthammer, Karl Ove karl.ove.hufthammer at helse-bergen.no
Mon May 11 17:03:35 CEST 2015


Also, the results don't seem to agree with the results from the old nlme package.

Example (using a slightly more sensible model than from the help file): Here's what the data looks like:

  xyplot(distance~age|Sex, groups=Subject, data=Orthodont, type="l", col="darkblue")

Let's fit lme4 and nlme models:

  l_new1 = lmer(distance ~ age + Sex + (age|Subject), data = Orthodont)
  l_new2 = update(l_new1, . ~ . + (0+dummy(Sex, "Female")|Subject))
  summary(l_new1)
  summary(l_new2)

  l_old1 = lme(distance ~ age + Sex, random=~age|Subject, data = Orthodont)
  l_old2 = update(l_old1, weights=varIdent(form=~1|Sex))
  summary(l_old1)
  summary(l_old2)

While the fixed effects estimates for the two models assuming homoscedasticity are very ~identical, the estimates from the heteroscedasticity models differ:

> fixef(l_new1)
(Intercept)         age   SexFemale 
 17.6351989   0.6601852  -2.1454883 
> fixef(l_new2)
(Intercept)         age   SexFemale 
 17.6487393   0.6601852  -2.1787238
> fixef(l_old1)
(Intercept)         age   SexFemale 
 17.6352002   0.6601852  -2.1454915
> fixef(l_old2)
(Intercept)         age   SexFemale 
 18.4385143   0.5795744  -1.9407706

It's a clever trick, but doesn't seem to work properly (or is it nlme that doesn't work properly?). (Also, it assumes that one *knows* which group has the lowest residual variance.)


Regards,
Karl Ove Hufthammer


-----Opphavleg melding-----
From: R-sig-mixed-models [mailto:r-sig-mixed-models-bounces at r-project.org] På vegner av Aaron Mackey
Sent: 11. mai 2015 15:50
To: Ben Bolker
Kopi: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] heterogeneous variance models

Thanks for this. What's confusing to me in the example is that the residual std. dev is still 1.3100 whether the dummy term is included or not; I would have thought that the residual sd without the dummy term would be intermediate between the two sex-specific sd's seen in the full model (reported as 1.0521 and 1.311).  The reported residuals also don't seem to change between the two models.

thanks for any clarification you can provide, -Aaron

On Sun, May 10, 2015 at 4:22 PM, Ben Bolker <bbolker at gmail.com> wrote:

> -----BEGIN PGP SIGNED MESSAGE-----
> Hash: SHA1
>
>
>    Someone sent me an e-mail recently asking how/whether heterogeneous 
> variances can be implemented in lme4.  I can't find it any more, so 
> I'm answering here in hopes that they'll see it. The short answer is 
> that they can, at least for heterogeneity among categorical groups 
> (similar to what varIdent does in the nlme package), although it's a 
> bit of a nuisance for categorical variables with lots of levels, as 
> the variance terms all need to be specified separately.
>
>    The trick is to set up a dummy variable for each level: there is a 
> helper function dummy() for this purpose, and example("dummy") gives 
> an example showing how to estimate heterogeneous residual variances 
> for females vs. males in the standard sleep study example.
>
>   Ben Bolker
>
>
> -----BEGIN PGP SIGNATURE-----
> Version: GnuPG v1.4.11 (GNU/Linux)
>
> iQEcBAEBAgAGBQJVT73oAAoJEOCV5YRblxUHlzQIAMuoYv9K06NGG+3mEcvHcmPW
> Xp7rLyxRfiTtRkmr/UieR52WpFDapxNgMurLT1aESD5TWeSBTXmeev2ltajPM5It
> A5Mz1cNob1Tl67zsBCUrKMVVclRicO/Ymbpab1QJYcXxGAVZ03gkaguPWJpw3aHQ
> kt2GZW6SGidbI17PwZXTyVb69cJ8p/mq459LwMAFNusD4MCiIt1/K2HTL8NnekVn
> G7XPMDBg/JoyHmnqWpJbhzccjxZetNMtK6GyXaamlnyPJaJI0i9j6t/X0jKiAqO0
> k1BdwNkLo03sneHodyhdRqAS5KMzeNSo8FOinDIUCLtVikz5k0IbQhHO6JWChz4=
> =KPfA
> -----END PGP SIGNATURE-----
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list 
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>

	[[alternative HTML version deleted]]

_______________________________________________
R-sig-mixed-models at 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