|
Kevin,
Forgot to include my second example with GLIMMIX and dist=poisson for comparison .. that is, invoke the robust standard error option with a random statement, even if there is only one obs per subject (identified here with id). In GLIMMIX, the robust error variances can be estimated by using the RANDOM statement as indicated with the subject identifier.
Assuming you have an id variable:
ods select estimates;
proc GLIMMIX data = indat;
CLASS gender id;
MODEL resp = gender / dist = poisson link = log ;
RANDOM residual / subject = id ;
estimate 'Gender' gender 1 -1 / exp cl;
title 'Relative risk estimation by Poisson regression with robust error variance';
run;
Relative risk estimation by Poisson regression with robust error variance (close to previous example)
The GLIMMIX Procedure
Standard Exponentiated Exponentiated Exponentiated
Label Estimate Error Estimate Lower Upper
Gender 0.4612 0.1945 1.5860 1.0782 2.333
Robin High
UNMC
-----Original Message-----
From: SAS(r) Discussion [mailto:SAS-L@LISTSERV.UGA.EDU] On Behalf Of Kevin F. Spratt
Sent: Wednesday, March 28, 2012 8:09 AM
To: SAS-L@LISTSERV.UGA.EDU
Subject: Problem with poisson regression in GLIMMIX
I have run a glimmix procedure using poisson regression to estimate
relative risks
associated with various binary predictors. I show two models below.
* full model;
PROC GLIMMIX DATA=save.demos METHOD=RSPL;
CLAsS pT_age_c2 pt_sex pt_emp pt_ret smoke
diabe hp_card hp_renu hp_psyc ;
MODEL dv = pt_age_c2 pt_sex pt_emp pt_ret
smoke diabe hp_card hp_renu hp_psyc / LINK=LOG
S DIST=POISSON;
LSMEANS pt_age_c2 / ILINK CL OM PDIFF;
LSMEANS pt_sex / ILINK CL OM PDIFF;
LSMEANS pt_emp / ILINK CL OM PDIFF;
LSMEANS pt_ret / ILINK CL OM PDIFF;
LSMEANS smoke / ILINK CL OM PDIFF;
LSMEANS diabe / ILINK CL OM PDIFF;
LSMEANS hp_card / ILINK CL OM PDIFF;
LSMEANS hp_renu / ILINK CL OM PDIFF;
LSMEANS hp_psyc / ILINK CL OM PDIFF;
run;
* one variable model;
PROC GLIMMIX DATA=save.demos METHOD=RSPL;
CLAsS pT_age_c2 ;
MODEL dv = pt_age_c2 / LINK=LOG S DIST=POISSON;
LSMEANS pt_age_c2 / ILINK CL OM PDIFF;
run;
My problem:
For the one variable model, the results produce the RR that is
consistent with what
I get with proc freq, but the confidence intervals for the glimmix and freq
procedures are not even close. This is probably due to the fact that the
dependent variable (dv), although 0/1, is not necessarily a rare event.
So, for the one variable models I can get my results using proc freq but for
the multi-variable model I plan to figure out a way to produce
non-parametric confidence
intervals using bootstrapping. Before I go to this trouble, I'm wondering if
there is a more straightforward way to do this within glimmix. Can I
changes some
options in the glimmix model that will "magically" solve this problem
or are the
model coefficients in the multi-variable Poisson model not able to
provide appropriate
adjusted relative risk estimates?
I have run the models using logistic regression and get odds ratios,
but, as we know
odds ratios and relative risks can be quite different when the 0/1
dependent variable
is not relative "rare." The differences in these data, at the
univariate level
suggest that ORs are markedly over-estimating variable effects
compared to the RRs.
As always, any help greatly appreciated.
PS: if the RR estimates are correct but bootstrapping to obtain
accurate non-parametric
confidence intervals for the values is needed, would this be a
reasonable paper to submit
for next year's SGF? If so, I would be willing to collaborate. I'm
sure I can figure
out how to do it, I'm less sure that I can produce a "really cool"
macro that would be
of use to the community.
______________________________________________________________________
Kevin F. Spratt, Ph.D.
Department of Orthopaedic Surgery
Dartmouth Medical School
One Medical Center Drive
DHMC
Lebanon, NH USA 03756
(603) 653-6012 (voice)
(603) 653-6013 (fax)
Kevin.F.Spratt@Dartmouth.Edu (e-mail)
Data is not information;
Information is not knowledge;
Knowledge is not understanding;
Understanding is not wisdom.
- Cliff Stoll and Gary Schubert
_______________________________________________________________________
| 
