Биномиальный блеск с категориальной переменной с полным успехом

Я использую glmm с биномиальной переменной ответа и категориальным предиктором. Случайный эффект дается вложенным дизайном, используемым для сбора данных. Данные выглядят так:

m.gen1$treatment
 [1] sucrose      control      protein      control      no_injection .....
Levels: no_injection control sucrose protein
m.gen1$emergence 
 [1]  1  0  0  1  0  1  1  1  1  1  1  0  0....
> m.gen1$nest
 [1] 1  1  1  2  2  3  3  3  3  4  4  4  .....
Levels: 1 2 3 4 5 6 8 10 11 13 15 16 17 18 20 22 24

Первая модель, которую я запускаю, выглядит так

m.glmm.em.<-glmer(emergence~treatment + (1|nest),family=binomial,data=m.gen1)

Я получаю два предупреждения, которые выглядят так:

Warning messages:
1: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
  Model failed to converge with max|grad| = 0.0240654 (tol = 0.001, component 4)
2: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
  Model is nearly unidentifiable: large eigenvalue ratio
 - Rescale variables?

Краткое описание модели показывает, что одна из обработок имеет необычно большую стандартную ошибку, которую вы можете увидеть здесь:

Fixed effects:
                 Estimate Std. Error z value Pr(>|z|)  
(Intercept)         2.565      1.038   2.472   0.0134 *
treatmentcontrol   -1.718      1.246  -1.378   0.1681  
treatmentsucrose   16.863   2048.000   0.008   0.9934  
treatmentprotein   -1.718      1.246  -1.378   0.1681 

Я попробовал разные оптимизаторы из glmer control и функции из других пакетов, и я получил похожий вывод. Я запустил модель, используя glm, игнорируя случайный эффект, и проблема сохраняется. Изучая данные, я понял, что лечение с высоким Std. Ошибка имеет только успехи в переменной ответа. Просто чтобы проверить, может ли это быть причиной проблемы, я добавил ложную точку данных с «неудачей» для этой обработки, и модель работает без сбоев и выдает разумную стандартную ошибку. Вы можете увидеть это здесь:

Fixed effects:
                 Estimate Std. Error z value Pr(>|z|)  
(Intercept)        3.4090     1.6712   2.040   0.0414 *
treatmentcontrol  -1.8405     1.4290  -1.288   0.1978  
treatmentsucrose  -0.2582     1.6263  -0.159   0.8738  
treatmentprotein  -2.6530     1.5904  -1.668   0.0953 .

Мне было интересно, верна ли моя интуиция насчет отсутствия неудач в этом лечении, мешающих хорошей оценке, и как я могу обойти эту проблему.

Заранее спасибо!

Ваша интуиция абсолютно правильна. Это явление называется полным разделением . Вы можете найти довольно много (теперь, когда вы знаете его имя), гуглят вокруг ... Это довольно подробно обсуждается здесь, в общем контексте , и здесь, в контексте GLMM . Стандартным решением этой проблемы является добавление небольшого члена, который сдвигает параметры обратно к нулю - в частых контекстах это называется штрафным или исправленным смещением методом. Стандартный алгоритм принадлежит Ферту (1993, «Уменьшение смещения оценок максимального правдоподобия», Biometrika 80, 27-38) и реализован в пакете logistf.на КРАН. В байесовском контексте это оформлено как добавление слабого перед параметрами с фиксированным эффектом.

Насколько мне известно, алгоритм Фёрта не был распространен на GLMM, но вы можете использовать байесовский трюк с помощью пакета blme , который накладывает тонкий байесовский слой поверх lme4пакета. Вот пример из вышеупомянутого обсуждения GLMM:

cmod_blme_L2 <- bglmer(predation~ttt+(1|block),data=newdat,
                   family=binomial,
                   fixef.prior = normal(cov = diag(9,4)))

Первые две строки в этом примере точно такие же, как мы использовали бы в стандартной glmerмодели; последнее указывает, что априором для фиксированных эффектов является многомерное нормальное распределение с диагональной матрицей дисперсии-ковариации. Матрица имеет размер 4x4 (поскольку в этом примере у нас есть 4 параметра с фиксированным эффектом), а предыдущая дисперсия каждого параметра равна 9 (что соответствует стандартному отклонению 3, что довольно мало - это означает, что +/- 2SD равно ( -6,6), что является очень большим диапазоном по шкале логита).

Очень большие стандартные ошибки параметров в вашем примере являются примером явления, тесно связанного с полным разделением (это происходит всякий раз, когда мы получаем экстремальные значения параметров в логистической модели), называемого эффектом Хаука-Доннера .

Еще две потенциально полезные ссылки (я еще не копался в них):

• Гельман А., Якулин А., Питтау М.Г. и Су Т.С. (2008) Слабоинформативное предварительное распределение по умолчанию для логистических и других регрессионных моделей. Летопись прикладной статистики , 2, 1360–383.
• José Cortiñas Abrahantes и Marc Aerts (2012) Решение для разделения для кластеризованных двоичных данных Статистическое моделирование 12 (1): 3–27 doi: 10.1177 / 1471082X1001200102

Более поздний поиск ученого Google для "bglmer 'полное разделение'" находит:

• Quiñones, AE и WT Wcislo. «Криптический расширенный уход за потомством в факультативно эусоциальной потовой пчеле Megalopta genalis ». Insectes Sociaux 62,3 (2015): 307–313.