I have count data (seconds of behaviour, n=145) with one explanatory categorical variable (Treatment) and 2 interaction terms (Sex and Time). I began by running a GLM with Poisson distribution:
model <- glm(Vigilance ~ Treatment * Sex * Time, family = "poisson", data = fallow)
The output gave significant results for all coefficients but one, and two were not defined because of singularities (so I will remove). Residual deviance = 2412.8 on 135 df and AIC was 2849.4. I checked for overdispersion with:
theta <- model$deviance / model$df.residual
which showed model is overdispersed - theta = ~17.8.
When I ran the same using quasipoisson, none of the coefficients were significant, yet when I ran the same with negative binomial, 2 coefficients were significant, residual deviance reduced to 161.7 on 135 df, AIC dropped to 950, and there was overdispersion - theta = 1.18.
Other options outlined by Zuur et al. (Mixed Effects Models and Extensions in Ecology with R), namely using drop1, are unavailable because my variables are included as interaction terms.
Is there anything else I should try, or can I safely use the NB distribution?