[R-sig-ME] no df of freedom to test the effect of an interaction in (lmer) mixed-model

Andrew Robinson A.Robinson at ms.unimelb.edu.au
Thu Mar 1 22:00:22 CET 2007


Hi Ahimsa,

I think that the problem is that model red1, as you have expressed it,
doesn't resepct heredity (in the context of linear models).  It
contains a three-way interaction term (S:C:t) but omits one of the
underlying two-way interaction terms (S:t).

Generally, bad things can happen when you do this, including a lack of
invariance to location and scale transformations.

In this case I'm guessing that lmer is just changing the
parameterization.  Here's an analogy.  (S:t) previously played the
role of one of the margins for (S:C:t), so lmer only reported one
level of (S:C:t), as the other level is aliased with the margin.  If
you eliminate the margin then lmer has to report two levels.  It's
like the difference between:

> x <- factor(c(1,1,2,2))
> y <- 1:4     
> lm(y ~ x)

Call:
lm(formula = y ~ x)

Coefficients:
(Intercept)           x2  
        1.5          2.0  

> lm(y ~ x - 1)

Call:
lm(formula = y ~ x - 1)

Coefficients:
 x1   x2  
1.5  3.5  

> 

Notice that even though we removed the intercept, R still used two
degrees of freedom to express the model.

I hope that this helps,

Andrew




On Thu, Mar 01, 2007 at 06:18:06PM +0900, ahimsa campos-arceiz wrote:
> Dear sig-mixed-models list-members,
> 
> I placed this question before on R-help main list, but I think this forum is
> more appropriate. I will reformulate the question though:
> 
> I have a model such as:
> full <- lmer(y ~ S*C*t + (S*t | id), method="ML")   # both S and C are
> factors with two levels each
> 
> First I want to test the effect of the interaction S:t. I fit the following
> reduced models:
> red1 <- lmer(y ~ S*C*t - S:t + (S*t | id), method="ML")   # excluding S:t
> from the fixed effects, and
> red2 <- lmer(y ~ S*C*t + (S + t | id), method="ML")   # excluding S:t:id
> from the random effects
> 
> then:
> anova(full, red1, red2)
> 
>             Df     AIC     BIC  logLik  Chisq Chi Df Pr(>Chisq)
> red2     14 10525.0 10614.5 -5248.5
> full       18 10493.7 10608.8 -5228.8 39.346  4  5.908e-08 ***
> red1     18 10493.7 10608.8 -5228.8  0.000   0
> 
> The problem is that there is no degree of freedom to see differences between
> full and red1
> 
> When looking at the fixed effects of both models:
> 
> 1. full model
> Fixed effects:
>                                   Estimate Std. Error t value
> (Intercept)                6.7822267  0.1475245   45.97
> S                           -0.4902659  0.0814779   -6.02
> C                           -0.2428243  0.2142228   -1.13
> t                             0.0122836  0.0015711    7.82
> S:C                         0.3716986   0.1232596    3.02
> S:t                         -0.0001985  0.0028627   -0.07
> C:t                         -0.0056198  0.0023020   -2.44
> S:C:t                      0.0022598  0.0040779    0.55
> 
> 2. red1 model
> Fixed effects:
>                               Estimate Std. Error t value
> (Intercept)              6.7822722  0.1483886   45.71
> S                         -0.4902454  0.0816564   -6.00
> C                         - 0.2428498  0.2154232   -1.13
> t                           0.0122829  0.0015684    7.83
> S:C                      0.3716582  0.1235135    3.01
> C:t                       -0.0056161  0.0022982   -2.44
> S:C1:t                  - 0.0001986  0.0028619   -0.07
> S:C2:t                   0.0020583  0.0029032    0.71
> 
> red1 model is analyzing the interaction S:C:t using both levels of C so that
> the result is one more factor, and therefore no degree of freedom left to
> compare with the full model.
> 
> Can anybody suggest how should I test the effect of S:t??
> 
> 
> Thank you very much for any feedback!!
> 
> Ahimsa
> 
> 
> For a more comprehensible explanation of the data you can check my previous
> post ( http://tolstoy.newcastle.edu.au/R/e2/help/07/03/11468.html). Below
> there is toi data and script.
> 
> **********************************************************************************************************
> 
> library(lme4)
> 
> # dataset
> A <- as.factor(rep(1:2, each=600))
> id <- as.factor(rep(1:6, each=200))
> S <- as.factor(rep(1:2, each=100, times=6))
> R <- as.factor(rep(1:10, each = 10, times = 12))
> t <- rep(c(2,4,6,8,10,12,14,16,18,20), times=120)
> y <- rnorm(1200, mean=100, sd=25)
> dummy.df <- data.frame(A,id,S,R,t,y)
> summary(dummy.df)
> str(dummy.df)
> 
> full.model <- lmer(y ~ A*S*t + (S*t|id), dummy.df, method="ML")
> summary(full.model)
> 
> red.model1 <- lmer(y ~ A*S*t - S:t + (S*t|id), dummy.df, method="ML")
> summary(red.model1)
> 
> red.model2 <- lmer(y ~ A*S*t + (S+t|id), dummy.df , method="ML")
> summary(red.model2)
> anova (full.model, red.model1, red.model2)
> 
> -- 
> ahimsa campos-arceiz
> www.camposarceiz.com
> 
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

-- 
Andrew Robinson  
Department of Mathematics and Statistics            Tel: +61-3-8344-9763
University of Melbourne, VIC 3010 Australia         Fax: +61-3-8344-4599
http://www.ms.unimelb.edu.au/~andrewpr
http://blogs.mbs.edu/fishing-in-the-bay/




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