Spline plot after proc phreg

[TomHsiung]

proc phreg data=work.d202_raw plot(cl) = surv;
class genderC AF Hypertension CHF akiC ddiC indicationC CYP2C9code VKORC1Code / desc;
effect spl = spline(serumAlb / naturalcubic basis = bspline knotmethod = percentilelist(25 50 75));
model timeToFirstInrAttain*firstInrAttain(0) = bsa ageY genderC AF Hypertension CHF 
spl akiC ddiC indicationC CYP2C9code VKORC1Code / rl = both;
hazardratio AF;
hazardratio Hypertension;
run;

We know the method to integrate the spline effect of the variable serumAlb into the Cox model and we get the regression coefficients of the spl variables.

But we don't know how we plot the spline plot, which has a y-axis for hazard ratios and x-axis for the original variable of serumAlb.

I know there is an example but it's for another procedure (which outputs the spl variable of each observation to a new dataset) instead of proc phreg.

I'm looking forward to your reply and I appreciate your help.

[Ksharp]

Cann't you yield these HazzardRatio by using hazardratio statement and plot it by PROC SGPLOT？

libname x v9 'C:\Users\xiakeshan\Downloads';

ods output  HazardRatios= HazardRatios;
proc phreg data = x.whas500;
 effect Spl=spline(age/naturalcubic);
class gender;
model lenfol*fstat(0) = gender Spl ;
hazardratio 'Effect of 5-unit change in age across age' age / at(age = (30 to 100 by 5)) units=5;
/*
model lenfol*fstat(0) = gender|Spl ;
hazardratio 'Effect of gender across ages' gender / at(age=(30 to 100 by 5));
*/
run;
data HazardRatios;
 set HazardRatios;
 age=input(scan(Description,-1,'='),best.);
run;

proc sgplot data=HazardRatios noautolegend;
        band upper=WaldUpper lower=WaldLower x=age;
        series y=HazardRatio x=age/markers;
		yaxis label='Hazard Ratio';
        refline 1;
        title "Hazzard Ratio of 5-unit change in age across age";
        run;

[TomHsiung]

Thanks. Better to know each observation's spline format of the continuous variable.

Stata achieves this goal by (another dataset but it use the Cox model too):

mkspline ages = ageY, cubic nknots(4) displayknots

stcox ages* BSA genderC AF Hypertension CHF serumAlb akiC ddiC indicationC CYP2C9code VKORC1code AFandWarfarinHistory

list ageY ages* if ageY == 60

predictnl xb = _b[ages1]*(ages1 - 60) + _b[ages2]*(ages2 - 13.19037) + _b[ages3]*(ages3 - 0.3942672) + _b[CYP2C9code]*(1 - 0), ci(lo hi)

twoway (line uci ageY, sort lc(black black) lp(- -)) (line lci ageY, sort lp(dash_dot)) (line hr ageY, sort lc(black) lp(l)), scheme(sj) ytitle(Hazard ratio (CYP2C9)) xtitle(Age (years)) yline(1.00, lpattern(dot))

[Ksharp]

I don't understand. You can get DESIGN MATRIX by my code :
https://communities.sas.com/t5/SAS-Programming/how-to-convert-a-continuous-variable-age-as-cubic-spl...

And can score it for each obs by design matrix and estimated coefficient with @Rick_SAS blog:
https://blogs.sas.com/content/iml/2019/10/16/visualize-regression-splines.html

And this:

https://support.sas.com/kb/57/682.html

[TomHsiung]

Thanks for the feedback.

We think we have a two-step approach.

First, we use the procedure to transfer continuous to spl variables, e.g., from the variable of age to spl variables age1, age2, and age3. Later, the regression procedure will estimate the coefficients of age1, age2, and age3.

Next, we want to map each age with its corresponding spl variables age1, age2, and age3. If there are n observations for age, there should be n for age1, age2, and age3 too.

The second step is that we use the data of age1, age2, and age3, along with their coefficients to estimate the outcome (e.g., HRs) for each observation by hand. Finally, we plot the relationship of the variable age and the outcome via proc sgplot.

So far, we don't know how to get the age1, age2, and age3 for each observation for the proc phreg.

[Ksharp]

Maybe @SteveDenham  @StatDave  @Rick_SAS  could give you a hand.

[SteveDenham]

Sorry @Ksharp and @TomHsiung, I don't know off hand how to get the estimates at the knots for each subject. I would start with the OUTEST= option in the PROC PHREG statement, and go from there.

SteveDenham

[Season]

Take a look at the STORE statement of the PHREG procedure as well as the PLM procedure that utilizes things generated from the STORE statement. Some statistical graphs can be generated from the PLM procedure, but delicate modifications of them still necessitates the employment of the SG procedures like PROC SGPLOT.