No standard error for the overdispersion factor in a simple Poisson regression with Proc Glimmix

[WillTheKiwi]

I have previously used Proc Glimmix with repeated measurement, hence I have used random SubjectID; or similar, and I have allowed for overdispersion with random _residual_;. Previously the overdispersion factor has appeared in the covariance parameters with a standard error (and confidence limits), but now I am just running a simple Poisson regression with one observation per subject, so there are no random effects. The overdispersion factor now appears in the parameter estimates as "Residual", along with the fixed effects, but it does not have a standard error. How come, and what can I do to get a standard error (and of course, confidence limits? Here's my code:

proc glimmix data=dat2;
class Sex;
model AntiN = Sex/s noint link=log dist=poisson;
estimate "Mean Males" Sex 1 0 /exp alpha=0.1;
estimate "Mean Females" Sex 0 1 /exp alpha=0.1;
estimate "Effect of Sex, F/M" Sex -1 1 /exp alpha=0.1;
random _residual_; 
run;

I seem to recall encountering this problem many years ago, but I can't find anything in my files about it.

Thank you!

Will

[SteveDenham]

I am not sure what is going to happen with this code, but it might trigger something on your side. What happens if you change the RANDOM _residual_ to RANDOM _residual_/group=sex ? I would expect to see two estimates for residual error in that case. To actually see them and their confidence bounds you need to a COVTEST statement. In this case it would look like:

covtest homogeneity /cl;

I don't see a good method in GLIMMIX to get the SE for the residual variance, but that doesn't mean there isn't one. GENMOD does give an SE for the scale parameter in a Poisson regression, that may or may not have an SE>0, depending on the data.

SteveDenham

[WillTheKiwi · Accepted solution]

Thanks for your feedback, Steve. I meant to say that I had already tried random _residual_/group=Sex, and weirdly, Glimmix ignored the statement altogether by not estimating the inflation factors and by making no warning or other statement in the LOG about ignoring the statement. However, your suggestion of adding covtest did the trick. So I used:

random _residual_/group=Sex; 
covtest homogeneity/cl(alpha=0.1);

ods output covparms=cov;

The cov dataset had the inflation factors and their standard errors and confidence limits.

I tried to get a single inflation factor by dropping off the group=Sex, but it didn't work, and again there was no indication of why in the LOG. So I tricked it by making a dummy class variable that had the same value for every observation in the data set, then using random _residual_/group=Dummy. It worked perfectly. But it seems to me that something needs to be done to Proc Glimmix and to its documentation so others don't have this problem.

I also tried the analysis with Proc Genmod, as you suggested, but it gave a different value for the scale parameter, way too large, and squaring it didn't make it right. Maybe it was because I actually have under-dispersion with these data, but the values given by Glimmix were correct, because they were exactly the same as the variance divided by the mean for the males and for the females.

So again, THANK YOU!