[R-sig-eco] heteroscedasticity and d.f.

[Jasper]

Dear all,
I am analyzing an experiment where the treatments cause heteroscedasticty in
the residuals. Following Pinheiro & Bates (2000 Mixed effects models in S
and S-plus) and Zuur et al (2009 Mixed effects models and extensions in
ecology with R) I try to include seperate treatment variances using varIdent
structure. This improves thing according to my residual plots and a LR test.
However, what I don't understand is that when I test for the effect of my
treatment, I don't seem to loose any degrees of freedom for modelling these
seperate variances! Does anyone know why this is?

Below a simple example dataset with 4 observations in each of 5 treatments
(the sample size is insufficient to show the differences in variance so the
LR is not significant, but that's not my point). Comparing both models I see
a difference in d.f. among the models for the LR test of variance structures
(6 vs. 10), but when I look at the tests of fixed effects for both models,
both still have 15 residual d.f. eventhough the second model needs to
estimate 5 variances instead of one.

I am curious about your thoughts!

Best wishes,
Jasper Wubs

> library(nlme)
> m1 <- gls(SM~Treat,weights=NULL,data=data)
> vf1 <- varIdent(form=~1|Treat) #estimate a variance per treatment
> m1v <- gls(SM~Treat,weights=vf1,data=data)
> anova(m1,m1v)

    Model df       AIC       BIC   logLik   Test  L.Ratio p-value
m1      1  6 -85.76002 -81.51172 48.88001                        
m1v     2 10 -84.85906 -77.77855 52.42953 1 vs 2 7.099036  0.1307

> anova(m1)
Denom. DF: 15 
            numDF   F-value p-value
(Intercept)     1 10246.756  <.0001
Treat           4     1.937  0.1564

> anova(m1v)
Denom. DF: 15 
            numDF   F-value p-value
(Intercept)     1 27538.274  <.0001
Treat           4     7.733  0.0014



--
View this message in context: http://r-sig-ecology.471788.n2.nabble.com/heteroscedasticity-and-d-f-tp7578507.html
Sent from the r-sig-ecology mailing list archive at Nabble.com.