Что такое R-структура G-структура в глмм?

16

Я недавно использовал MCMCglmmпакет. Меня смущает то, что в документации упоминается как R-структура и G-структура. Похоже, что они связаны со случайными эффектами - в частности, указанием параметров для предварительного распределения по ним, но обсуждение в документации, похоже, предполагает, что читатель знает, что это за термины. Например:

необязательный список предыдущих спецификаций, имеющий 3 возможных элемента: R (R-структура) G (G-структура) и B (фиксированные эффекты) ............ Приоритеты для структур дисперсии (R и G ) являются списками с ожидаемыми (со) дисперсиями (V) и параметром степени достоверности (nu) для обратного Вишарта

... взяты отсюда .

РЕДАКТИРОВАТЬ: Обратите внимание, что я переписал остальную часть вопроса после комментариев от Стефана.

Может ли кто-нибудь пролить свет на то, что такое R-структура и G-структура, в контексте простой модели компонентов дисперсии, где линейный предиктор равен с e 0 i jN ( 0 , σ 2 0 e ) и u 0 jN ( 0 , σ 2 0 u )

β0+e0ij+u0j
e0ijN(0,σ0e2)u0jN(0,σ0u2)

Я сделал следующий пример с некоторыми данными, которые идут с MCMCglmm

> require(MCMCglmm)
> require(lme4)
> data(PlodiaRB)
> prior1 = list(R = list(V = 1, fix=1), G = list(G1 = list(V = 1, nu = 0.002)))
> m1 <- MCMCglmm(Pupated ~1, random = ~FSfamily, family = "categorical", 
+ data = PlodiaRB, prior = prior1, verbose = FALSE)
> summary(m1)


 G-structure:  ~FSfamily

         post.mean l-95% CI u-95% CI eff.samp
FSfamily    0.8529   0.2951    1.455      160

 R-structure:  ~units

      post.mean l-95% CI u-95% CI eff.samp
units         1        1        1        0

 Location effects: Pupated ~ 1 

            post.mean l-95% CI u-95% CI eff.samp  pMCMC    
(Intercept)   -1.1630  -1.4558  -0.8119    463.1 <0.001 ***
---

> prior2 = list(R = list(V = 1, nu = 0), G = list(G1 = list(V = 1, nu = 0.002)))
> m2 <- MCMCglmm(Pupated ~1, random = ~FSfamily, family = "categorical", 
+ data = PlodiaRB, prior = prior2, verbose = FALSE)
> summary(m2)


 G-structure:  ~FSfamily

         post.mean l-95% CI u-95% CI eff.samp
FSfamily    0.8325   0.3101    1.438    79.25

 R-structure:  ~units

      post.mean l-95% CI u-95% CI eff.samp
units    0.7212  0.04808    2.427    3.125

 Location effects: Pupated ~ 1 

            post.mean l-95% CI u-95% CI eff.samp  pMCMC    
(Intercept)   -1.1042  -1.5191  -0.7078    20.99 <0.001 ***
---

> m2 <- glmer(Pupated ~ 1+ (1|FSfamily), family="binomial",data=PlodiaRB)
> summary(m2)
Generalized linear mixed model fit by the Laplace approximation 
Formula: Pupated ~ 1 + (1 | FSfamily) 
   Data: PlodiaRB 
  AIC  BIC logLik deviance
 1020 1029   -508     1016
Random effects:
 Groups   Name        Variance Std.Dev.
 FSfamily (Intercept) 0.56023  0.74849 
Number of obs: 874, groups: FSfamily, 49

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -0.9861     0.1344  -7.336  2.2e-13 ***

Поэтому, основываясь на комментариях Стефана, я думаю, что структура G для . Но в комментариях также говорится, что структура R предназначена для σ 2 0 e, однако, похоже, что она не появляется в выходных данных.σ0u2σ0e2lme4

Обратите внимание, что результаты lme4/glmer()согласуются с обоими примерами из MCMC MCMCglmm.

Итак, является ли структура R для и почему она не появляется в выходных данных для ?σ0e2lme4/glmer()

Джо Кинг
источник
1
С помощью терминологии SAS (но, возможно, это более распространенная терминология), матрица G является матрицей дисперсии случайных эффектов, а матрица R является матрицей дисперсии «слагаемых ошибок» (в вашем случае, возможно, это предполагаемый остаток дисперсия ?)σ0e2
Стефан Лоран
@ StéphaneLaurent спасибо. Я задавался вопросом, можно ли оценить но когда я впервые узнал об обобщенной линейной модели, я помню, что σ 2 0 e не оценивается - вычисляется только «отклонение» (как с ). Может я что-то упустил? σ0e2σ0e2lme4
Джо Кинг,
1
может быть, смысл остаточной дисперсии неясен, когда семейство распределений не гауссово
Стефан Лоран
1
@ Стефан Лоран Да! Пожалуйста, смотрите мой комментарий к ответу Майкла минуту назад - для двоичного результата, это должно быть исправлено (как в моих моделях в моем OP)
Джо Кинг
1
Если у вас есть модель ME / Multilevel, есть несколько отклонений. Представьте простейший случай: . Существует разница в пересечениях b i и в ошибочном члене ε i . G часто используется для матрицы вар-коваров случайных эффектов (в данном случае скаляр σ 2 b ), а R i - для матрицы вар-коваров остаточных дисперсий ε i.Yi=β0+β1X+bi+εibiεiGσb2Riεiпосле учета фиксированных и случайных эффектов этого кластера. Это, как правило , понимается как диагональная матрица «s. Кроме того, обе стороны рассматриваются как многовариантные нормальные значения w / mean = 0. σ2
gung - Восстановить Монику

