Предупреждение GLMER: матрица дисперсии-ковариации [] не является положительно определенной или содержит значения NA

Иногда я обнаруживаю, что мои GLMM из glmer, пакет lme4, показывают следующие предупреждающие сообщения, когда вызывается их сводка:

Warning messages:
1: In vcov.merMod(object, use.hessian = use.hessian) :
  variance-covariance matrix computed from finite-difference Hessian is
not positive definite or contains NA values: falling back to var-cov estimated from RX
2: In vcov.merMod(object, correlation = correlation, sigm = sig) :
  variance-covariance matrix computed from finite-difference Hessian is
not positive definite or contains NA values: falling back to var-cov estimated from RX

Подобные вопросы, которые я нашел здесь, в Stackoverflow, относятся к другим функциям, а не к glmer, и LME4 Wiki не уточняет на том либо. В этом вопросе проблема была решена до появления таких сообщений об ошибках были рассмотрены, и здесь обсуждение сосредоточено на конкретной модели, а не на значении предупреждающего сообщения.

Итак, вопрос: должен ли я беспокоиться об этом сообщении, или это нормально, потому что это просто предупреждение, а не ошибка, и, как говорится, оно «откатывается к var-cov, оцененному из RX» (независимо от RX) тем не мение.

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

Вот (минимальный) набор данных:

