首页
学习
活动
专区
圈层
工具
发布
社区首页 >问答首页 >GLMER警告:方差协方差矩阵[.]不为正定或包含NA值。

GLMER警告:方差协方差矩阵[.]不为正定或包含NA值。
EN

Stack Overflow用户
提问于 2016-08-17 13:00:01
回答 1查看 10.5K关注 0票数 9

有时,我发现来自glmer的GLMMs包lme4显示了以下警告消息,当它们的摘要被调用时:

代码语言:javascript
复制
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维基也没有详细说明这一点。在问题中,在解决这类错误消息之前已经解决了这个问题,这里讨论的重点是一个特定的模型,而不是警告消息的含义。

所以问题是:我应该担心这条信息吗,还是因为它只是一个警告而不是一个错误,正如它所说的,它是“回到了RX估计的var-cov”(不管RX是什么)。

有趣的是,虽然总结指出模型无法收敛,但我没有得到通常的红色收敛警告。

下面是一个(最小的)数据集:

代码语言:javascript
复制
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))

下面是导致警告的代码:

代码语言:javascript
复制
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.nb估计的θ对glmer的调用确实给了我警告。后者是使用dredge (MuMIn)在相应的glmer.nb完整模型上生成的模型调用。

EN

回答 1

Stack Overflow用户

回答已采纳

发布于 2016-08-17 21:58:04

此警告表示您的标准错误估计可能不太准确。但是和所有的警告一样,很难确定,最好的办法就是尽可能的反复检查。

在本例中,我从glmer.nbglmer中将您的两次匹配保存为g1g2。你可以看到估计(点估计,SE,Z值.)已经改变了一点点,但不是很多,所以至少这应该能让你安心。

代码语言:javascript
复制
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的开发版本( test_mods分支,希望很快集成到主分支中:如果您想安装它,可以使用devtools::install_github("lme4/lme4",ref="test_mods")),它允许您为标准错误选择更精确(但更慢)的计算:这使我们回到(几乎)与glmer.nb相同的标准错误。

代码语言:javascript
复制
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)也给出了几乎相同的结果:

代码语言:javascript
复制
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 ** 
票数 10
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/38997371

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档