This is an extension of the difference in difference analysis for longitudinal binary data discussed in this note. With the MONTH variable as a continuous predictor, one way to estimate the slopes is to compute average marginal effects within the pre- and post-intervention periods. This can be done with the Margins macro, which is also illustrated in the above note for difference in difference analysis. These statements fit a GEE model allowing for repeated binary measures on each ID and plot the fitted values on the mean (probability) scale and on the link (log odds) scale. proc gee data=mydata; class id iv; model y(event="1") = iv|month / dist=bin; repeated subject=id; effectplot slicefit(x=month sliceby=iv); effectplot slicefit(x=month sliceby=iv)/link; run; The Margins macro can also fit the model and request estimates of the mean (event probability) at each month. %margins(data=mydata, class=id iv, response=y, roptions=event="1", model=iv|month, dist=binomial, geesubject=id, margins=month iv) You can also use the Margins macro to compute the slope in each period as the average marginal effect over the months in the period.  %margins(data=mydata, class=id iv, response=y, roptions=event="1", model=iv|month, dist=binomial, geesubject=id, effect=month, within=iv=0) %margins(data=mydata, class=id iv, response=y, roptions=event="1", model=iv|month, dist=binomial, geesubject=id, effect=month, within=iv=1) If "PRE/POST GAP" means that you want to estimate the change over the two months immediately on either side of the intervention, then this can be done by estimating a contrast on the monthly margins. This code assumes there are 6 months with intervention between months 3 and 4. data c; length label f $32767; infile datalines delimiter='|'; input label f; datalines; gap: month4-month3 | 0 0 -1 1 0 0 ; %margins(data=mydata, class=id iv, response=y, roptions=event="1", model=iv|month, dist=binomial, geesubject=id, margins=iv month, contrasts=c) ... View more

Use the LSMEANS statement on the 3-way interaction to show the estimates and order of the means. Then use an LSMESTIMATE statement to test the joint hypothesis of the pairwise differences of the ones you want to compare.  For example, the following compares the first three and abbreviation treatment as t, visit as v, and age as a. It will be a 2 df test. For your case, you'll have a lot more zeros in each part to pick out the right means. Be careful to match your coefficients to the order of from the LSMEANS statement. lsmeans t*v*a; lsmestimate t*v*a '111=112' 1 -1 0, '111=113' 1 0 -1, '112=113' 0 1 -1 / joint; ... View more

Producing a 3D plot for two continuous variables is very easy - in fact, it is the default of the EFFECTPLOT statement in that situation. The following produces a contour plot that should be pretty straightforward for any audience to understand.  proc logistic data=expi2; model criterion (event= '1') = predictor1 predictor2 ; effectplot; run; ... View more

In general when you think "adjusted" predictions, you are talking about either LS-means, as available from the LSMEANS statement, or predictive margins, as from the Margins macro. For LS-means, you can add the LSMEANS statement in your PROC LOGISTIC step. Note that least squares means are simply linear combinations of the model parameters. The adjustment is in the coefficients used on the predictors other than the one of interest. You can see these coefficients by adding the E option, and you can adjust these, if needed, with the OM= option. The ILINK option applies the inverse of the logit link to get predicted probabilities (what you are calling a "rate").   lsmeans sex / ilink e; The problem with the LSMEANS statement is that it can only provide adjusted predicted probabilities for the first k-1 levels of a response with k levels (since there are only k-1 logits). An obvious way to get adjusted estimates for all response levels is to simply average the predicted probabilities from the model for each response level. You can do that by adding an OUTPUT statement to save the predicted probabilities for all observations.   output out=preds predprobs=individual; and then average them     proc means data=preds mean; class sex; var ip:; run; The other possibility is predictive margins. Unfortunately, the Margins macro cannot be used with a multinomial response model. However, you can easily compute point estimates of the predictive margins since they are simply averages of predicted probabilities when all observations are fixed at one level of the predictor. You can do that by creating versions of your data with all observations set to Men or Women. data m; set taxes.analyse; sex='Men'; run; data w; set taxes.analyse; sex='Women'; run; then add SCORE statements to apply the fitted model to each of these data sets. score data=m out=mpreds; score data=w out=wpreds; and then average each proc means data=mpreds mean; var p:; run; proc means data=wpreds mean; var p:; run;     ... View more