MCMCglmm мультиномиальная модель R

Question

Я пытаюсь создать модель с помощью MCMCglmm пакета в R.

Данные структурированы следующим образом, где двойка, очаговые, другие являются все случайные эффекты, predict1-2 являются предиктором переменные и ответа 1-5 являются итоговыми переменными, которые захватывают # наблюдаемого поведение различных подтипов:

dyad focal other r present village resp1 resp2 resp3 resp4 resp5 
1 10101 14302 0.5 3  1  0  0  4  0  5 
2 10405 11301 0.0 5  0  0  0  1  0  1 
… 

так модель только один исход (обучение) выглядит следующим образом:

prior_overdisp_i <- list(R=list(V=diag(2),nu=0.08,fix=2), 
G=list(G1=list(V=1,nu=0.08), G2=list(V=1,nu=0.08), G3=list(V=1,nu=0.08), G4=list(V=1,nu=0.08))) 

m1 <- MCMCglmm(teaching ~ trait-1 + at.level(trait,1):r + at.level(trait,1):present, 
random= ~idh(at.level(trait,1)):focal + idh(at.level(trait,1)):other + 
idh(at.level(trait,1)):X + idh(at.level(trait,1)):village, 
rcov=~idh(trait):units, family = "zipoisson", prior=prior_overdisp_i, 
data = data, nitt = nitt.1, thin = 50, burnin = 15000, pr = TRUE, pl = TRUE, verbose = TRUE, DIC = TRUE) 

Замечания курса Хадфилда (Ch 5) приводят пример мультиномиальной модели, которая использует только одну переменную результата с тремя уровнями (овечья рожки 3-х типов). Подобное лечение можно найти здесь: http://hlplab.wordpress.com/2009/05/07/multinomial-random-effects-models-in-r/ Это не совсем так, как я делаю, но содержит полезную справочную информацию.

Другая ссылка (Hadfield 2010) дает пример MCMCglmm с несколькими ответами, который следует за одним и тем же форматом, но использует cbind() для прогнозирования вектора ответов, а не для одного результата. Та же модель с множественными ответами будет выглядеть следующим образом:

m1 <- MCMCglmm(cbind(resp1, resp2, resp3, resp4, resp5) ~ trait-1 + 
at.level(trait,1):r + at.level(trait,1):present, 
random= ~idh(at.level(trait,1)):focal + idh(at.level(trait,1)):other + 
idh(at.level(trait,1)):X + idh(at.level(trait,1)):village, 
rcov=~idh(trait):units, 
family = cbind("zipoisson","zipoisson","zipoisson","zipoisson","zipoisson"), 
prior=prior_overdisp_i, 
data = data, nitt = nitt.1, thin = 50, burnin = 15000, pr = TRUE, pl = TRUE, verbose = TRUE, DIC = TRUE) 

У меня есть два вопроса здесь: Программирование Радио

1. Как указать предшествующий для этой модели? Я просмотрел материалы, упомянутые в этом посте, но просто не могу понять это.

2. У меня есть аналогичная версия с двумя переменными ответа, но я получаю только один наклон - там, где я думал, что должен получить другой уклон для каждой переменной resp. Где я ошибаюсь, или я неправильно понял модель?

Answer 1

Ответ на мой первый вопрос, на основе HLP пост и некоторую помощь от Colleage/статистики консультанта:

# values for prior 
k <- 5 # originally: length(levels(dative$SemanticClass)), so k = # of outcomes for SemanticClass  aka categorical outcomes 
I <- diag(k-1) #should make matrix of 0's with diagonal of 1's, dimensions k-1 rows and k-1 columns 
J <- matrix(rep(1, (k-1)^2), c(k-1, k-1)) # should make k-1 x k-1 matrix of 1's 

И для моей модели, используя multinomial5 семьи и 5 итоговых переменных до это:

prior = list(
      R = list(fix=1, V=0.5 * (I + J), n = 4), 
      G = list(
       G1 = list(V = diag(4), n = 4)) 

на мой второй вопрос, мне нужно добавить член взаимодействия с фиксированными эффектами в этой модели:

m <- MCMCglmm(cbind(Resp1, Resp2...) ~ -1 + trait*predictorvariable, 
... 

Результат дает как основные эффекты для переменных ответа, так и задние оценки для взаимодействия Response/Predictor (влияние переменной-предиктора на каждую переменную ответа).