How to visualization the interaction and odds ratio/ risk ratio?

[Diels_O]

Dear all, I have a question of how to create the plot using continuous in the model as the independent variable and the odds ratio/ predicted probability as the dependent variable.

I try to use the effectplot statement in proc genmod, however, it seems to not work

proc genmod data = ITT1;
class enroll_order group(desc) score_0to1(desc)/PARAM = GLM;
model score_0to1 = group age/ dist = poisson link = log covB corrb;
effectplot fit(x = age);
repeated subject = enroll_order/type = un;
estimate "beta" group 1 -1/exp;
lsmeans group/diff exp cl;
run;

The score_0to1 is a binary responce(0 and 1 ), why this get negtive value?

How can I get plot like the this? (This plot shows the odds ratio for group changing with time)

[Ksharp]

/*
Check this URL:

https://support.sas.com/kb/24/188.html
http://support.sas.com/kb/37/344.html

*/

data Seizure;
input ID Count Visit Trt Age Weeks;
datalines;
104 11 0 0 31 8
104 5 1 0 31 2
104 3 2 0 31 2
104 3 3 0 31 2
104 3 4 0 31 2
106 11 0 0 21 8
106 15 1 0 21 2
106 13 2 0 21 2
106 30 3 0 21 2
106 32 4 0 21 2
234 2 0 1 35 8
234 1 1 1 35 2
234 4 2 1 35 2
234 3 3 1 35 2
234 12 4 1 35 2
235 2 0 1 30 8
235 11 1 1 30 2
235 40 2 1 30 2
235 32 3 1 30 2
235 22 4 1 30 2
236 12 0 1 30 8
236 1 1 1 37 2
236 4 2 1 37 2
236 3 3 1 37 2
236 2 4 1 37 2
;
data Seizure;
set Seizure;
if ID ne 207;
if Visit = 0 then do;
X1=0;
Ltime = log(8);
end;
else do;
X1=1;
Ltime=log(2);
end;
run;

proc gee data = Seizure;
class ID Visit trt;
model Count = age Trt  / dist=poisson link=log offset= Ltime;
repeated subject = ID /within=visit type=unstr covb corrw;
effectplot slicefit(x=age sliceby=trt) /  limits ;
run;

[Diels_O]

Thanks for your help! However, how could I plot the ODDS ratio or risk ratio of a categorical changing with the continous variable?

[Ksharp · Accepted solution]

You can get risk ratio from URL I post.
Or you could use
lsmeans A / at x1=1.2 ODDSRATIO diff cl ;
lsmeans A / at x1=1.4 ODDSRATIO diff cl ;
lsmeans A / at x1=1.6 ODDSRATIO diff cl ;
...........
to get oddsratio, after that use PROC SGPLOT to plot it.

and maybe @StatDave know more things.

[Diels_O]

Thank you very much, and I hope @StatDave and @Rick_SAS have other ways to visualize odds/risk ratio changing with continuous variable

[StatDave]

If you want a plot of how the predicted probability, relative risk, risk difference, or odds ratio changes over the levels of age, then you need to include the group*age interaction in the model so that the model allows for change over age. Since your response is binary, then this is really much easier to do with a logistic model, which is the only model you can use if you want odds ratios. And while PROC LOGISTIC would normally be the best choice, it is more convenient to use PROC NLMIXED to fit the model since you can then use its PREDICT statement to save data sets of the predicted probabilities and the estimated relative risks and odds ratios. If you aggregate the data for the group-age combinations so that each observation has a count of the events (r) and a count of the total number of cases (n), then the following code fits the model, saves the needed data sets, and produces each of the plots against age. The difference of probabilities could similarly be computed and plotted.

proc nlmixed;
   p=logistic(b0 + b1*(group='A') + b2*age + b3*age*(group='A'));
   model r ~ binomial(n,p);
   predict p out=pred;
   predict  logistic(b0+b1+b2*age+b3*age) / logistic(b0+b2*age) out=rr;
   predict exp(b0+b1+b2*age+b3*age)/exp(b0+b2*age) out=or;
   run;
proc sort data=pred; by group age; run;
proc sgplot data=pred;
   band upper=upper lower=lower x=age/group=group transparency=.8;
   series y=pred x=age/group=group;
   yaxis label="Predicted probability";
   run;
proc sort data=rr; by age; run;
proc sgplot data=rr noautolegend;
   band upper=upper lower=lower x=age/transparency=.8;
   series y=pred x=age;
   yaxis label="Relative risk";
   run;
proc sort data=or; by age; run;
proc sgplot data=or noautolegend;
   band upper=upper lower=lower x=age/transparency=.8;
   series y=pred x=age;
   yaxis label="Odds ratio";
   run;

[Diels_O]

Thank you Dave , yet I met problems when preparing my data to proc nlmixed, and the score_0to1 is a binary response, n is the total number of events and r is the number of predicted event. How shoud I transform my data? Thank you very much!

data q2; 
set q1; 
n = _n_;
r =(score_0to1 = 1);  
output;  
run;

[StatDave]

The NLMIXED code I gave is appropriate for aggregated data rather than individual data. If your data is individual level data, then you can either aggregate it with a procedure like MEANS or FREQ to save the n and r counts, or you could change the NLMIXED MODEL statement to be appropriate for individual level data. To change NLMIXED, your response variable needs to be numeric and coded 1 for the event and 0 for the nonevent. Then change the MODEL statement to: model y ~ binary(p);

[Diels_O]

Thank you Dave! Wish you a luck day!