testdata=structure(list(Site = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 
5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 
6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L), .Label = c("EO1", "EO2", 
"EO3", "EO4", "EO5", "EO6"), class = "factor"), Treatment = structure(c(1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L), .Label = c("control", 
"no ants", "no birds", "no birds no ants"), class = "factor"), 
    Tree = structure(c(2L, 3L, 4L, 16L, 12L, 13L, 14L, 15L, 5L, 
    6L, 7L, 8L, 1L, 9L, 10L, 11L, 28L, 29L, 30L, 31L, 17L, 25L, 
    26L, 27L, 18L, 19L, 20L, 32L, 21L, 22L, 23L, 24L, 33L, 41L, 
    42L, 43L, 37L, 38L, 39L, 40L, 44L, 45L, 46L, 47L, 34L, 35L, 
    36L, 48L, 49L, 57L, 58L, 59L, 50L, 51L, 52L, 64L, 53L, 54L, 
    55L, 56L, 60L, 61L, 62L, 63L, 66L, 67L, 68L, 80L, 69L, 70L, 
    71L, 72L, 76L, 77L, 78L, 79L, 65L, 73L, 74L, 75L, 82L, 83L, 
    84L, 96L, 92L, 93L, 94L, 95L, 85L, 86L, 87L, 88L, 81L, 89L, 
    90L, 91L), .Label = c("EO1 1", "EO1 10", "EO1 11", "EO1 12", 
    "EO1 13", "EO1 14", "EO1 15", "EO1 16", "EO1 2", "EO1 3", 
    "EO1 4", "EO1 5", "EO1 6", "EO1 7", "EO1 8", "EO1 9", "EO2 1", 
    "EO2 10", "EO2 11", "EO2 12", "EO2 13", "EO2 14", "EO2 15", 
    "EO2 16", "EO2 2", "EO2 3", "EO2 4", "EO2 5", "EO2 6", "EO2 7", 
    "EO2 8", "EO2 9", "EO3 1", "EO3 10", "EO3 11", "EO3 12", 
    "EO3 13", "EO3 14", "EO3 15", "EO3 16", "EO3 2", "EO3 3", 
    "EO3 4", "EO3 5", "EO3 6", "EO3 7", "EO3 8", "EO3 9", "EO4 1", 
    "EO4 10", "EO4 11", "EO4 12", "EO4 13", "EO4 14", "EO4 15", 
    "EO4 16", "EO4 2", "EO4 3", "EO4 4", "EO4 5", "EO4 6", "EO4 7", 
    "EO4 8", "EO4 9", "EO5 1", "EO5 10", "EO5 11", "EO5 12", 
    "EO5 13", "EO5 14", "EO5 15", "EO5 16", "EO5 2", "EO5 3", 
    "EO5 4", "EO5 5", "EO5 6", "EO5 7", "EO5 8", "EO5 9", "EO6 1", 
    "EO6 10", "EO6 11", "EO6 12", "EO6 13", "EO6 14", "EO6 15", 
    "EO6 16", "EO6 2", "EO6 3", "EO6 4", "EO6 5", "EO6 6", "EO6 7", 
    "EO6 8", "EO6 9"), class = "factor"), predators_trunk = c(7L, 
    10L, 9L, 15L, 18L, 11L, 5L, 7L, 15L, 12L, 6L, 12L, 7L, 13L, 
    24L, 17L, 3L, 0L, 0L, 2L, 4L, 3L, 0L, 6L, 2L, 3L, 5L, 1L, 
    5L, 12L, 18L, 15L, 7L, 0L, 5L, 1L, 17L, 7L, 13L, 19L, 7L, 
    3L, 5L, 10L, 11L, 7L, 13L, 7L, 7L, 0L, 4L, 2L, 5L, 7L, 4L, 
    7L, 8L, 7L, 9L, 20L, 13L, 2L, 12L, 7L, 0L, 7L, 2L, 2L, 2L, 
    4L, 17L, 2L, 3L, 1L, 1L, 1L, 11L, 1L, 1L, 8L, 8L, 18L, 5L, 
    6L, 6L, 5L, 6L, 5L, 9L, 2L, 8L, 13L, 13L, 5L, 3L, 5L), pH_H2O = c(4.145, 
    4.145, 4.145, 4.145, 4.1825, 4.1825, 4.1825, 4.1825, 4.1325, 
    4.1325, 4.1325, 4.1325, 4.14125, 4.14125, 4.14125, 4.14125, 
    4.265, 4.265, 4.265, 4.265, 4.21, 4.21, 4.21, 4.21, 4.18375, 
    4.18375, 4.18375, 4.18375, 4.09625, 4.09625, 4.09625, 4.09625, 
    4.1575, 4.1575, 4.1575, 4.1575, 4.1125, 4.1125, 4.1125, 4.1125, 
    4.20875, 4.20875, 4.20875, 4.20875, 3.97125, 3.97125, 3.97125, 
    3.97125, 4.025, 4.025, 4.025, 4.025, 4.005, 4.005, 4.005, 
    4.005, 4.04, 4.04, 4.04, 4.04, 4.03125, 4.03125, 4.03125, 
    4.03125, 4.4575, 4.4575, 4.4575, 4.4575, 4.52, 4.52, 4.52, 
    4.52, 4.505, 4.505, 4.505, 4.505, 4.34875, 4.34875, 4.34875, 
    4.34875, 4.305, 4.305, 4.305, 4.305, 4.32, 4.32, 4.32, 4.32, 
    4.35, 4.35, 4.35, 4.35, 4.445, 4.445, 4.445, 4.445), ant_mean_abundance = c(53.85714, 
    53.85714, 53.85714, 53.85714, 24.28571, 24.28571, 24.28571, 
    24.28571, 45.5, 45.5, 45.5, 45.5, 51.14286, 51.14286, 51.14286, 
    51.14286, 66.28571, 66.28571, 66.28571, 66.28571, 76.5, 76.5, 
    76.5, 76.5, 65.71429, 65.71429, 65.71429, 65.71429, 8.642857, 
    8.642857, 8.642857, 8.642857, 109.3571, 109.3571, 109.3571, 
    109.3571, 25.14286, 25.14286, 25.14286, 25.14286, 101.3571, 
    101.3571, 101.3571, 101.3571, 31.78571, 31.78571, 31.78571, 
    31.78571, 78.64286, 78.64286, 78.64286, 78.64286, 93.28571, 
    93.28571, 93.28571, 93.28571, 63.14286, 63.14286, 63.14286, 
    63.14286, 67.14286, 67.14286, 67.14286, 67.14286, 44.0625, 
    44.0625, 44.0625, 44.0625, 23.875, 23.875, 23.875, 23.875, 
    95.8125, 95.8125, 95.8125, 95.8125, 49.125, 49.125, 49.125, 
    49.125, 57, 57, 57, 57, 38.125, 38.125, 38.125, 38.125, 40.6875, 
    40.6875, 40.6875, 40.6875, 22, 22, 22, 22), bird_activity = c(153.24, 
    153.24, 153.24, 153.24, 153.24, 153.24, 153.24, 153.24, 0, 
    0, 0, 0, 0, 0, 0, 0, 240.96, 240.96, 240.96, 240.96, 240.96, 
    240.96, 240.96, 240.96, 0, 0, 0, 0, 0, 0, 0, 0, 154.54, 154.54, 
    154.54, 154.54, 154.54, 154.54, 154.54, 154.54, 0, 0, 0, 
    0, 0, 0, 0, 0, 107.68, 107.68, 107.68, 107.68, 107.68, 107.68, 
    107.68, 107.68, 0, 0, 0, 0, 0, 0, 0, 0, 172.42, 172.42, 172.42, 
    172.42, 172.42, 172.42, 172.42, 172.42, 0, 0, 0, 0, 0, 0, 
    0, 0, 113.8, 113.8, 113.8, 113.8, 113.8, 113.8, 113.8, 113.8, 
    0, 0, 0, 0, 0, 0, 0, 0)), .Names = c("Site", "Treatment", 
"Tree", "predators_trunk", "pH_H2O", "ant_mean_abundance", "bird_activity"
), class = "data.frame", row.names = c(NA, -96L))

