Эквивалентность (0 + фактор | группа) и (1 | группа) + (1 | группа: фактор) характеристики случайных эффектов в случае составной симметрии

13

Дуглас Бейтс утверждает, что следующие модели эквивалентны «если матрица дисперсии-ковариации для векторно-значных случайных эффектов имеет особую форму, называемую составной симметрией» ( слайд 91 в этой презентации ):

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

В частности, Бейтс использует этот пример:

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

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

с соответствующими выходами:

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

Кто-нибудь может объяснить разницу между моделями и как m1сводится к m2(с учетом составной симметрии) интуитивно понятным способом?

statmerkur
источник
6
+1 и, имхо, это абсолютно по теме. Проголосовать, чтобы открыть.
говорит амеба, восстанови Монику
2
@ Питер Флом, почему вы считаете этот вопрос не по теме?
statmerkur
3
Возможно, не было ясно, что вы спрашивали о моделях, а не о lme4синтаксисе. Было бы полезно - и расширить круг потенциальных ответчиков - если бы вы объяснили их людям, с которыми незнакомы lme4.
Scortchi - Восстановить Монику
Похоже, речь идет о кодировании.
Питер Флом - Восстановить Монику
1
Если это полезно, вот два хороших поста о том, что делает синтаксис lme4, и что такое составная симметрия в контексте смешанных моделей (см. Принятые ответы на оба вопроса). stats.stackexchange.com/questions/13166/rs-lmer-cheat-sheet и stats.stackexchange.com/questions/15102/…
Джейкоб Соколар

Ответы:

11

В этом примере есть три наблюдения для каждой комбинации трех машин (A, B, C) и шести рабочих. Я буду использовать для обозначения n- мерной единичной матрицы и 1 n для обозначения n- мерного вектора единиц. Скажем, у - вектор наблюдений, который, как я предполагаю, упорядочен рабочим, затем машиной, а затем воспроизведен. Пусть µ будет соответствующими ожидаемыми значениями (например, фиксированными эффектами), и пусть γ будет вектором групповых отклонений от ожидаемых значений (например, случайных эффектов). Условно на γ модель для y можно записать так:Inn1nnyμγγy

yN(μ+γ,σy2I54)

где - «остаточная» дисперсия.σy2

Чтобы понять, как ковариационная структура случайных эффектов индуцирует ковариационную структуру среди наблюдений, более интуитивно понятно работать с эквивалентным «маргинальным» представлением , которое интегрируется по случайным эффектам . Предельная форма этой модели,γ

yN(μ,σy2I54+Σ)

Здесь - ковариационная матрица, которая зависит от структуры γ (например, «компоненты дисперсии», лежащие в основе случайных эффектов). Я буду называть Σ «маргинальной» ковариацией.ΣγΣ

По вашему m1, случайные эффекты разлагаются на:

γ=Zθ

Там , где является дизайн матрица , которая отображает случайные коэффициенты на наблюдения, а θ Т = [ θ 1 , , θ 1 , B , θ 1 , С ... θ 6 , , & thetas ; 6 , B , θ 6 , с ] представляет собой 18-мерный вектор случайных коэффициентов упорядоченных по рабочим затем машина, и распределяются следующим образом: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]

(Другими словами, корреляционная матрица, лежащая в основе имеет все элементы на диагонали, установленные на одно и то же значение.)Λ

Предельная структура ковариации , вызванная этими случайных эффектов является , так что дисперсия данного наблюдения является σ 2 & thetas ; + т 2 + сг 2 г и ковариации между двумя (отдельными) наблюдениями от работников i , j и машин u , v is: c o v ( y i , u , y j , v ) =Σ=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

Для вас m2случайные эффекты разлагаются на:

γ=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

σθ2σω2τ2ση2 m1

Краткость не является моей сильной стороной: это всего лишь длинный, запутанный способ сказать, что каждая модель имеет два параметра дисперсии для случайных эффектов, и это просто два разных способа написания одной и той же «маргинальной» модели.

В коде ...

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
Нейт Папа
источник
1
Очень хороший ответ! Но я думаю, что фраза «машина, вложенная в работника» может вводить в заблуждение, поскольку одни и те же три машины появляются на более чем одном (фактически, каждом) уровне работника.
statmerkur
@statmerkur Спасибо, я попытался уточнить эту строку. Дайте мне знать, если у вас есть другое предложение.
Нейт Папа
1
Должен Икс определяться как Иксзнак равноя619?
С. Каттералл восстановил Монику
1
@S.Catterall Yup, that's a typo -- thanks for catching it! I've fixed in my answer.
Nate Pope
2
@statmerkur Вы можете уточнить, что вы имеете в виду? Здесь нет непрерывных ковариат, поэтому не уверен, что вы подразумеваете под «уклоном». То, как я думаю о модели, заключается в том, что существуют систематические различия в среднем отклике между машинами (фиксированные эффекты); затем случайное отклонение для каждого работника (случайные перехваты / работник); затем случайное отклонение для каждой комбинации машина-рабочий; и, наконец, случайное отклонение за наблюдение. Чем больше дисперсия случайных отклонений на одного работника, тем больше будет коррелированных наблюдений от данного работника и т. Д.
Нейт Папа