Dear list, 

I'm just starting to learn how to use glmmPQL and I have a very basic question. Unfortunately no one around me works with glmmPQL and the answers I found online were not directly related to my question. 

I work with habitat use data (summarized binomial data) which I fitted in a mixed model of the binomial family with a random factor accounting for repeated sampling of the same individual. I should specify that this is not a habitat selection study, and this particular analysis is peripheral to the whole project, which is a wildlife toxicology study. We simply want basic descriptions of the data - effects of sex, year, etc on the proportion of time spent in each habitat.

The code is as follow:
Colony.glmmSexYear <- glmmPQL(cbind(ColoYes.allnoF24, ColoNo.allnoF24) ~ year.coded + sex, random= ~1| BirdID, family=binomial(), data=mydata)
summary(Colony.glmmSexYear)

Here is the fixed effects section of my output: 

Fixed effects: cbind(ColoYes.allnoF24, ColoNo.allnoF24) ~ year.coded + sex
 
                     Value Std.Error DF   t-value p-value
(Intercept)      1.4442270 0.1824173 74  7.917161  0.0000
year.codedY2011 -0.1733713 0.1864744 74 -0.929732  0.3555
year.codedY2012 -0.3004284 0.2027411 74 -1.481833  0.1426
sexM             0.4403108 0.1507471 74  2.920857  0.0046


I simply want to report the fact that year (categorical, three levels) is not significant. In traditional linear models I would report the p value for treatment "year" from an anova table. But my understanding is that producing an anova from glmmPQL is not possible (?). Reporting those individual p values of my three years separately seems clumsy and out of place...

What is the proper way to do this ? 

Thank you very much for your help,

Marie
	[[alternative HTML version deleted]]

