Erin McMullen Jonaitis <jonaitis <at> wisc.edu> writes:
> Hello, all, I am running some simulations to estimate power for a
> complicated epidemiological study, and am using lme and lmer to get
> these estimates. I have to run a few thousand iterations, and once
> in a great while, an iteration will create fake data such that the
> model won't converge. I see from Google searches that this is not
> an uncommon situation.
> My question: is there a way to extract the convergence value from an
> lme or lmer model? It prints on the screen, but how can I get hold
> of it to evaluate it with some function? What I'd like to do is
> build a failsafe into my program so that if a particular model in an
> iteration doesn't converge, I call a redo on that iteration. This
> way the program will keep running and not stop in a fit of pique in
> the middle of my long simulation.
> If I can't do this, my fallback will be to try setting lmeControl
> options such that even "bad" models return parameter estimates etc
> -- once or twice in 10,000 iterations should not ruin things too
> badly -- but I'd like to try it the cleaner way first.
There's a somewhat hack-ish solution, which is to use options(warn=2)
to 'upgrade' warnings to errors, and then use try() or tryCatch() to
catch them.
More fancily, I used code that looked something like this to save
warnings as I went along (sorry about the <<- ) in a recent simulation
study. You could also check w$message to do different things in the
case of different warnings.
withCallingHandlers(tryCatch(fun(n=nvec[j],tau=tauvec[i],...),
error = function(e) {
warn[k,i,j] <<- paste("ERROR:",e$message)
NA_ans}),
warning = function(w) {
warn[k,i,j] <<- w$message
invokeRestart("muffleWarning")
})
In the slightly longer run, we are working on getting the development
(lme4Eigen) version of lme4 to save convergence warnings in an
accessible slot. I don't know if I would hold my breath for this
to be back-ported to nlme, though ...
Some of these discussions might be better suited for r-sig-mixed-models
<at> r-project.org
Ben Bolker
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
