[R-sig-ME] problem with nested mixed model

[Seth Bigelow]

Well, if it were me, I would begin with the simplest reasonable model:

M1 <- lme(altura~ vegetation*bromeliad, random=~1|site)

(You don't have to put vegetation + bromeliad + vegetation*bromeliad, just
vegetation*bromeliad will do)

Then also do the more complicated model:

M2 <- lme(altura ~ vegetation*bromeliad, random=~1|vegetation/site)...

(I think this is better syntax than the way you had it before)

...and do a likelihood ratio test to see if there is an improvement in fit
by nesting site within vegetation.

anova(M1,M2)

You should be checking residuals for these models -- Chapter 1 in Pinheiro &
Bates has all this. The quantiles of your residuals look rather unbalanced,
some kind of transformation may be necessary

Boa sorte

-Seth



-----Original Message-----
From: r-sig-mixed-models-bounces at r-project.org
[mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf Of Érika Tsuda
Sent: Friday, April 19, 2013 11:53 AM
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] problem with nested mixed model





Hi guys, 


I am trying to fit a nested mixed model using the lme function but not sure
on how to state a random factor is nested within another factor. 

Here is my sampling design: 

I measured height of seedlings (my response variable called height) in two
different vegetation strata (herbaceous and shrubby), this is my factor
vegetation, which I believe is fixed. 


At each strata I collected in three different site, this is my factor site,
which I believe is random and nested whithin vegetation. 


Also, at each vegetation strata I collected seedlings inside and outside
bromeliads, this is my factor bromeliad, which I believe is fixed and
orthogonal to vegetation. 


So, I tried to model this using the function lme as follows: 


model<-
lme(height~vegetation+bromeliad+vegetation*bromeliad+(vegetation/site),
~1|site, data=altura)


Resulting in: 


Linear mixed-effects model fit by REML Data: altura  AIC      BIC    logLik
13290.98 13342.81 -6635.491 Random effects: Formula: ~1 | site (Intercept)
Residual
StdDev:    4.682701  36.8902 Fixed effects: height ~ vegetation + bromeliad
+ vegetation * bromeliad + (vegetation/site)  Value Std.Error   DF   t-value
p-value
(Intercept)                  57.37983  5.856448 1314  9.797718  0.0000
vegetationHERB              -21.48592  4.363718 1314 -4.923764  0.0000
bromeliadBRO                  6.73168  2.441297 1314  2.757421  0.0059
vegetationHERB:bromeliadBRO  12.72138  5.274164 1314  2.412019  0.0160
vegetationARB:sitea2          5.41959  7.776715 1314  0.696900  0.4860
vegetationHERB:sitea2         0.05141  9.265842 1314  0.005549  0.9956
vegetationARB:sitea3          1.69237  7.570281 1314  0.223555  0.8231
vegetationHERB:sitea3       -14.65097  7.692126 1314 -1.904672  0.0570
Correlation:  (Intr) vgHERB brmBRO vHERB:B vARB:2 vHERB:2 vARB:3
vegetationHERB              -0.484

bromeliadBRO                -0.202  0.271

vegetationHERB:bromeliadBRO  0.093 -0.282 -0.463

vegetationARB:sitea2        -0.722  0.323 -0.002  0.001

vegetationHERB:sitea2       -0.404 -0.144  0.000 -0.032   0.609

vegetationARB:sitea3        -0.742  0.332 -0.001  0.000   0.559  0.313

vegetationHERB:sitea3       -0.487 -0.176  0.000 -0.024   0.367  0.394
0.753 Standardized Within-Group Residuals: Min         Q1        Med
Q3        Max 
-1.7492750 -0.5854722 -0.2114832  0.3335703  7.7013967  Number of
Observations: 1324 Number of Groups: 3 

When I ask the anova table, I get: 

numDF denDF  F-value p-value
(Intercept)              1  1314 350.5013  <.0001
vegetation               1  1314 155.3877  <.0001
bromeliad                1  1314  18.7275  <.0001
vegetation:bromeliad     1  1314   5.3803  0.0205
vegetation:site          4  1314   2.6398  0.0325

Based on this results I guess the analysis is right, however I am not sure
if I assigned the random and fixed factors properly. 
Also, I am not sure if the factor sites is nested whithin vegetation.
Another thing I dont know is what does the number 1 in ~1|sitemeans, since
when I changed for 0 the results is the same, but the models fails if I
remove this number.

Sorry by the big message, and thank you very much for your attention. 



----Irika Tiemi Tsuda
Bisloga, MSc. em Ecologia
Laboratsrio de Ecologia Vegetal

Universidade Federal de Santa Catarina
(48) 9165-8822/ (48) 3721-5520
	[[alternative HTML version deleted]]