Equivalência das especificações de efeito aleatório (0 + fator | grupo) e (1 | grupo) + (1 | grupo: fator) em caso de simetria composta

13

Douglas Bates afirma que os seguintes modelos são equivalentes "se a matriz de variância-covariância para efeitos aleatórios com valor vetorial tiver uma forma especial, chamada simetria composta" ( slide 91 desta apresentação ):

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

Especificamente, Bates usa este exemplo:

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

m1a <- lmer(score ~ Machine + (0 + Machine|Worker), Machines)
m2a <- lmer(score ~ Machine + (1|Worker) + (1|Worker:Machine), Machines)

com as saídas correspondentes:

print(m1a, corr = FALSE)

Linear mixed model fit by REML ['lmerMod']
Formula: score ~ Machine + (0 + Machine | Worker)
   Data: Machines
REML criterion at convergence: 208.3112
Random effects:
 Groups   Name     Std.Dev. Corr     
 Worker   MachineA 4.0793            
          MachineB 8.6253   0.80     
          MachineC 4.3895   0.62 0.77
 Residual          0.9616            
Number of obs: 54, groups:  Worker, 6
Fixed Effects:
(Intercept)     MachineB     MachineC  
     52.356        7.967       13.917  

print(m2a, corr = FALSE)

Linear mixed model fit by REML ['lmerMod']
Formula: score ~ Machine + (1 | Worker) + (1 | Worker:Machine)
   Data: Machines
REML criterion at convergence: 215.6876
Random effects:
 Groups         Name        Std.Dev.
 Worker:Machine (Intercept) 3.7295  
 Worker         (Intercept) 4.7811  
 Residual                   0.9616  
Number of obs: 54, groups:  Worker:Machine, 18; Worker, 6
Fixed Effects:
(Intercept)     MachineB     MachineC  
     52.356        7.967       13.917

Alguém pode explicar a diferença entre os modelos e como se m1reduz a m2(dada simetria composta) de maneira intuitiva?

statmerkur
fonte
6
+1 e, imho, isso é absolutamente sobre tópico. Vote para reabrir.
Ameba diz Reinstate Monica
2
@ Peter Flom, por que você considera essa pergunta fora de tópico?
Statmerkur
3
Talvez não estivesse claro que você estava perguntando sobre os modelos, e não sobre a lme4sintaxe. Seria útil - e ampliaria o grupo de possíveis respondedores - se você os explicasse para pessoas não familiarizadas lme4.
Scortchi - Restabelece Monica
Parece que é sobre codificação.
Peter Flom - Restabelece Monica
1
Se for útil, aqui estão duas boas postagens sobre o que a sintaxe lme4 está fazendo e qual simetria composta está no contexto de modelos mistos (consulte as respostas aceitas nas duas perguntas). stats.stackexchange.com/questions/13166/rs-lmer-cheat-sheet e stats.stackexchange.com/questions/15102/…
Jacob Socolar 28/17

Respostas:

11

Neste exemplo, há três observações para cada combinação das três máquinas (A, B, C) e os seis trabalhadores. Usarei para denotar uma matriz de identidade n- dimensional e 1 n para denotar um vetor n- dimensional de uns. Digamos que y é o vetor de observações, que assumirei que é ordenado pelo trabalhador e depois pela máquina e depois pela replicação. Seja μ os valores esperados correspondentes (por exemplo, os efeitos fixos) e γ seja um vetor de desvios específicos do grupo aos valores esperados (por exemplo, os efeitos aleatórios). Condicional em γ , o modelo para y pode ser escrito:Inn1nnyμγγy

yN(μ+γ,σy2I54)

onde é a variação "residual".σy2

Para entender como a estrutura de covariância dos efeitos aleatórios induz uma estrutura de covariância entre as observações, é mais intuitivo trabalhar com a representação "marginal" equivalente , que se integra aos efeitos aleatórios . A forma marginal deste modelo é,γ

yN(μ,σy2I54+Σ)

Aqui, é uma matriz de covariância que depende da estrutura de γ (por exemplo, os "componentes de variância" subjacentes aos efeitos aleatórios). Vou me referir a Σ como a covariância "marginal".ΣγΣ

No seu caso m1, os efeitos aleatórios se decompõem como:

γ=Zθ

Z=I1813θT=[θ1,A,θ1,B,θ1,Cθ6,A,θ6,B,θ6,C]

θN(0,I6Λ)

ΛΛσθτ

Λ=[σθ2+τ2τ2τ2τ2σθ2+τ2τ2τ2τ2σθ2+τ2]

Λ

Σ=Z(I6Λ)ZTσθ2+τ2+σy2i,ju,v

cov(yi,u,yj,v)={0if ijτ2if i=j,uvσθ2+τ2if i=j,u=v

Para você m2, os efeitos aleatórios se decompõem em:

γ=Zω+Xη

X=I619ωT=[ω1,A,ω1,B,ω1,C,,ω6,A,ω6,B,ω6,C]ηT=[η1,,η6]

ηN(0,ση2I6)
ωN(0,σω2I18)
ση2,σω2

m2Σ=σω2ZZT+ση2XXTσω2+ση2+σy2i,ju,v

cov(yi,u,yj,v)={0if ijση2if i=j,uvσω2+ση2if i=j,u=v

So ... σθ2σω2 and τ2ση2. If m1 assumed compound symmetry (which it doesn't with your call to lmer, because the random effects covariance is unstructured).

Brevity is not my strong point: this is all just a long, convoluted way of saying that each model has two variance parameters for the random effects, and are just two different ways of writing of the same "marginal" model.

In code ...

sigma_theta <- 1.8
tau         <- 0.5
sigma_eta   <- tau
sigma_omega <- sigma_theta
Z <- kronecker(diag(18), rep(1,3))
rownames(Z) <- paste(paste0("worker", rep(1:6, each=9)), 
                     rep(paste0("machine", rep(1:3, each=3)),6))
X <- kronecker(diag(6), rep(1,9))
rownames(X) <- rownames(Z)
Lambda <- diag(3)*sigma_theta^2 + tau^2

# marginal covariance for m1:
Z%*%kronecker(diag(6), Lambda)%*%t(Z)
# for m2:
X%*%t(X)*sigma_eta^2 + Z%*%t(Z)*sigma_omega^2
Nate Pope
fonte
1
Very nice answer! But I think the phrase "machine nested within worker" could be misleading as the same three machines appear in more than one (in fact every) level of worker.
statmerkur
@statmerkur Thanks, I've tried to clarify this line. Let me know if you have another suggestion.
Nate Pope
1
Should X be defined as X=I619?
S. Catterall Reinstate Monica
1
@S.Catterall Yup, that's a typo -- thanks for catching it! I've fixed in my answer.
Nate Pope
2
@statmerkur can you clarify what you mean? There's no continuous covariates here, so not sure what you mean by "slope". The way I think of the model is that there are systematic differences in the mean of the response between machines (the fixed effects); then a random deviation for each worker (random intercepts/worker); then a random deviation for each machine-worker combination; and finally a random deviation per observation. The greater the variance of the random deviations per worker, the more correlated observations from a given worker would be, etc.
Nate Pope