А вот код, ведущий к предупреждениям:

library(lme4)
summary(glmer.nb(predators_trunk ~ scale(ant_mean_abundance) + scale(bird_activity) + scale(pH_H2O) + (1 | Site/Treatment), testdata, na.action = na.fail))
summary(glmer(predators_trunk ~ scale(ant_mean_abundance) + scale(bird_activity) + scale(pH_H2O) + (1 | Site/Treatment), testdata, family = negative.binomial(theta = 4.06643400243645), na.action = na.fail))

Мне интересно, что сводка glmer.nb не дает никаких предупреждений, но вызов glmer с использованием тэты, оцененной glmer.nb, дает мне предупреждения. Последний представляет собой вызов модели, который генерируется с использованием dredge (MuMIn) в соответствующей glmer.nb полной модели.


r lme4 mixed-models
person kdarras    schedule 17.08.2016    source источник
comment
Не могли бы вы включить данные и/или код, который предоставит нам воспроизводимый пример ? Это говорит о том, что ваши оценки стандартной ошибки могут быть менее точными, но их трудно понять без примера.   -  person Ben Bolker    schedule 17.08.2016
comment
Я сделал воспроизводимый пример сейчас. Мой минимальный набор данных был довольно большим, поэтому я хотел спросить об общем значении и беспокойстве предупреждения.   -  person kdarras    schedule 17.08.2016

Ответы (1)


arrow_upward
6
arrow_downward

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

В этом случае я сохранил два ваших соответствия, от glmer.nb и glmer, как g1 и g2. Вы можете видеть, что оценки (точечные оценки, SE, значения Z...) немного изменились, но не очень сильно, так что, по крайней мере, это должно вас успокоить.

printCoefmat(coef(summary(g1)),digits=2)
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)                  1.844      0.111    16.7   <2e-16 ***
scale(ant_mean_abundance)   -0.347      0.077    -4.5    7e-06 ***
scale(bird_activity)        -0.122      0.076    -1.6    0.107    
scale(pH_H2O)               -0.275      0.104    -2.6    0.008 ** 

> printCoefmat(coef(summary(g2)),digits=2)
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)                  1.846      0.108    17.1   <2e-16 ***
scale(ant_mean_abundance)   -0.347      0.077    -4.5    6e-06 ***
scale(bird_activity)        -0.122      0.075    -1.6    0.102    
scale(pH_H2O)               -0.275      0.102    -2.7    0.007 **

У меня есть разрабатываемая версия lme4 на Github (ветка test_mods, надеюсь, скоро будет интегрирована в основную ветку: если вы хотите установить ее, вы можете использовать devtools::install_github("lme4/lme4",ref="test_mods")), которая позволяет вам выбрать более точный (но более медленный) расчет для стандартного ошибки: это возвращает нас (почти) к тем же стандартным ошибкам, что и glmer.nb.

g3 <- update(g2, control=glmerControl(deriv.method="Richardson"))
printCoefmat(coef(summary(g3)),digits=2)
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)                  1.846      0.111    16.7   <2e-16 ***
scale(ant_mean_abundance)   -0.347      0.077    -4.5    6e-06 ***
scale(bird_activity)        -0.122      0.076    -1.6    0.106    
scale(pH_H2O)               -0.275      0.104    -2.6    0.008 ** 

all.equal(coef(summary(g1))[,"Std. Error"],
          coef(summary(g3))[,"Std. Error"])
[1] "Mean relative difference: 0.001597978"

Пакет glmmTMB (на Github) тоже дает почти такие же результаты:

library(glmmTMB)
g5 <- glmmTMB(predators_trunk ~ scale(ant_mean_abundance) +
                     scale(bird_activity) + scale(pH_H2O) +
                     (1 | Site/Treatment), testdata,
              family=nbinom2)
printCoefmat(coef(summary(g5))[["cond"]],digits=2)
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)                  1.852      0.110    16.8   <2e-16 ***
scale(ant_mean_abundance)   -0.348      0.077    -4.5    7e-06 ***
scale(bird_activity)        -0.123      0.076    -1.6    0.106    
scale(pH_H2O)               -0.276      0.105    -2.6    0.008 **
person Ben Bolker    schedule 17.08.2016
comment
Я получаю то же предупреждающее сообщение, что и OP, и пытаюсь реализовать ответ @Ben. Я установил из git, используя devtools::install_github("lme4/lme4",ref="test_mods"), а также (в случае, если ветка test_mods была объединена с основной веткой) с помощью install.packages('lme4'). Независимо от метода установки я получаю следующую ошибку при вызове choose = glmerControl(deriv.method=c("Richardson")): Error in glmerControl(deriv.method = c("Richardson")) : unused argument (deriv.method = c("Richardson")). Любые указатели? - person Cole Robertson; 05.12.2018