Понимание и исправление ошибок ложной сходимости в glmmTMB + lme4 - PullRequest
Новая
       118

Понимание и исправление ошибок ложной сходимости в glmmTMB + lme4

1 голос
/ 20 июня 2020

Я пытаюсь использовать glmmTMB для выполнения нескольких итераций моей модели, но продолжаю получать одну и ту же постоянную ошибку. Я попытался объяснить свой эксперимент ниже и вставил полную модель, которую пытаюсь запустить.

Предпосылки эксперимента

Зависимой переменной, которую я пытаюсь смоделировать, является количество копий бактериального гена 16S, в данном случае используется в качестве прокси для бактериальной биомассы.

Экспериментальный план состоит в том, что у меня есть отложения из 8 потоков, которые падают по градиенту загрязнения (до нетронутого). (Фактор 1 = поток, с 8 уровнями).

Для каждого из 8 потоков было выполнено следующее: Осадок был добавлен в 6 лотков. 3 из этих лотков были помещены в искусственный канал, нагретый до 13 ° C, в то время как остальные 3 были нагреты до 17 ° C (фактор 2 = обработка с подогревом, с 2 уровнями). Всего имеется 16 каналов, и обработка Warming была назначена каналу случайным образом.

Затем я неоднократно измерял 3 лотка в каждом канале потока в течение четырех временных точек (фактор 3 = день, с 4 уровнями).

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

Итак, чтобы подвести итог: условия модели (все указаны как факторы):

  • Обработка тепла (13 против 17o C)
  • StreamID ( 1,2,3,4,5,6,7,8)
  • День (T1, T4, T7, T14)

Полная модель, которую я предлагал, была, 

X4_tmb.nb2<-glmmTMB(CopyNo~Treatment*Stream*Time, family=nbinom2, data=qPCR)

Несмотря на то, что эта версия модели не включает случайный эффект, я хотел использовать пакет glmmTMB, а не запускать его с помощью lme4, потому что я хотел изучить идею добавления компонент модели для учета дисперсии, а также изучить возможность добавления лотка в качестве случайного эффекта (не уверен, правильно ли это). Запустив все версии модели в glmmTMB, я могу с уверенностью сравнить их оценки AI C. Я бы не смог этого сделать, если бы я запустил полную модель без компонента дисперсии в lme4, а остальные с glmmTMB.

К сожалению, для большинства итераций полной модели при использовании glmmTMB (я имею в виду последовательное удаление членов модели) я получаю одно и то же постоянное предупреждение:

Предупреждение: In fitTMB (TMBStru c): проблема сходимости модели; ложная сходимость (8). См. Виньетку («устранение неполадок»)

Я пытался понять ошибку, но мне трудно понять, потому что сбивает с толку то, что когда я запускаю полную модель с использованием lme4, она работает без ошибок.

Это версия полной модели, которая работает в lme4, 

X4 = glm.nb(CopyNo~Treatment*Stream*Time, data = qPCR

Насколько я понимаю из чтения https://www.biorxiv.org/content/10.1101/132753v1.full.pdf @ строка 225, это можно использовать этот пакет для перекрестного сравнения GLM и GLMM. Вы знаете, правильно ли я это понял?

Я также использовал пакет DHARMa, чтобы помочь проверить модели и версию, которая не смогла сойтись, используя glmmTMB, пройти KStest, тест дисперсии, тест выбросов и комбинированный тест скорректированного квантиля, но в идеале Я не хочу ошибки сходимости.

Любая помощь будет принята с благодарностью.

1 Ответ

0 голосов
/ 20 июня 2020

Здесь куча.

предупреждающее сообщение

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

 control = glmmTMBControl(optimizer = optim, optArgs = list(method="BFGS"))

в свой вызов. Надеюсь, это даст очень похожий ответ (в этом случае вы сделаете вывод, что предупреждение о конвергенции, вероятно, было ложным, поскольку разные оптимизаторы вряд ли потерпят неудачу одинаково) без предупреждения. (Вы можете попробовать method="CG" выше в качестве третьей альтернативы.) (Обратите внимание, что есть небольшая ошибка с печатью и суммированием вывода при использовании альтернативных оптимизаторов, которая была недавно исправлена ​​; вы можете установить версию для разработки, если вы работаете над этим до того, как исправление будет распространено на CRAN.)

модель «lme4»

Функция glm.nb() - это не из пакета lme4, это из пакета MASS. Если бы у вас были случайные эффекты в модели, вы бы использовали glmer.nb(), что равно в пакете lme4 ... как и в тестах переключения оптимизатора выше, если вы получите аналогичные ответы с glmmTMB и glm.nb вы можете сделать вывод, что предупреждение от glmmTMB (на самом деле, это от оптимизатора nlminb(), который glmmTMB вызывает внутренне), вероятно, является ложным срабатыванием.

Самый простой способ проверить соизмеримость правдоподобия / AIC из разных пакетов - по возможности уместить одну и ту же модель в оба пакета, например, 

library(MASS)
library(glmmTMB)

quine.nb1 <- glm.nb(Days ~ Sex/(Age + Eth*Lrn), data = quine)
quine.nb2 <- glmmTMB(Days ~ Sex/(Age + Eth*Lrn), data = quine,
                     family=nbinom2)
all.equal(AIC(quine.nb1),AIC(quine.nb2))  ## TRUE

другие детали

Одна из возможных проблем вашей модели заключается в том, что, подбирая полное трехстороннее взаимодействие трех категориальных переменных, вы пытаетесь оценить (2 * 4 * 8 =) 64 параметра до 64 * 3 = 192 наблюдений (если я правильно понять ваш экспериментальный план). Это более вероятно приведет к численным проблемам (как указано выше) и может дать неточные результаты. Хотя некоторые рекомендуют это, я не лично рекомендую подход к выбору модели (все подмножества или последовательный, на основе AI C или на основе p-значения); Я бы предложил преобразовать Stream в случайный эффект, т.е. 

CopyNo ~ Treatment + (Treatment|StreamID) + (1|Time/StreamID)

Это соответствует (1) общему лечебному эффекту; (2) изменение по потокам, изменение эффекта обработки по потокам; (3) изменение во времени и между потоками в определенные моменты времени. Здесь используется только 2 (фиксировано: перехват + обработка) + 3 (дисперсия перехвата, дисперсия обработки и их ковариация по потокам) + 2 (дисперсия во времени и между потоками во времени). Это не совсем «максимальная» модель; поскольку обе обработки измеряются каждый раз в каждом потоке, последний термин может быть (Обработка | Время / ID потока), но это значительно усложнит модель. Поскольку 4 группы - это не так много для случайного эффекта, вы можете обнаружить, что хотите превратить последний член в Time + (1|Time:StreamID), который будет соответствовать времени как фиксированному, а (потоки во времени) как случайному ...

...