Я использую биномиальную логистику c с использованием lme4, которая построена следующим образом:
m2_DL <- glmer(Dest ~ Pb_on * Pb_type + Dest_tree + ZF_Dest + (1|Site), data = DF_DL_3, family = binomial(link = "logit"))
Все переменные являются коэффициентами, кроме ZF_Dest, который является числовым. Полная модель (m1) включает в себя одну дополнительную числовую переменную и работает нормально. Эта вторая модель, однако, не сходится с предупреждающим сообщением:
fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
Warning message:
In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00144561 (tol = 0.001, component 1)
В моем случае ожидается недостаток ранга и (я полагаю) не имеет отношения к проблеме здесь. Посмотрев онлайн, я обнаружил, что могу попытаться добавить следующий аргумент к своему блеску: control = glmerControl(optimizer = "bobyqa"), что действительно решает проблему.
Однако я хотел бы лучше понять, что делает эта функция, и если бы я все еще мог сравнивать свои модели, если только моя вторая модель имеет этот дополнительный аргумент. Это привело меня к этому вопросу на SO . В одном из комментариев к этому вопросу говорится: «Вы не использовали другой оптимизатор (по умолчанию для glmer используется bobyqa)». Но если это так, то почему добавление аргумента разрешило предупреждение, поскольку я не увеличивал количество итераций, как это иногда делается?
Из ответа на этот вопрос я понял, что если разность вероятностей обеих моделей близка, я указываю «ложный положительный результат» ошибки сходимости, что опять-таки имеет место для моих моделей. Я также прочитал help(converge) в соответствии с рекомендациями и использовал указанную там функцию allFit. Это выводит:
bobyqa : fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
[OK]
Nelder_Mead : fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
[OK]
nlminbwrap : fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
[OK]
nmkbw : fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
[OK]
optimx.L-BFGS-B : fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
[OK]
nloptwrap.NLOPT_LN_NELDERMEAD : fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
[OK]
nloptwrap.NLOPT_LN_BOBYQA : fixed-effect model matrix is rank deficient so dropping 6 columns / coefficients
[OK]
Warning messages:
1: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0013785 (tol = 0.001, component 1)
2: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00904036 (tol = 0.001, component 1)
3: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00904036 (tol = 0.001, component 1)
, и объект, содержащий результаты allFit, содержит:
original model:
Dest ~ Pb_on * Pb_type + Dest_tree + ZF_Dest + (1 | Site)
data: DF_DL_3
optimizers (7): bobyqa, Nelder_Mead, nlminbwrap, nmkbw, optimx.L-BFGS-B, nloptwrap.NLOPT_LN_N...
1 optimizer(s) failed
differences in negative log-likelihoods:
max= 6.65e-05 ; std dev= 3.41e-05
Здесь вы видите, что 1 оптимизатор не прошел, и я получаю несколько предупреждений при использовании allFit. Теперь я понимаю, что не все оптимизаторы сделаны равными, и некоторые не будут работать для определенных данных. Проблема в том, что я не понимаю, как интерпретировать результат этой функции allFit, и я не уверен, какой оптимизатор дал сбой и стоит ли мне беспокоиться. Итак, в заключение, мои вопросы:
Как оптимизатор BOBYQA мог решить проблему, если это действительно оптимизатор по умолчанию?
Как мне интерпретировать результат allFit.
Каков наилучший способ устранения предупреждения для второй модели. Является ли просто добавление control = glmerControl(optimizer = "bobyqa") подходящим?
Если я добавлю этот аргумент к этой модели, могу ли я тогда еще честно сравнить эту модель с другими? (Примечание: я использую функцию anova для сравнения моделей).
Спасибо
РЕДАКТИРОВАТЬ 6-3-2020: В соответствии с запросом Роланда вот сводка объекта allFit:
$which.OK
bobyqa Nelder_Mead nlminbwrap nmkbw
TRUE TRUE TRUE TRUE
optimx.L-BFGS-B nloptwrap.NLOPT_LN_NELDERMEAD nloptwrap.NLOPT_LN_BOBYQA
FALSE TRUE TRUE
$msgs
$msgs$bobyqa
NULL
$msgs$Nelder_Mead
NULL
$msgs$nlminbwrap
NULL
$msgs$nmkbw
[1] "Model failed to converge with max|grad| = 0.0013785 (tol = 0.001, component 1)"
$msgs$nloptwrap.NLOPT_LN_NELDERMEAD
[1] "Model failed to converge with max|grad| = 0.00904036 (tol = 0.001, component 1)"
$msgs$nloptwrap.NLOPT_LN_BOBYQA
[1] "Model failed to converge with max|grad| = 0.00904036 (tol = 0.001, component 1)"
$fixef
(Intercept) Pb_onL Pb_onN Pb_onS Pb_typeNG Pb_typeSong Dest_treeL ZF_Dest
bobyqa 0.4421569 0.5633434 0.1788602 -0.2046165 0.1777216 -0.1115736 0.5845617 0.01578428
Nelder_Mead 0.4420709 0.5633071 0.1788488 -0.2046262 0.1777301 -0.1115980 0.5845859 0.01578563
nlminbwrap 0.4421020 0.5633365 0.1788603 -0.2046102 0.1777225 -0.1115746 0.5845682 0.01578510
nmkbw 0.4423417 0.5635447 0.1791150 -0.2045143 0.1781367 -0.1114576 0.5844313 0.01578639
nloptwrap.NLOPT_LN_NELDERMEAD 0.4434285 0.5612392 0.1765273 -0.2063090 0.1773619 -0.1140108 0.5852130 0.01579032
nloptwrap.NLOPT_LN_BOBYQA 0.4434285 0.5612392 0.1765273 -0.2063090 0.1773619 -0.1140108 0.5852130 0.01579032
Pb_onL:Pb_typeNG Pb_onN:Pb_typeNG Pb_onL:Pb_typeSong Pb_onN:Pb_typeSong
bobyqa -1.140386 -0.1635209 -0.5446501 -0.2497627
Nelder_Mead -1.140380 -0.1635345 -0.5445915 -0.2497404
nlminbwrap -1.140389 -0.1635264 -0.5446236 -0.2497616
nmkbw -1.140884 -0.1640907 -0.5448729 -0.2500731
nloptwrap.NLOPT_LN_NELDERMEAD -1.138619 -0.1624441 -0.5413400 -0.2461967
nloptwrap.NLOPT_LN_BOBYQA -1.138619 -0.1624441 -0.5413400 -0.2461967
$llik
bobyqa Nelder_Mead nlminbwrap nmkbw
-826.3143 -826.3143 -826.3143 -826.3143
nloptwrap.NLOPT_LN_NELDERMEAD nloptwrap.NLOPT_LN_BOBYQA
-826.3143 -826.3143
$sdcor
Site.(Intercept)
bobyqa 1.537530
Nelder_Mead 1.537517
nlminbwrap 1.537501
nmkbw 1.537453
nloptwrap.NLOPT_LN_NELDERMEAD 1.538042
nloptwrap.NLOPT_LN_BOBYQA 1.538042
$theta
Site.(Intercept)
bobyqa 1.537530
Nelder_Mead 1.537517
nlminbwrap 1.537501
nmkbw 1.537453
nloptwrap.NLOPT_LN_NELDERMEAD 1.538042
nloptwrap.NLOPT_LN_BOBYQA 1.538042
$times
user.self sys.self elapsed user.child sys.child
bobyqa 6.51 0 6.55 NA NA
Nelder_Mead 3.06 0 3.06 NA NA
nlminbwrap 2.22 0 2.22 NA NA
nmkbw 2.17 0 2.17 NA NA
nloptwrap.NLOPT_LN_NELDERMEAD 0.94 0 0.94 NA NA
nloptwrap.NLOPT_LN_BOBYQA 0.86 0 0.86 NA NA
$feval
bobyqa Nelder_Mead nlminbwrap nmkbw
3098 1274 NA 672
nloptwrap.NLOPT_LN_NELDERMEAD nloptwrap.NLOPT_LN_BOBYQA
97 97
attr(,"class")
[1] "summary.allFit"
Использование car::vif(full_model) для проверки коллинеарности приводит к:
> car::vif(m1_DD)
GVIF Df GVIF^(1/(2*Df))
Pb_on 17.606990 3 1.612925
Pb_type 20.180212 2 2.119490
ZF_Dest 1.015237 1 1.007589
ZF_Other 1.006738 1 1.003364
Pb_on:Pb_type 96.910527 4 1.771317
Поскольку я не уверен на 100%, как их интерпретировать, я использовал эту ссылку, которая подсказывает, что я не должен не слишком переживай. В одном ответе говорится: «Правило GVIF (1 / (2 × Df)) <2 применяется в некоторых публикациях, что равняется обычному VIF 4 для однофакторных переменных». Таким образом, Pb_type может создать проблему, но я часто видел, что VIF 10 считается «приемлемым», поэтому я думаю, что это нормально. </p>
РЕДАКТИРОВАТЬ 10-3-2020 После еще немного поиска (теперь, когда я знаю, что на самом деле делает bobyqa, благодаря meriops !), я наткнулся на эту тему . Читая ответ там, я считаю, что мои ответы теперь выглядят следующим образом:
В соответствии с тем, что сказал meriops: «по умолчанию glmer использует двухэтапную оптимизацию с bobyqa для первой стадии. и Nelder-Mead для второго этапа. Поэтому, когда вы меняете оптимизатор на bobyqa, вы используете bobyqa для обоих этапов. "
Основано на ответе Бена Болкера в этой теме , и мои собственные вероятности отличаются <0,1, это может означать, что предупреждение было ложно положительным с самого начала. На данный момент я все еще хочу настроить оптимизаторы. </p>
Исходя из вышеупомянутой темы, я теперь полагаю, что это действительно хорошо.
Аналогично, вышеупомянутая тема предлагает этого не делать быть проблемой в моем случае. Так как мои модели действительно вложенные.
Если кто-либо из вас sh разработает ваши комментарии или внесет их в полные ответы, сделайте это, тогда я приму их.