Ответы:

8

Я бы предпочел опубликовать свои комментарии ниже как комментарий, но этого будет недостаточно. Это скорее вопросы, чем ответ (по сравнению с @gung я не достаточно силен в этой теме).

У меня сложилось впечатление, что MCMCglmm не реализует «настоящий» байесовский блеск. Истинная байесовская модель описана в разделе 2 этой статьи . Как и в случае модели с частотой, каждый имеет и для параметра дисперсии ϕ 1 требуется априор в дополнение к фиксированным параметрам β и дисперсии «G» случайным образом эффект у .g(E(yu))=Xβ+Zuϕ1βu

Но в соответствии с этой MCMCglmm виньетка , модель реализована в MCMCglmm задается , и он не связан с параметром дисперсии φ 1 . Это не похоже на классическую модель частых.g(E(yu,e))=Xβ+Zu+eϕ1

Поэтому я бы не удивился , что нет аналога с glmer.σe

Пожалуйста, извинитесь за эти грубые комментарии, я просто взглянул на это.

Стефан Лоран
источник
σeglmerMCMCglmmMCMCglmmϕ1
Извините, мои слова были не совсем уместны. MCMCglmm действительно байесовский, но он не совсем реализует классический glmm (я думаю). Кроме того, вы должны знать, что трудно установить априорные значения, дающие вывод по компонентам дисперсии, близким к выводу по частоте.
Стефан Лоран
Еще раз спасибо. В моем исследовании я обнаружил, что могу использовать стандартное распределение обратных желаний для компонентов дисперсии при MCMCglmmиспользовании множества параметров, и доверительные интервалы 95% всегда содержат значение дисперсии для оценки случайных эффектов, glmerпоэтому я чувствовал, что это разумно , но как мне интерпретировать этот случай, который может быть нетипичным, в результате чего MCMCglmmинтервалы не очень чувствительны к выбору предшествующего? Может быть, я должен задать новый вопрос по этому поводу?
Джо Кинг,
σe=0σe0
σeσeσuglmer
11

Я опоздал на игру, но несколько заметок. рструктура - это остаточная структура. В вашем случае «структура» имеет только один элемент (но это не обязательно так). Для гауссовой переменной отклика остаточная дисперсияσе2обычно оценивается. Для бинарных результатов он остается постоянным. Из-за того, как настроен MCMCglmm , вы не можете исправить это в нуле, но это относительно стандартно, чтобы исправить это в1(также верно для пробной модели). Для данных подсчета (например, с распределением Пуассона) вы не фиксируете это, и это автоматически автоматически оценивает параметр избыточной дисперсии.

граммструктура - это структура случайных эффектов. Опять же, в вашем случае, просто случайный перехват, но если бы у вас было несколько случайных эффектов, они бы сформировали матрицу дисперсии-ковариации,грамм,

И последнее замечание: поскольку остаточная дисперсия не зафиксирована на нуле, оценки не будут совпадать с оценками из glmer. Вы должны изменить их масштаб. Вот небольшой пример (не использующий случайные эффекты, но обобщающий). Обратите внимание, как дисперсия структуры R фиксируется на 1.

# example showing how close the match is to ML without separation
m2 <- MCMCglmm(vs ~ mpg, data = mtcars, family = "categorical",
  prior = list(
    B = list(mu = c(0, 0), V = diag(2) * 1e10),
    R = list(V = 1, fix = 1)),
  nitt = 1e6, thin = 500, burnin = 10000)
summary(m2)

Вот константа масштабирования для биномиального семейства:

k <- ((16*sqrt(3))/(15*pi))^2

Теперь разделите решение на него и получите задние моды

posterior.mode(m2$Sol/(sqrt(1 + k)))

Что должно быть достаточно близко к тому, что мы получаем от glm

summary(glm(vs ~mpg, data = mtcars, family = binomial))
Джошуа
источник
Вы случайно не знаете, как указать гетероскедастичность на первом уровне в MCMCglmm? Это структура R? Какой синтаксис тогда?
Maxim.K
@ Джошуа, ты можешь объяснить «константу пересчета для биномиального семейства»? PS: для семян 123я получаю (с поправкой) от m2значений -8.164и 0.421; и из glmзначений -8.833и 0.430.
Касвед
Константу масштабирования можно найти в Diggle et. и др. ( amazon.de/Analysis-Longtial-Oxford-Statistical-Science/dp/… ) - в соответствии с cran.r-project.org/web/packages/MCMCglmm/vignettes/… экв. 2.14 на стр. 47.
Ксвед