Post-glmer 的临时测试

Question

我正在使用广义线性混合模型（glmer、lme4-package）用 R 分析我的二项式数据集。我想使用 Tukey 的 post-hoc 测试（glht，multcomp-package）对某个固定效应（"Sound"）进行成对比较。

大部分工作正常，但我的一个固定效应变量 ("SoundC") 完全没有变化（96 次是“1”，零次是“0”），看起来Tukey 的测试无法解决这个问题。所有与此 "SoundC" 的成对比较给出的 p 值为 1.000，而有些则明显显着。

作为验证，我将 96 个“1”中的一个更改为“0”，之后我再次获得了正常的 p 值和我预期的显着差异，而差异实际上在之后变小了我的手动更改。

有人有解决办法吗？如果不是，是否可以使用我修改后的数据集的结果并报告我的手动更改？

可重现的例子：

Response <- c(1,0,1,1,0,1,1,0,1,1,0,1,1,0,1,1,1,1,0,
              0,1,1,0,1,1,0,1,1,0,1,1,0,1,1,0,1,1,0,1,1,0,
              1,1,0,1,1,0,1,1,1,1,0,0,1,1,0,1,1,0,1,1,0,1,1,0,1)    
Data <- data.frame(Sound=rep(paste0('Sound',c('A','B','C')),22),
                   Response,
                   Individual=rep(rep(c('A','B'),2),rep(c(18,15),2)))


# Visual
boxplot(Response ~ Sound,Data)

# Mixed model
library (lme4)
model10 <- glmer(Response~Sound + (1|Individual), Data, family=binomial)

# Post-hoc analysis
library (multcomp)
summary(glht(model10, mcp(Sound="Tukey")))

Answer 1

这接近 CrossValidated question; you are definitely seeing complete separation, where there is a perfect division of your response into 0 vs 1 results. This leads to (1) infinite values of the parameters (they're only listed as non-infinite due to computational imperfections) and (2) crazy/useless values of the Wald standard errors and corresponding $p$ values (which is what you're seeing here). Discussion and solutions are given here, here, and here，但我将在下面说明更多。

暂时成为统计上的抱怨者：无论如何，您真的不应该尝试只用 3 个水平来拟合随机效应（参见 http://glmm.wikidot.com/faq）...

Firth 校正逻辑回归：

library(logistf)
L1 <- logistf(Response~Sound*Individual,data=Data,
        contrasts.arg=list(Sound="contr.treatment",
         Individual="contr.sum"))

                                 coef se(coef)            p
(Intercept)              3.218876e+00 1.501111 2.051613e-04 
SoundSoundB             -4.653960e+00 1.670282 1.736123e-05 
SoundSoundC             -1.753527e-15 2.122891 1.000000e+00 
IndividualB             -1.995100e+00 1.680103 1.516838e-01 
SoundSoundB:IndividualB  3.856625e-01 2.379919 8.657348e-01 
SoundSoundC:IndividualB  1.820747e+00 2.716770 4.824847e-01

标准误差和 p 值现在是合理的（A 与 C 比较的 p 值为 1，因为实际上没有区别...）

具有弱先验的混合贝叶斯模型：

library(blme)
model20 <- bglmer(Response~Sound + (1|Individual), Data, family=binomial,
                  fixef.prior = normal(cov = diag(9,3)))

##              Estimate Std. Error    z value     Pr(>|z|)
## (Intercept)  1.711485   2.233667  0.7662221 4.435441e-01
## SoundSoundB -5.088002   1.248969 -4.0737620 4.625976e-05
## SoundSoundC  2.453988   1.701674  1.4421024 1.492735e-01

固定效应方差-协方差矩阵的规范diag(9,3)产生

$$ \剩下（ \begin{数组}{ccc} 9 & 0 & 0 \ 0 & 9 & 0 \ 0 & 0 & 9 \end{数组} \正确的） $$

换句话说，3指定矩阵的维度（等于固定效应参数的数量），9指定方差——这对应于3的标准差或95%的范围大约 $\pm 6$，这对于逻辑尺度响应来说相当 large/weak/uninformative。

这些大致一致（模型非常不同）

library(multcomp)
summary(glht(model20, mcp(Sound="Tukey")))

##                     Estimate Std. Error z value Pr(>|z|)    
## SoundB - SoundA == 0   -5.088      1.249  -4.074 0.000124 ***
## SoundC - SoundA == 0    2.454      1.702   1.442 0.309216    
## SoundC - SoundB == 0    7.542      1.997   3.776 0.000397 ***

正如我上面所说，无论如何，在这种情况下我不会推荐混合模型...

Post-glmer 的临时测试

Post-hoc test for glmer

r

lme4