[R] How to make try to catch warnings in logistic glm

Fredrik Nilsson laf.nilsson at gmail.com
Mon Jun 22 14:55:32 CEST 2009


Dear list,

>From an earlier post I got the impression that one could promote
warnings from a glm to errors (presumably by putting
options(warn=1)?), then try() would flag them as errors. I’ve spent
half the day trying to do this, but no luck. Do you have an explicit
solution?

My problems is that I am trying to figure out during what conditions
one may find 5 significant parameters in a logistic regression with
only 21 bad and 28 outcomes (this appears in a published article) and
I would like to have  my.glm<-try(glm(outcome~ald + fat+ hstop + btun
+ time, data=DF, family=binomial)) to indicate if warnings
(convergence / probabilities equal to 0/1) occurs (see attempt to code
below if my explanation is terse).

Best regards,

Fredrik Nilsson

------------------------------------------------------------------------------
Fredrik Nilsson, PhD
Competence Centre for Clinical Research
Lund University Hospital


Warnings I'd like to catch:
Warning in glm.fit(x = X, y = Y, weights = weights, start = start,
etastart = etastart,  :
 algorithm did not converge
Warning in glm.fit(x = X, y = Y, weights = weights, start = start,
etastart = etastart,  :
 fitted probabilities numerically 0 or 1 occurred

script
ngood<-28
nbad<-21
resmat<-"DFr"
j<-0
nsim<-100
res<-matrix(NA, nrow=nsim, ncol=6)
ow <- options("warn")
options(warn= 1)
outcome<-rep(c("G","B"), times=c(ngood, nbad))
outcome<-as.factor(outcome)
agood<-rep(c(0,1), times=c(27,1))
abad<-rep(c(0,1), times=c(14,7))
thgood<-rep(c(0,1), times=c(24,4))
thbad<-rep(c(0,1), times=c(12,9))

for (i in 1:nsim)
{
 tgood<-sample(agood, ngood)
 tbad<-sample(abad, nbad)
 hstop<-c(tgood, tbad)
 hstop<-as.factor(hstop)
 aldgood<-rnorm(ngood, mean=54, sd=8/0.675)
 aldbad<-rnorm(nbad, mean=64, sd=8/0.675)
 ald<-c(aldgood, aldbad)
 fatgood<-exp(rnorm(ngood, mean=log(23.9)-0.06^2/2,sd=0.06))
 fatbad<-exp(rnorm(nbad, mean=log(27.4)-0.09^2/2, sd=0.09))
 fat<-c(fatgood, fatbad)
 thgood<-sample(thgood, ngood)
 thbad<-sample(thbad, nbad)
 btun<-c(thgood, thbad)
 btun<-as.factor(btun)
 timegood<-exp(rnorm(ngood, mean=log(443)-0.5^2/2, sd=0.5))
 timebad<- exp(rnorm(nbad, mean=log(555)-0.6^2/2, sd=0.6))
 time<-c(timegood, timebad)

 DF<-data.frame(outcome, ald, fat, hstop, btun, time)
 my.glm<-try(glm(outcome~ald + fat+ hstop + btun + time, data=DF,
family=binomial))
 test<-!is.null(attr(my.glm, "try-error"))
 if(test | !my.glm$converged)
 {
   res[i, 1:6]<-NA
 } else
 {
   W<-summary(my.glm)$coefficients
   lres<-as.numeric(W[,4]<=0.05)
   if (sum(lres)==5)
   {
     j<-j+1
     assign(paste(resmat, j, sep=".", collapse=""), DF)
   }
   res[i,1:6]<-lres
 }
}

options(ow)



From: Dieter Menne <dieter.menne_at_menne-biomed.de>
Date: Fri, 14 Dec 2007 20:49:20 +0000 (UTC) Caroline Paulsen <cpaulsen
<at> u.washington.edu> writes:
>
> I'm attempting to run 250 permutations of a negative binomial GLM
> model for data on fish counts. Many of the models are fit
> appropriately, but some issue warnings such as "convergence not
> reached" or "step size truncated due to divergence."
....
> I'd like to figure out a way to flag the models that have warnings and
> output them into either a table or into R console. Then I could
> evaluate the problems associated with these models and know not choose
> them as the best models for fitting the data.
You could promote the options to errors (?warning) and use try() Dieter




More information about the R-help mailing list