[R-sig-ME] Help with anova result for model comparison

Ben Bolker bbolker at gmail.com
Fri Jul 8 20:41:36 CEST 2011


-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

On 07/08/2011 12:24 PM, Catarina Miranda wrote:
> Hello;
> 
> When I run an anova to compare two lmer models I have in some situations
> (always when I compare a model with interactions with a model where I
> removed one simple fixed effect) an output that I don't fully understand. In
> the anova output the models have the same degrees of freedom, as well as AIC
> and BIC; the Chisq and its degrees of freedom are 0, and the p value is 1.
> This makes me think that I cannot make this model comparison, but I would
> like to understand why and to know how could I correctly compare two models
> of this type.
> The R output is below, any help or suggestions would be greatly appreciated!
> 


  Removing the main effect of "Origin" is just reparameterizing the
model, not actually changing it.

  The easiest way to figure this stuff out is to play with design
matrices for simple examples:

> d <- expand.grid(a=factor(1:2),b=factor(1:2))
> model.matrix(~(a+b)^2,data=d)
  (Intercept) a2 b2 a2:b2
1           1  0  0     0
2           1  1  0     0
3           1  0  1     0
4           1  1  1     1

> model.matrix(~(a+b)^2-a,data=d)
  (Intercept) b2 a2:b1 a2:b2
1           1  0     0     0
2           1  0     1     0
3           1  1     0     0
4           1  1     0     1

  (some details omitted)

  It is **occasionally** sensible (but often not: see Venables'
_Exegeses on Linear Models_, unpublished but Googlable) to really drop
the main effect term.  This is actually quite difficult to do in R
because model.matrix() (which is used internally by most if not all
modeling functions) isn't quite flexible enough. Some modeling functions
will take a raw design matrix; the only way I know to do it in lme4
would be to generate the design matrix yourself, drop the undesired
columns, and use the result as the data for a new fit ...

  Ben



> Catarina
> 
> R version 2.13.0 (2011-04-13)
> 
>>
> m2_Exp<-lmer(lat_perch~(Sex+Origin+trial)^2+SMI_Neo+(1|Nest_ID)+(1|ID),data=full,REML=F)
> 
> # lat_perch and SMI are numeric, but the other variables are factors
> 
>>
> m2_Exp_Origin<-lmer(lat_perch~(Sex+Origin+trial)^2+SMI_Neo-Origin+(1|Nest_ID)+(1|ID),data=full,REML=F)
> 
>> anova(m2_Exp,m2_Exp_Origin)  # NS but sething wierd is happening
> Data: full
> Models:
> m2_Exp: lat_perch ~ (Sex + Origin + trial)^2 + SMI_Neo + (1 | Nest_ID) +
> m2_Exp:     (1 | ID)
> m2_Exp_Origin: lat_perch ~ (Sex + Origin + trial)^2 + SMI_Neo - Origin + (1
> |
> m2_Exp_Origin:     Nest_ID) + (1 | ID)
>               Df  AIC    BIC logLik Chisq Chi Df Pr(>Chisq)
> m2_Exp        14 3056 3098.9  -1514
> m2_Exp_Origin 14 3056 3098.9  -1514     0      0          1
> 
>> summary(m2_Exp)
> Linear mixed model fit by maximum likelihood
> Formula: lat_perch ~ (Sex + Origin + trial)^2 + SMI_Neo + (1 | Nest_ID)
> +      (1 | ID)
>    Data: full
>   AIC  BIC logLik deviance REMLdev
>  3056 3099  -1514     3028    2865
> Random effects:
>  Groups   Name        Variance Std.Dev.
>  ID       (Intercept) 2516808  1586.45
>  Nest_ID  (Intercept)  481927   694.21
>  Residual             9961255  3156.15
> Number of obs: 158, groups: ID, 53; Nest_ID, 18
> 
> Fixed effects:
>                     Estimate Std. Error t value
> (Intercept)         -3995.97    3731.88  -1.071
> SexMale             -1754.35    1225.74  -1.431
> OriginUrban          1618.35    1268.25   1.276
> trialT2             -2989.72    1128.92  -2.648
> trialT3             -1457.01    1145.84  -1.272
> SMI_Neo                98.49      43.34   2.272
> SexMale:OriginUrban    64.10    1390.64   0.046
> SexMale:trialT2      1364.19    1230.84   1.108
> SexMale:trialT3      2570.80    1237.91   2.077
> OriginUrban:trialT2  -336.72    1231.60  -0.273
> OriginUrban:trialT3   821.13    1238.67   0.663
> 
> Correlation of Fixed Effects:
>             (Intr) SexMal OrgnUr trilT2 trilT3 SMI_Ne SxM:OU SxM:T2 SxM:T3
> OrU:T2
> SexMale
> -0.066
> OriginUrban -0.127
> 0.380
> trialT2     -0.106  0.312
> 0.299
> trialT3      0.009  0.324  0.304
> 0.499
> SMI_Neo     -0.962 -0.120 -0.063 -0.047
> -0.165
> SxMl:OrgnUr -0.047 -0.593 -0.578 -0.007 -0.030
> 0.159
> SexMl:trlT2  0.118 -0.498 -0.027 -0.609 -0.296 -0.027
> -0.004
> SexMl:trlT3  0.077 -0.505 -0.033 -0.304 -0.607  0.016  0.009
> 0.496
> OrgnUrbn:T2  0.135 -0.025 -0.482 -0.608 -0.293 -0.044 -0.007  0.061
> 0.029
> OrgnUrbn:T3  0.130 -0.029 -0.483 -0.301 -0.598 -0.039  0.001  0.031  0.070
> 0.498
> 
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.10 (GNU/Linux)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org/

iEYEARECAAYFAk4XT2AACgkQc5UpGjwzenMCygCcC8xoejtQEqnInkr7DowLbR+f
PzEAn0wgctfVY+vauvMfkBaFuRxTPaof
=pazH
-----END PGP SIGNATURE-----




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