Plotting PROC PHREG after Multiple Imputation

fikel — Thu, 03 Jun 2021 20:12:01 GMT

Hello, I had previously asked about plotting Proc PHREG for restricted splines analysis and was able to use this code:

proc phreg data = have;
title "Colorectal-cancer incidence attempt";
	
effect resid_total_polys = spline(resid_total_poly / basis=tpf(noint) NATURALCUBIC details knotmethod=Percentilelist(5 35 65 95));

class menopause_with_males (ref="0") race (ref="1") hhincome (ref="1")
		enrollment_source (ref="C") fh_colorectalcancer (ref="0")
		BMI_category_expanded (ref="0") alcohol_category (ref="0") smokestatus_packyear (ref="0") hei10_category (ref="0");
 
model eof_age_CRC*Inc_CRC(0) = resid_total_polys 
		menopause_with_males race hhincome enrollment_source 
		BMI_category_expanded totalactivitymethr smokestatus_packyear Comorbidity_Index 
		ffq_kcal hei10_category alcohol_category fh_colorectalcancer
		/ entry=enrollment_agemonths rl=wald;
		
store work.coxr;
run;



%macro est(ref=212, start=0, end=3000, by=1);

%Do i = 1 %To %eval(%SysFunc( Ceil( %SysEvalF( ( &End - &Start ) / &By ) ) ) +1) ;

   %Let value=%SysEvalF( ( &Start - &By ) + ( &By * &I ) ) ;

    estimate "&value." resid_total_polys [-1, &ref] [1, &value] / exp cl;

    %end;

%mend est;

ods exclude all;

ods dataset Estimates=Estimates;

proc plm restore=work.coxr;

    %est(ref=212, start=0, end=3000, by=1);

run;

ods exclude none;

data estimates;

    set estimates;

    resid_total_poly=label*1;

run;

proc sgplot data=estimates NOAUTOLEGEND Noborder;

    Series y=ExpEstimate x=resid_total_poly / LINEATTRS=(color=black thickness=3px);

    Series y=LowerExp x=resid_total_poly / LINEATTRS=(pattern=ShortDash color=Black THICKNESS=1);

    Series y=UpperExp x=resid_total_poly / LINEATTRS=(pattern=ShortDash color=Black THICKNESS=1);
    
    band x=resid_total_poly upper=UpperExp lower=LowerExp / fillattrs=(color=graydd) transparency=.50;

    REFLINE 1 / axis=y;

    REFLINE 3000 / Axis=X LINEATTRS=(pattern=ThinDot color=Black THICKNESS=1);;

    yaxis Values=(0.5 0.8 1 1.2 1.5) Label="Hazard ratio" Type=LOG LABELATTRS=(weight=BOLD);

    xaxis min=0 VALUES=(0 to 3000 by 100) LABELATTRS=(weight=BOLD);

run;

To produce this graph:

Now, however, I am using multiple imputation to handle missing values, but I am not sure how to adapt the code to be able to handle this. This is the code I have:

proc phreg data = mult_imp.mi_mvn namelen=50;
title "All-cancer incidence attempt";

effect resid_total_polys = spline(resid_total_poly / basis=tpf(noint) NATURALCUBIC details knotmethod=Percentilelist(5 35 65 95));

model eof_age_CRC*Inc_CRC(0) = resid_total_polys
		menopause_males_dum2 menopause_males_dum3 racedum2 racedum3 incomedum2 incomedum3 incomedum4 
		enrolldum2 BMI_dum2 BMI_dum3 totalactivitymethr packyears_dum2 packyears_dum3 packyears_dum4 packyears_dum5 Comorbidity_Index 
		ffq_kcal hei10_category alcohol_dum2 alcohol_dum3 fh_colorectalcancer 
		/ entry=enrollment_agemonths rl=wald;
by _imputation_;

store work.coxr;

run;
quit;

%macro est(ref=212, start=0, end=3000, by=1);
%Do i = 1 %To %eval(%SysFunc( Ceil( %SysEvalF( ( &End - &Start ) / &By ) ) ) +1) ;
   %Let value=%SysEvalF( ( &Start - &By ) + ( &By * &I ) ) ;
    estimate "&value." resid_total_polys [-1, &ref] [1, &value] / exp cl;
    %end;
%mend est;

ods exclude all;
ods dataset Estimates=Estimates;

proc plm restore=work.coxr;
    %est(ref=212, start=0, end=3000, by=1);
run;

ods exclude none;
data estimates;
    set estimates;
    resid_total_poly=label*1;
run;

proc sgplot data=estimates NOAUTOLEGEND Noborder;
    Series y=ExpEstimate x=resid_total_poly / LINEATTRS=(color=black thickness=3px);
    Series y=LowerExp x=resid_total_poly / LINEATTRS=(pattern=ShortDash color=Black THICKNESS=1);
    Series y=UpperExp x=resid_total_poly / LINEATTRS=(pattern=ShortDash color=Black THICKNESS=1); 
    band x=resid_total_poly upper=UpperExp lower=LowerExp / fillattrs=(color=graydd) transparency=.50;
    REFLINE 1 / axis=y;
    REFLINE 3000 / Axis=X LINEATTRS=(pattern=ThinDot color=Black THICKNESS=1);;
    yaxis Values=(0.5 0.8 1 1.2 1.5) Label="Hazard ratio" Type=LOG LABELATTRS=(weight=BOLD);
    xaxis min=0 VALUES=(0 to 3000 by 100) LABELATTRS=(weight=BOLD);

run;

Where the set "mult_imp.mi_mvn" is the dataset obtained after running Proc MI with 5 imputations. I am stuck trying to figure out how to get a single graph for this model that combines the 5 imputed sets. However, this produces this plot:

I was thinking that I somehow need to use proc mianalyze to combine the estimates for the models from each imputation but was not sure how to accomplish this.

topic Plotting PROC PHREG after Multiple Imputation in Statistical Procedures

Plotting PROC PHREG after Multiple Imputation