Equivalence of (0 + factor|group) and (1|group) + (1|group:factor) random effect specifications in case of compound symmetry

Sat, Sep 23, 2017 10:26 AM

Thanks, this really helped me understand the models!

There is one thing I still don't understand: I think that compound symmetry
implies equal variances within the factor levels and equal
covariances/correlations between the levels. If there are equal variances
within the factor levels and no correlations between the levels I
understand the equivalence of m1, m.zcp and m2.

But what if there are correlations between the levels? you explained that
m2 forces the correlation parameters to zero. But Douglas Bates states that
the models are equivalent given compound symmetry (which I think doesn't
imply zero correlations).

Just to double check: This m1 -> m.zcp -> m2 -> m.vio is a reasonable model
reduction approach, right?

2017-09-23 13:51 GMT+02:00 Reinhold Kliegl <reinhold.kliegl at gmail.com>:

This should always be "k(k-1)/2" correlation parameters, of course.

Also perhaps you may want to read this post to the list by Douglas Bates:
https://stat.ethz.ch/pipermail/r-sig-mixed-models/2009q1/001736.html

On Sat, Sep 23, 2017 at 1:41 PM, Reinhold Kliegl <
reinhold.kliegl at gmail.com> wrote:

The models are on a continuum of complexity wrt to the random-effects
structure. Specifically:

m1 estimates k variance components and k(k-1) correlation parameters for
the k levels of factor f

m2 estimates 1 variance component for the intercept and a 1 variance
component for the k-1 contrasts defined for the k levels of factor f, that
it constrains the k-1 contrasts for factor f to the same value. The
correlation parameters are forced to zero.

The continuum becomes transparent if one re-parameterizes m1 as a m1a
(see below) with 1 variance component for the intercept and k-1 variance
components for the k levels of factor f, and k(k-1) correlation parameters.
 m1 and m1a have the same number of parameters and identical deviance.

m1a <- lmer(y ~ factor + (factor|group))

There is an additional model specification on this continuum between m1
and m2.  If contrasts are converted to numeric covariates, one can force
correlation parameters to zero, but estimate different variance components
for the k-1 contrasts of factor f. We call this the zero-correlation
parameter model.

cB_A <- model.matrix(m1)[,2]
cC_A <- model.matrix(m1)[,3]

m.zcp <- lmer(score ~ 1 + Machine + (1 + cB_A + cC_A  ||  Worker),
 data=Machines, REML=FALSE)
print(summary(m.zcp), corr=FALSE)

Note that m.zcp without the double-bar is equivalent to m1 and m1b.
m.1b <- lmer(score ~ 1 + Machine + (1 + cB_A + cC_A  |  Worker),
 data=Machines, REML=FALSE)
print(summary(m.zcp), corr=FALSE)

Of course, there is also a simpler model than m2 - the
varying-intercept-only model:

m.vio <- lmer(score ~ 1 + Machine + (1 | Worker),  data=Machines,
REML=FALSE)
print(summary(m.zcp), corr=FALSE)

Here is some R code demonstrating all this for the Machines data.

####
library(lme4)
#library(RePsychLing)

data("Machines", package = "MEMSS")

# OP m1
m1 <- lmer(score ~ 1 + Machine + (0 + Machine | Worker), data=Machines,
REML=FALSE)
print(summary(m1), corr=FALSE)

# re-parameterization of m1
m1a <- lmer(score ~ 1 + Machine + (1 + Machine | Worker), data=Machines,
REML=FALSE)
print(summary(m1a), corr=FALSE)

# alternative specification of m1a
cB_A <- model.matrix(m1)[,2]
cC_A <- model.matrix(m1)[,3]

m1b <- lmer(score ~ 1 + Machine + (1 + cB_A + cC_A  | Worker),
 data=Machines, REML=FALSE)
print(summary(m1b), corr=FALSE)

anova(m1, m1a, m1b)

# zero-correlation parameter LMM
m.zcp <- lmer(score ~ 1 + Machine + (1 + cB_A + cC_A  || Worker),
 data=Machines, REML=FALSE)
print(summary(m.zcp), corr=FALSE)

# OP m2
m2 <- lmer(score ~ 1 + Machine + (1 | Worker) + (1 | Machine:Worker),
data=Machines, REML=FALSE)
print(summary(m2), corr=FALSE)

# varying-intercept-only LMM
m.vio <- lmer(score ~ 1 + Machine + (1 | Worker), data=Machines,
REML=FALSE)
print(summary(m.vio), corr=FALSE)

anova(m1, m.zcp, m2, m.vio)

sessionInfo()
####

You may also want to look this RPub:
http://www.rpubs.com/Reinhold/22193

On Sat, Sep 23, 2017 at 11:43 AM, Maarten Jung <
Maarten.Jung at mailbox.tu-dresden.de> wrote:

Hello everyone,

I have a question regarding the equivalence of the following models:

m1 <- lmer(y ~ factor + (0 + factor|group))
m2 <- lmer(y ~ factor + (1|group) + (1|group:factor))

Douglas Bates states (slide 91 in this presentation [1])  that these
models
are equivalent in case of compound symmetry.

1. I realized that I don't really understand the random slope by factor
model (m1) and espacially why it reduces to m2 given compound symmetry.
Also, why is there no random intercept in m1?
Can anyone explain the difference between the models and how m1 reduces
to
m2 in an intuitive way.

2. If m1 is a special case of m2 ? this could be an interesting option
for
model reduction but I?ve never seen something like m2 in papers. Instead
they suggest dropping the random slope and thus the interaction
completely
(e.g. Matuschek et al. 2017 [3]).
What do you think about it?

Please note that I asked the question on Stack Exchange [2] but some
consider it as off-topic. So, I hope one of you can help me out.


Best regards,
Maarten

[1] http://www.stat.wisc.edu/~bates/UseR2008/WorkshopD.pdf
[2] https://stats.stackexchange.com/q/304374/136579
[3] https://doi.org/10.1016/j.jml.2017.01.001

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

Equivalence of (0 + factor|group) and (1|group) + (1|group:factor) random effect specifications in case of compound symmetry

Thread (8 messages)