[R] Mixed-effects models: question about the syntax to introduce interactions

From: Anaid Diaz <syadp_at_yahoo.com.mx>
Date: Tue, 25 Mar 2008 08:42:34 -0500 (CDT)


hello everyone,

I would like to as for advice for the use of “lmer”
(package ‘lme4’) and writing the proper syntax to best
describe my data using a mixed-effects model.

I have just started to use these models, and although I have read some good examples (Extending the Linear Model with R, Faraway 2005; and the R book, Crawley 2007), I am still not sure of the syntax to test my hypothesis.

Thanks in advance for reading me.

Briefly, I describe the data and the situation:

I want to describe the age-specific fecundity of the ith individual from the jth replicate (or line) from the kth strain.

Variables:

Categorical factors:

A[a] = Age (1,2,3…n=8) #Because the fecundity is not linear, I decided to include it in the model as a factor
s[k] =strain (A and B, n=2) # for the moment two, but it’s likely to increase as the work progresses l[j] = line (1,2,..n=10)
i[i] = Ind(1,2… n=50)

(Note: I use capital letters for fixed factors and low
case for random effects)

Because the experimental design, the data follows a hierarchical structure: where the ith individual is nested within the jth line, and line within the kth strain

Continuous (response) variable:
Y =Age-specific fecundity (362 observations)

Models:

Because I was (I am still) not sure of how to include all the variables in a single model, I started by splitting up the data and assessing which is the best model for each strain, therefore the “Simplest” model for each strain is:

Linear model.

Y[aij] = A[a] + error[aij]

R code:

m1 <- lm(fecudnity ~ Age)

Reduced mixed-effect model:

Y[aij] = A[a] + l[j] + i[i] + error[aij]

R code:

m2 <- lmer(fecudnity ~ Age + (1 | line/ind), method=”ML”)

And a “Full mixed-effects model” model (looking for interactions between Age and line/ind)

Y[aij] = A[a] + A[a]*l[j] + A[a]*i[i] + error[aij]

R code:

m3 <- lmer(fecudnity ~ Age + (Age | line/ind), method=”ML”)

I have used Likelihood test ratio (LTR) to compare between models, and I have found that for strain A the best model is m3 (X^2 [36 d.f] =164.8, p-value= 4.73e-13), whereas for strain B, the best one is m1
(X^2 [2 d.f] =1.47, p-value= 0.473). Therefore, I
interpret these results as follow:

These results, on its own, are interesting, but I would like to have a model where I include both strains (and still can make some interpretations)

My first guess is

m4 <- lmer(fecudnity ~ Age + (1 | strain/line/ind), method=”ML”)
m5 <- lmer(fecudnity ~ Age + (Age | strain/line/ind), method=”ML”)

Using LTR, I find that m5 describes better the data
(X^2 [105 d.f] = 347.15, p-value < 2.2 e-16), but I
feel like I can not say much of which strain has more individual variance (or perhaps I am wrong and not looking in the right place).

Then I though about using strain as a fixed factor, because I am now interested in the differences between strains

m6 <- lmer(fecudnity ~ Age * Strain + (Age |strain/line/ind), method=”ML”)

or perhaps include it in the random interaction?

m7 <- lmer(fecudnity ~ Age + (Age * Strain |strain/line/ind), method=”ML”)

I have to be honest, at this point, I am just not sure of how to write the model to describe the age-specific fecundity and test the hypothesis of whether one strain shows more variance between individuals and lines or not. I hope some one could give some advise.

Thanks in advance

Anaid Diaz



ˇCapacidad ilimitada de almacenamiento en tu correo! No te preocupes más por el espacio de tu cuenta con Correo Yahoo!:

R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Tue 25 Mar 2008 - 14:04:19 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Tue 25 Mar 2008 - 14:30:23 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive