[R-sig-ME] Odd ANOVA degrees of freedom with ZI component of glmmTMB model

Ben Bolker bbo|ker @end|ng |rom gm@||@com
Thu Jan 12 00:27:00 CET 2023


   It would be great if you can send along a reproducible example, 
privately if necessary.  It's always alarming to have undiagnosed 
weirdness happening, even if you can make it go away by limiting the 
problem ...

On 2023-01-11 6:08 p.m., Elliot Johnston wrote:
> Thanks for getting back to me Ben and John. As I was making a 
> reproducible example for this thread, I dropped all of the dataframe 
> columns not used in the analysis and ended up trimming the dataframe 
> down from 35 columns to 6. The zero-inflated ANOVA output now appears 
> more sensical:
> 
>  > car::Anova(m1, component = "zi")
> 
> Analysis of Deviance Table (Type II Wald chisquare tests)
> 
> Response: Count
> 
>                                   Chisq   Df  Pr(>Chisq)
> Time_Period            1.9271  2    0.38154
> Assignment             1.5043  1    0.22001
> Time_Period:Assignment 7.9605  2    0.01868 *
> 
> The statistics for the interaction are the same as before, but the Time 
> Period and Assignment terms now make more sense. No observations were 
> dropped and the model specification remained the same, I just dropped a 
> number of columns. Seems strange, but I have the desired output so I 
> don't feel the need to troubleshoot the source of the error at this 
> point. Hope that's alright with you both. Thanks again for getting back 
> to me.
> 
> -Elliot
> 
> On Wed, Jan 11, 2023 at 2:52 PM John Fox <jfox using mcmaster.ca 
> <mailto:jfox using mcmaster.ca>> wrote:
> 
>     Dear Elliot and Ben,
> 
>     Yes, something definitely seems wrong here, and as usual a reproducible
>     example would help. Given that, and as soon as I have some time, I'll
>     try to see what went wrong, but I won't be able to do that this week.
> 
>     My apologies for the problem,
>        John
> 
>     John Fox, Professor Emeritus
>     McMaster University
>     Hamilton, Ontario, Canada
>     web: https://socialsciences.mcmaster.ca/jfox/
>     <https://socialsciences.mcmaster.ca/jfox/>
> 
>     On 2023-01-11 1:12 p.m., Ben Bolker wrote:
>      >       The difference in 'Df' between the two components, which
>     appear to
>      > have the same fixed-effect model specification, is definitely
>     surprising.
>      >
>      >       It's not surprising that chisq=9.46 with 2 df could have a
>     lower
>      > p-value than chisq=9.89 with 3 df; the larger the df (i.e. the
>     larger
>      > the difference in the number of parameters/complexity between the
>     two
>      > models implicitly being compared), the more dispersed the null
>      > distribution of the deviance difference (=='chisq').
>      >
>      >     To troubleshoot I would look at the guts of
>     glmmTMB:::Anova.glmmTMB
>      > and see what's going on. I'm not claiming that will be obvious:
>     if you
>      > can post a *reproducible* example to the glmmTMB issues list I'd be
>      > happy to take a look.
>      >
>      > On 2023-01-11 12:14 PM, Elliot Johnston wrote:
>      >> Hi all,
>      >>
>      >> I am using the car package to run ANOVAs (type II Wald chi square
>      >> tests) on
>      >> the following model:
>      >>
>      >> m1 <- glmmTMB(Count ~ Time_Period*Assignment + (1|Region/Site_ID),
>      >>                    ziformula = ~ Time_Period*Assignment +
>      >> (1|Region/Site_ID),
>      >>                    data = allbirds, family = poisson)
>      >>
>      >> Time Period has three levels and Assignment has two levels. When
>     running
>      >> the ANOVA on the conditional component -- car::Anova(m1, component =
>      >> "cond") -- the degrees of freedom in the output is as I would
>     expected
>      >> (n-1):
>      >>
>      >>                                  Chisq    Df  Pr(>Chisq)
>      >> Time_Period            0.9105  2    0.63429
>      >> Assignment             2.1043  1    0.14689
>      >> Time_Period:Assignment 6.8486  2    0.03257 *
>      >>
>      >> But when I run the ANOVA for the zero-inflated component --
>      >> car::Anova(m1,
>      >> component = "zi") -- the output looks strange:
>      >>
>      >>                                  Chisq    Df Pr(>Chisq)
>      >> Time_Period            9.8876  3   0.019546 *
>      >> Assignment             9.4648  2   0.008805 **
>      >> Time_Period:Assignment 7.9605  2   0.018681 *
>      >>
>      >> Why would the degrees of freedom change? FWIW this df
>     discrepancy between
>      >> the conditional and ZI ANOVAs does *not* happen when running the
>     above
>      >> glmmTMB model with subsetted data frames based on different bird
>      >> guilds. It
>      >> also seems strange that between the Time Period and Assignment
>     terms the
>      >> smaller chi square value leads to greater statistical
>     significance. Do
>      >> you
>      >> agree that something seems wrong here or am I misunderstanding
>     what is
>      >> going on under the hood? Any ideas on how to troubleshoot?
>      >>
>      >> Thank you!
>      >>
>      >> -Elliot
>      >>
>      >>     [[alternative HTML version deleted]]
>      >>
>      >> _______________________________________________
>      >> R-sig-mixed-models using r-project.org
>     <mailto:R-sig-mixed-models using r-project.org> mailing list
>      >> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>     <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>
>      >
>      > _______________________________________________
>      > R-sig-mixed-models using r-project.org
>     <mailto:R-sig-mixed-models using r-project.org> mailing list
>      > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>     <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>
> 

-- 
Dr. Benjamin Bolker
Professor, Mathematics & Statistics and Biology, McMaster University
Director, School of Computational Science and Engineering
(Acting) Graduate chair, Mathematics & Statistics
> E-mail is sent at my convenience; I don't expect replies outside of 
working hours.



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