[R-sig-ME] mixed model negative bionomial

Paul Johnson pauljohn32 at gmail.com
Fri Apr 5 08:09:00 CEST 2013


Hi

I refitted with his data, making the change in the formula suggested.

It appears "underdispersed" rather than over, if the relevant
benchmark is 120, you expect deviance to be about that same value.
Isn't that how you would read it?

> sick2.glmm <- glmer(cbind(NSick,(Ntest - NSick)) ~ Breed + (1 | Farm), data = dat,
+ family = binomial)
> summary(sick2.glmm)
Generalized linear mixed model fit by the Laplace approximation
Formula: cbind(NSick, (Ntest - NSick)) ~ Breed + (1 | Farm)
   Data: dat
   AIC   BIC logLik deviance
 88.99 100.1  -40.5    80.99
Random effects:
 Groups Name        Variance Std.Dev.
 Farm   (Intercept)  0        0
Number of obs: 120, groups: Farm, 4

Fixed effects:
               Estimate Std. Error z value Pr(>|z|)
(Intercept)     -3.9404     0.2524 -15.610   <2e-16 ***
BreedPerindale  -0.7601     0.5153  -1.475   0.1402
BreedRomney     -0.9463     0.4813  -1.966   0.0493 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
            (Intr) BrdPrn
BreedPerndl -0.490
BreedRomney -0.524  0.257


And if one of these is underdispersed, is one supposed to change something?

Ever see this article? It will make you re-think. Maybe even
second-guess things you thought you knew.

Skrondal, A. and Rabe-Hesketh, S. (2007). Redundant overdispersion
parameters in multilevel models for categorical responses. Journal of
Educational and Behavioral Statistics 32, 419-430.

which you can retrieve from the GLLAAM website
(http://www.gllamm.org/sophia.html).

On Fri, Mar 29, 2013 at 5:33 PM, Philippi, Tom <tom_philippi at nps.gov> wrote:
> Andy--
> Also, check your glmer() formula.  I believe that for family=binomial, the
> response is specified as cbind(success,failure) and not
> cbind(success,trials) as in SAS.  Your example code
> glmer(cbind(dat$NSick,dat$Ntest)) might need to be cbind(NSick,NWell) or
> cbind(NSick,Ntest-NSick).  I don't suggest that this will eliminate all of
> your overdispersion, but if your Ntest is really the number of individuals
> at risk in each group, your coding will inflate the number of
> non-responses.
>
> Tom
>
>
> On Wed, Mar 27, 2013 at 6:58 PM, Ben Bolker <bbolker at gmail.com> wrote:
>
>> Andrew McFadden (Andy <Andrew.McFadden at ...> writes:
>>
>> >
>> > Hi all
>>
>> > Really appreciate a hand here. I am trying to model some data with a
>> > bionomial outcome. I believe that I need to use a negative bionomial
>> > distribution as there were a lot of samples where a large number of
>> > zeros were present i.e. none sick and when modelled using a
>> > bionomial distribution in the lme4 package the residuals were
>> > extremely high. Hence the attempted use of gamlss package.
>>
>>   I can't help with gamlss at the moment, but: the negative binomial
>> is *not* an appropriate generalization of the binomial (unless you
>> have low probabilities everywhere and want to approximate the binomial
>> by a Poisson, in which case you would then get a NB).  Beta-binomial
>> is to binomial and NB is to Poisson.  You can model overdispersion
>> (which is *one* route to many zeros -- another is simply very low
>> overall prevalence, and a third is zero-inflation) in various ways:
>> see http://glmm.wikidot.com/faq ...
>>
>> > I have had difficulty coding the model for the gamlss package,
>> > perhaps I have done something wrong. Also I would like to include
>> > the denominator in the outcome as the sample size varied per group
>> > i.e. in the lme4 package I coded it as:
>> > glmer(cbind(dat$NSick,dat$Ntest); but couldn't seem to do this in
>> > gamlss.
>>
>> > The data below is "made up data", but reflects the analysis I am
>> > doing i.e. data clustered within farm (hence the need for a random
>> > effect), Lots of zero outcome (and the need to use a negative
>> > bionomial equation).
>>
>>  Overdispersion can be modeled in various packages; glmmADMB handles
>> zero-inflation (although it's not well tested for binomial models,
>> so check your results especially carefully).
>>
>>   Ben Bolker
>>
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>
>
>
> --
> -------------------------------------------
> Tom Philippi, Ph.D.
> Quantitative Ecologist & Data Therapist
> Inventory and Monitoring Program
> National Park Service
> (619) 523-4576
> Tom_Philippi at nps.gov
> http://science.nature.nps.gov/im/monitor
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models



-- 
Paul E. Johnson
Professor, Political Science      Assoc. Director
1541 Lilac Lane, Room 504      Center for Research Methods
University of Kansas                 University of Kansas
http://pj.freefaculty.org               http://quant.ku.edu



More information about the R-sig-mixed-models mailing list