[R-sig-ME] "bootstrap runs failed" in bootMer within confint.merMod
Paul Johnson
paul.johnson at glasgow.ac.uk
Tue Sep 9 10:46:29 CEST 2014
An update on this question: it's now an issue on github, at https://github.com/lme4/lme4/issues/231
Thanks to Ben Bolker for following this up - see the comments for updates on progress.
Paul Johnson
On 14 Aug 2014, at 14:15, Paul Johnson <Paul.Johnson at glasgow.ac.uk> wrote:
> Hi,
>
> I’ve been using confint.merMod to get 95% CIs by parametric bootstrapping, and have been getting the warnings of the type "In bootMer(object, bootFun, nsim = nsim, ...) : some bootstrap runs failed (30/100)”. (I realise than 100 samples are too few - this is just an example). The function still returns CIs using the bootstrap samples that didn’t fail, but I don’t think these CIs are valid as the failures are very unlikely to be a random sample. Here’s an example using the grouseticks data set that comes with lme4 (see ?grouseticks). I chose this data set because it’s a classic example data set so presumably isn’t pathological, but I’ve been encountering the same problem with other apparently healthy glmer fits.
>
>> library(lme4)
>> full_mod1 <- glmer(TICKS ~ YEAR+scale(HEIGHT)+(1|BROOD)+(1|INDEX)+(1|LOCATION), family="poisson", data=grouseticks)
>> set.seed(12345)
>> confint(full_mod1, method = "boot", boot.type = "basic", nsim = 100, parallel = "multicore", ncpus = 8)
> Computing bootstrap confidence intervals ...
> 2.5 % 97.5 %
> sd_(Intercept)|INDEX 0.4624700492 0.6550336
> sd_(Intercept)|BROOD 0.5740796527 1.0482547
> sd_(Intercept)|LOCATION 0.2497343230 1.0574176
> (Intercept) -0.0002611865 0.7208569
> YEAR96 0.7770837606 1.5857202
> YEAR97 -1.4648325648 -0.4484355
> scale(HEIGHT) -1.1335522182 -0.5879858
> Warning message:
> In bootMer(object, bootFun, nsim = nsim, ...) :
> some bootstrap runs failed (30/100)
>
> (NB I scaled height to prevent a warning but it makes no difference to the likelihood.)
>
> I tried to replicate the problem using bootMer directly to produce fixed effect estimates from 10 samples:
>
>> bootMer(full_mod1, fixef, nsim = 10)$t
> (Intercept) YEAR96 YEAR97 scale(HEIGHT)
> [1,] 0.6387534 0.7467317 -1.1291369 -0.8741711
> [2,] NA NA NA NA
> [3,] 0.1920928 1.3666922 -0.5821328 -0.6742498
> [4,] NA NA NA NA
> [5,] NA NA NA NA
> [6,] 0.5312426 0.9639612 -1.1041643 -0.7550492
> [7,] NA NA NA NA
> [8,] 0.2328530 1.1016917 -0.8010541 -0.8480263
> [9,] 0.2685594 1.2806690 -0.7452719 -0.8302675
> [10,] 0.1726928 1.1423894 -0.8989542 -1.0344128
> Warning message:
> In bootMer(full_mod1, fixef, nsim = 10) : some bootstrap runs failed (4/10)
>
> …and get the same problem, not surprisingly. Here’s a manual version of what (I believe) bootMer is doing:
>
>
>> simTICKS.tab <- simulate(full_mod1, nsim=10)
>> t(apply(simTICKS.tab, 2, function(simTICKS) fixef(glmer(simTICKS ~ YEAR+scale(HEIGHT)+(1|BROOD)+(1|INDEX)+(1|LOCATION), family="poisson", data=grouseticks))))
> (Intercept) YEAR96 YEAR97 scale(HEIGHT)
> sim_1 0.1793276 1.2265976 -0.5753168 -0.6435871
> sim_2 0.3918824 1.4583529 -0.8958396 -0.8513233
> sim_3 0.2565916 1.1230342 -0.9369728 -1.0817019
> sim_4 0.3507285 1.1865904 -1.0654895 -0.9193403
> sim_5 0.5688213 0.7291811 -0.9734654 -0.8575988
> sim_6 0.6160919 1.0552260 -1.5139770 -1.0050926
> sim_7 0.1594954 1.1342809 -1.0624134 -0.8107930
> sim_8 0.7956874 0.8290606 -1.2491060 -1.0024299
> sim_9 0.3081545 1.2677562 -1.0398100 -0.9458849
> sim_10 0.4325345 0.9390844 -0.9241206 -0.7341589
>
> …no errors, not even a warning. The discrepancy between bootMer and the DIY version above appears to arise from the bootMer fitting function, refit, being less tolerant of warning signs (non-convergence?) than glmer, which can be seen by putting refit into the DIY bootstrap:
>
>> t(apply(simTICKS.tab, 2, function(simTICKS) fixef(refit(full_mod1, simTICKS))))
> Warning messages:
> 1: In pwrssUpdate(pp, resp, tolPwrss, GQmat, compDev, fac, verbose) :
> Cholmod warning 'not positive definite' at file:../Cholesky/t_cholmod_rowfac.c, line 431
> 2: In pwrssUpdate(pp, resp, tolPwrss, GQmat, compDev, fac, verbose) :
> Cholmod warning 'not positive definite' at file:../Cholesky/t_cholmod_rowfac.c, line 431
> Error in t(apply(simTICKS.tab, 2, function(simTICKS) fixef(refit(full_mod1, :
> error in evaluating the argument 'x' in selecting a method for function 't': Error in t(apply(simTICKS.tab, 2, function(simTICKS) fixef(refit(full_mod1, :
> (maxstephalfit) PIRLS step-halvings failed to reduce deviance in pwrssUpdate
>
> My questions are:
>
> Should I trust the error- and warning-free DIY bootstrap results? I would say “yes”, on the (slightly vague) basis that glmer produced no warnings, and that the original model fits well to what is a pretty information-rich data set. So is bootMer being oversensitive? I’m aware that there have been some teething problems with the optimisers and convergence checks in lme4 1.0+ — are these ongoing? I’m writing a tutorial document in which I’d like to demonstrate parametric bootstrapping for CIs. I could use the DIY approach, but I’d much rather use confint, so it would be great if there was a simple fix to this problem.
>
> Thanks in advance,
> Paul Johnson
> (the Glasgow one, not the Kansas or Oxford ones)
>
>> sessionInfo()
> R version 3.1.1 (2014-07-10)
> Platform: x86_64-apple-darwin13.1.0 (64-bit)
>
> locale:
> [1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
>
> attached base packages:
> [1] parallel stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] lme4_1.1-7 Rcpp_0.11.2 Matrix_1.1-4
>
> loaded via a namespace (and not attached):
> [1] boot_1.3-11 grid_3.1.1 lattice_0.20-29 MASS_7.3-33 minqa_1.2.3 nlme_3.1-117 nloptr_1.0.0 splines_3.1.1 tools_3.1.1
>
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
More information about the R-sig-mixed-models
mailing list