PROC GENMOD - Probabilities for each "response level" by "treatment"

[Nick_Serao]

Hi,

I would like to get the probability of having any "response level" for each level of "treatment" in my data, and not only the "probabilities of levels of y having LOWER Ordered Values in the

response profile table." that SAS automatically outputs.

(I'm simplifying a lot the real scenario just to have an idea).

Let's say I would like to see the effect of 3 treatments (T1, T2 and T3) on a given health status (5 categories: low, mid-low, med, mid-high and high), correcting for the covariate Weight (W).

I'm running an ordinal multinomial model:

----------------------------------------------------------------------------------

proc genmod data=health;

class T;

model Status = T W/dist=multinomial link=clogit type3;

estimate 'T1 - T2' T 1 -1/exp;

run;

---------------------------------------------------------------------------

This (estimate statement) code gives that the odds of T1 being in lower health categories is 2.8 the odds of T2 being in lower health categories.

Although this is interesting, that's not what I want. What I really want are all the probabilities of being in each of the 5 Health Categories for each of the 3 Treatments:

P(low|T1)

P(mid-low|T1)

...

P(high|T3).

Is is possible to have this in Proc Genmod? If so, can anyone help me?

Thanks,

Nick

[SteveDenham · Accepted solution]

Not sure this will give what you want, but it might lead to it:

proc genmod data=health;

class T;

model Status = T W/dist=multinomial link=clogit type3 aggregate=T;

estimate 'P(low|T1)' intercept 1 T 1 0 0 / exp e;

estimate 'P(mid-low, low|T1) intercept 0 1  T 1 0 0 / exp e;

estimate 'P(med, mid-low, low|T1) intercept 0 0 1  T 1 0 0 / exp e;

<put more estimates in here, as needed, following these examples>

ods output estimates=estimates;

run;

I based it on the cumulative logit.  There would have to be a post-processing step on the estimates dataset, as the mean estimates are the cumulative probabilities.  The estimated probabilities per level would have to be calculated by difference, with the probability of the 'high' level as 1 - P(med-high, med, mid-low, low|T1).  Unfortunately, this will not give standard errors for the probabilities, but if you only need point estimates, this might be an approach you could use.

Steve Denham

[Nick_Serao]

Thanks Steve...

Sorry for the late reply, I was travelling...

This was actually the first approach I tried, but I was not having much success with these estimates because my model is a little bit more complex, but looking at your code, I could finally understand what I had to do.

But, as you said, I could only get the point estimates.

Thanks!

[SteveDenham]

Nick,

What about interval estimates?  You could get these by adding the cl option after the slash.  Approximate standard errors are a lot harder.  While the delta method is by far to be preferred, one approach would be to estimate the standard error as sqrt( estimate * (1 - estimate) / df)).  This means post-processing of the estimates dataset using ODS output estimates=estimates;

Steve Denham