# [R] R: lmer specification for random effects: contradictory reults

Benedetta Cesqui b.cesqui at hsantalucia.it
Mon Nov 25 15:45:30 CET 2013

```Dear Thierry,

thank you for the quick reply.
I have only one question about the approach you proposed.
As you suggested, imagine that the model we end up after the model selection
procedure is:

mod2.1 <- lmer(dT_purs ~ T + Z + (1 +T+Z| subject), data =x, REML= FALSE)

According to the common procedures specified in many manuals and recent
papers, if I want to compute the p_values relative to each term, I will
perform a likelihood test, in which the deviance of the (-2LL) of a model
containing the specific term is compared to another model without it.
In the case of the fixed effect terms I have no problem in the
interpretation of the results. Each comparison returns a significance
associated with the estimated coefficient of the term.
Thus in this case:

mod2.2 <- lmer(dT_purs ~ Z + (1 +T+Z|soggetto)  , data = x, REML = FALSE)
mod2.3 <- lmer(dT_purs ~ T + (1 +T+Z|soggetto)  , data = x, REML = FALSE)
anova(mod2.1, mod2.2)
p_valueT = 3.203e-05
anova(mod2.1, mod2.3)
p_valueZ = 0.001793

What about the p_value relative to the (1+T+Z|subject)?
One option is to compute:
mod2.4 <- lm(dT_purs ~ T + (1 +T+Z|soggetto)  , data = x)
and then execute the loklikelihood test as follows:

L0 <-logLik(mod2.4)
L1 <-logLik(mod2.1)
LR <--2*(L1-L0)
pv <- pchisq(LR,2,ncp = 0, lower.tai=FALSE,log.p = FALSE)

However, what can I conclude on the random slopeif it is significant?

With the previouse approach using the model:
mod2 <- lmer(dT_purs ~ T + Z + (1 +T| subject) + (1+ Z| subject), data =x)

The comparison among the models in which the different termd were
included/excluded provided me the following results:
p_valueT = 1.269e-07;
p_valueZ =0.00322
p_valueTS =  0.4277
p_valueZS = 0.005701

I interpreted the ones relative to the random effects as if the subjects
differed not only in their overall responses, but also in the nature of
their response dT_purse values in the different T conditions, but not in the
different Z conditions.

Benedetta

-----Messaggio originale-----
Da: ONKELINX, Thierry [mailto:Thierry.ONKELINX at inbo.be]
Inviato: lunedì 25 novembre 2013 14:48
A: Benedetta Cesqui; r-help at r-project.org
Cc: r-sig-mixed-models at r-project.org
Oggetto: RE: [R] lmer specification for random effects: contradictory reults

Dear Benedetta,

I think you might want (1+T+Z|subject) as random effects  rather than
(1+T|subject) + (1 + Z|subject). The latter has two random intercepts per
subject: a recipe for disaster.

Follow-up posts should only go to the mixed models mailing list which I'm
cc'ing.

Best regards,

ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
Forest team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium
+ 32 2 525 02 51
+ 32 54 43 61 85
Thierry.Onkelinx at inbo.be
www.inbo.be

To call in the statistician after the experiment is done may be no more than
asking him to perform a post-mortem examination: he may be able to say what
the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey

-----Oorspronkelijk bericht-----
Van: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
Namens Benedetta Cesqui
Verzonden: maandag 25 november 2013 11:13
Aan: r-help at r-project.org
Onderwerp: [R] lmer specification for random effects: contradictory reults

Hi All,

I was wondering if someone could help me to solve this issue with lmer.
In order to understand the best mixed effects model to fit my data, I
compared the following options according to the procedures specified in many
papers (i.e. Baayen
&url=http%3A%2F%2Fwww.ualberta.ca%2F~baayen%2Fpublications%2FbaayenDavidsonB
ates.pdf&ei=FhqTUoXuJKKV7Abds4GYBA&usg=AFQjCNFst7GT7mBX7w9lXItJTtELJSKWJg&si
g2=KGA5MHxOvEGwDxf-Gcqi6g&bvm>  R.H. et al 2008) Here, dT_purs is the
response variable, T and Z are the fixed effects, and subject is the random
effect. Random and fixed effects are crossed.:

mod0 <- lmer(dT_purs ~ T + Z + (1|subject), data = x)
mod1 <- lmer(dT_purs ~ T + Z + (1 +tempo| subject), data = x)
mod2 <- lmer(dT_purs ~ T + Z + (1 +tempo| subject) + (1+ Z| subject), data =
x)
mod3 <- lmer(dT_purs ~ T * Z + (1 +tempo| subject) + (1+ Z| subject), data =
x)
mod4 <- lmer(dT_purs ~ T * Z + (1| subject), data = x)

anova(mod0, mod1,mod2, mod3, mod4)

Data: x

Models:

mod0: dT_purs ~ T + Z + (1 | subject)

mod4: dT_purs ~ T * Z + (1 | subject )

mod1: dT_purs ~ T + Z + (1 + T| subject)

mod2: dT_purs ~ T + Z + (1 + T| subject ) + (1 + Z | subject)

mod3: dT_purs ~ T * Z + (1 + T| subject) + (1 + Z | subject)

Df     AIC     BIC logLik deviance   Chisq Chi Df Pr(>Chisq)

mod0  5 -689.81 -669.46 349.91  -699.81

mod4  6 -689.57 -665.14 350.78  -701.57  1.7532      1   0.185473

mod1  7 -689.12 -660.62 351.56  -703.12  1.5504      1   0.213070

mod2 10 -695.67 -654.97 357.84  -715.67 12.5563      3   0.005701 **

mod3 11 -695.83 -651.05 358.92  -717.83  2.1580      1   0.141825

---

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

It turns out that mod2 has the right level of complexity for this dataset.

However when I looked at its summary, I got a correlation of -0.87 for the
random effects relative to the T effect and -1 for the random effects
relatively to the Z.

summary(mod2)

Linear mixed model fit by maximum likelihood ['lmerMod']

Formula: dT_purs ~T + Z + (1 + T | subject) + (1 + Z | subject)

Data: x

AIC       BIC    logLik  deviance

-695.6729 -654.9655  357.8364 -715.6729

Random effects:

Groups     Name        Variance  Std.Dev. Corr

subject   (Intercept) 0.0032063 0.05662

T       0.0117204 0.10826  -0.87

subject.1 (Intercept) 0.0005673 0.02382

Z           0.0025859 0.05085  1.00

Residual               0.0104551 0.10225

Number of obs: 433, groups: soggetto, 7

Fixed effects:

Estimate Std. Error t value

(Intercept)  0.02489    0.03833   0.650

T        0.52010    0.05905   8.808

Z           -0.09019    0.02199  -4.101

Correlation of Fixed Effects:

(Intr) tempo

T -0.901

Z      0.218 -0.026

If I understand correctly what the correlation parameters reported in the
table are, the correlation of 1 means that, for the Z effects the random
intercept is perfectly collinear with the random slope. Thus, we fit the
wrong model. A random intercept only model would have sufficed.

Am I correct?

If so, should I take mod1 (mod1 <- dT_purs ~ T + Z + (1 + T | subject )
instead of mod2 to fit my data?

Finally is a correlation value of -0.87 a too high or an acceptable value ?

Thanks for help me in advance!

Best

Benedetta

---

Benedetta Cesqui, Ph.D.

Laboratory of Neuromotor Physiology

IRCCS Fondazione Santa Lucia

Via Ardeatina 306

00179 Rome, Italy

tel: (+39) 06-51501485

fax:(+39) 06-51501482

E_mail:  b.cesqui at hsantalucia.it

[[alternative HTML version deleted]]

______________________________________________
R-help at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help