计算的R中使用confint固定效应的CI

Question

我想进行引导，以获得我的固定效应的95％CI在二项式GLMM：

m <- glmer(cbind(df$Valid.detections, df$Missed.detections) ~ distance + 
       Habitat + Replicate + transmitter.depth + receiver.depth + 
       wind.speed + wtc + Transmitter + (1 | Unit) + 
       (1 | SUR.ID) + distance:Transmitter + 
       distance:Habitat + distance:transmitter.depth + distance:receiver.depth + 
       distance:wind.speed, data = df, family = binomial(link=logit),control=glmerControl(calc.derivs=F)) 

我发现confint()功能是能够做到这一点，所以我指定了功能：

confint(m, method = "boot", boot.type = "basic", seed = 123, nsim = 1000) 

该函数在我决定终止之前已运行了8个多小时。协议终止后，我得到了返回下面的警告消息（X10）：

Warning messages: 
1: In (function (fn, par, lower = rep.int(-Inf, n), upper = rep.int(Inf, : 
    failure to converge in 10000 evaluations 

我的问题是：1）我担心这些警告消息？如果是这样，我该如何处理它们？2）因为8小时后它仍在运行我不知道执行此功能需要多长时间。因此，在执行此功能时有一些进度条会很好。我读了confint（）可以从bootMer接受参数，所以我包括参数.progress =“TXT”，导致：

confint(m, method = "boot", boot.type = "basic", seed = 123, nsim = 1000, .progress = "txt") 

，但它似乎并没有工作。或者，如果有更好的方法来实现相同的目标，我愿意接受建议。

感谢所有帮助

Answer 1

下面是一个例子：

library("lme4") 
(t1 <- system.time(
    gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd), 
       data = cbpp, family = binomial))) 
## user system elapsed 
## 0.188 0.000 0.186 

nranpars <- length(getME(gm1,"theta")) 
nfixpars <- length(fixef(gm1)) 

(t2 <- system.time(c1 <- confint(gm1,method="boot", nsim=1000, 
        parm=(nranpars+1):(nranpars+nfixpars), 
        .progress="txt"))) 

## user system elapsed 
## 221.958 0.164 222.187 

注意，这个进度条长仅80个字符，所以它只是每个1000至80年= 12的自举迭代后递增。如果你的原始模型花了一个小时，以适应那么你不应该指望能看到进度条的活动直到12小时后...

(t3 <- system.time(c2 <- confint(gm1, 
        parm=(nranpars+1):(nranpars+nfixpars)))) 

## user system elapsed 
## 5.212 0.012 5.236 

1000引导销售代表可能是矫枉过正 - 如果你的模型拟合很慢，您可能可以从200个自举代表获得合理的配置项。

我试过这个也是optimizer="nloptwrap"，希望它能加快速度。它确实有，但有一个缺点（见下文）。

(t4 <- system.time(
    gm1B <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd), 
       data = cbpp, family = binomial, 
       control=glmerControl(optimizer="nloptwrap")))) 
## user system elapsed 
## 0.064 0.008 0.075 

(t5 <- system.time(c3 <- confint(gm1B,method="boot", nsim=1000, 
        parm=(nranpars+1):(nranpars+nfixpars), 
        .progress="txt",PBargs=list(style=3)))) 
## 
## user system elapsed 
## 65.264 2.160 67.504 

这是更快，但（在这种情况下37）给出警告有关 收敛。根据all.equal()，以这种方式计算的置信区间仅有约2％的差异。 （在包装本身中仍然存在一些皱纹...）

加速这一过程的最佳方法是并行化 - 不幸的是，这种方式会失去使用进度条的能力......

(t6 <- system.time(c4 <- confint(gm1,method="boot", nsim=1000, 
        parm=(nranpars+1):(nranpars+nfixpars), 
        parallel="multicore", ncpus=4))) 

## 
##  user system elapsed 
## 310.355 0.916 116.917 

这需要更多的用户时间（其依靠的所有内核使用的时间），但过去的时间减少了一半。 （用4个内核做得更好，但是两倍的速度仍然不错，它们是虚拟Linux机器上的虚拟内核，真实（非虚拟）内核可能会有更好的性能。）

With nloptwrap and multicore结合我可以把时间缩短到91秒（用户）/ 36秒（经过）。