[R] glm binomial with no successes

[Gavin Simpson]

[Following up for my own personal education so a bit OT!]

Naively, I would have thought that package multcomp would be of use
here. So I tried, for my own comprehension and education, to answer the
OP's question using multcomp. Here's what I got:

## make this reproducible (I hope)
set.seed(1234)
s <- c(rpois(8, 4), rep(0, 4))
f <- rpois(12, 30)
tr <- gl(3, 4)
sf <- cbind(s,f)

## fit the glm
mod <- glm(sf ~ tr, family=binomial)
summary(mod) ## tr2 and tr3 not different from reference level tr1
anova(mod, test = "Chisq") ## tr is signif

## multiple comparison of levels of tr
require(multcomp)
mod.glht <- glht(mod, linfct = mcp(tr = "Tukey"))
mod.glht
summary(mod.glht)

If I interpret this correctly, both summary(mod) and summary(mod.glht)
suggest that there are no significant differences between the 3 levels
of tr, but that tr, as a whole, is better than the Null model (as shown
by anova(mod) )?

Is my interpretation correct, for this specific example, or am I abusing
multcomp and statistics in this case?

Thanks for your time and indulgence of a more statistically-related than
R-related question

All the best,

G

On Wed, 2008-02-27 at 12:51 +0100, juli pausas wrote: 
> Thank you very much for your reply.
> Then I understand that would not be correct to perform the test in
> summary for testing the significance of the different levels of a
> factor in relation to the first level, including when there are more
> than 2 levels, as in my real case; at least for binomial regressions.
> So here a more close-to-real example, with a 3-level factor
> 
> s <- c(rpois(8, 4), rep(0, 4))
> f <- rpois(12, 30)
> tr <- gl(3, 4)
> sf <- cbind(s,f)
> drop1(glm(sf ~ tr, family="binomial"), test="Chisq") # significant
> summary(glm(sf ~ tr, family="binomial"))             # the 3rd level
> is not significant from the 1st
> 
> So I understand that I need to explite the data and perform the two
> tests separately:
> 
> drop1(glm(sf ~ tr, family="binomial", subset=(tr %in% c("1", "2"))),
> test="Chisq") # ns as expected
> 
> drop1(glm(sf ~ tr, family="binomial", subset=(tr %in% c("1", "3"))),
> test="Chisq") # significant, as expected
> 
> Is this the correct approach?
> Many thanks
> 
> Juli
> 
> On Wed, Feb 27, 2008 at 12:13 PM, Prof Brian Ripley
> <ripley at stats.ox.ac.uk> wrote:
> > On Wed, 27 Feb 2008, juli pausas wrote:
> >
> >  > Dear all,
> >  > I have a question on glm, family binomial. I do not see significant
> >  > differences between the levels of a factor (treatment) if all data for
> >  > a level is 0; and replacing a 0 for a 1 (in fact reducing the
> >  > difference), then I detect the significant difference that I expected.
> >
> >  This is because you are using the wrong test, one with negligible power.
> >  See MASS4 pp.197-8 -- you need to use the LRT, as in
> >
> >  > drop1(glm(sf ~ tr, family="binomial"), test="Chisq")
> >  Single term deletions
> >
> >  Model:
> >  sf ~ tr
> >         Df Deviance    AIC    LRT   Pr(Chi)
> >  <none>       1.595 17.730
> >  tr      1   24.244 38.379 22.649 1.944e-06
> >
> >  (and in your example you can replace 'drop1' by 'anova').
> >
> >
> >  > Is there a way to overcome this problem? or this is an expected
> >  > behaviour ?  Here is an example:
> >  >
> >  > s <- c(2,4,4,5,0,0,0,0)
> >  > f <- c(31,28,28,28,32,37,34,35)
> >  > tr <- gl(2, 4)
> >  > sf <- cbind(s,f)  # numbers of successes and failures
> >  > summary(glm(sf ~ tr, family="binomial"))  # tr ns
> >  >
> >  > sf[8,1] <- 1
> >  > summary(glm(sf ~ tr, family="binomial"))  # tr significative **
> >  >
> >  > Thanks for any suggestion
> >  >
> >  > Juli
> >  >
> >  > --
> >  > http://www.ceam.es/pausas
> >  >
> >  > ______________________________________________
> >  > R-help at 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.
> >  >
> >
> >  --
> >  Brian D. Ripley,                  ripley at stats.ox.ac.uk
> >  Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
> >  University of Oxford,             Tel:  +44 1865 272861 (self)
> >  1 South Parks Road,                     +44 1865 272866 (PA)
> >  Oxford OX1 3TG, UK                Fax:  +44 1865 272595
> >
> 
> 
> 
-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson             [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,          [f] +44 (0)20 7679 0565
 Pearson Building,             [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London          [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT.                 [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%