Hello SAS Community,
I am currently working on a survival analysis and would like to create a summary table that closely resembles the Kaplan-Meier analysis output one typically sees in publications. My goal is to have a table that shows, for each treatment group, the time in months, overall survival probability, number at risk, cumulative number of censored events, and number of events remaining.
Here's an example of the format I'm targeting (similar to an attached image):
+---------------+-------------------------+----------------+-------------------+-------------------+------------+ | Time (months) | Overall Survival Prob. | Number at risk | Number censored | Number of events | Treatment | +---------------+-------------------------+----------------+-------------------+-------------------+------------+ | 0 | 1 | 406 | 0 | ... | Treatment1 | | 3 | 0.83795 | 312 | 32 | ... | Treatment1 | | ... | ... | ... | ... | ... | ... | +---------------+-------------------------+----------------+-------------------+-------------------+------------+
Here's what I've done so far using PROC LIFETEST
to generate survival estimates:
ods output ProductLimitEstimates=SurvivalEstimates;
proc lifetest data=adtteqs method=KM plots=survival TIMELIST=(3, 6, 9, 12, 18, 24, 36, 48);
time avalM * CNSR(1);
strata trtp;
run;
ods trace off;
Following that, I've attempted to calculate the cumulative number of censored events:
/* Sort the SURVIVALESTIMATES by treatment group and time for the cumulative calculation */
proc sort data=SurvivalEstimates out=sorted_estimates;
by trtp TimeList;
run;
/* Calculate the cumulative count of censored subjects */
data cumulative_censored;
set sorted_estimates;
by trtp TimeList;
retain CumulativeCensored 0;
if first.trtp then CumulativeCensored=0;
CumulativeCensored + Censor;
run;
However, I am struggling with the correct calculation of cumulative censored counts and how to incorporate the total number of patients per treatment group into this table. Ideally, I would like the censored count to start at 0 and then accumulate as shown in the example table above.
Can anyone provide guidance on how to adjust my code to achieve this format, particularly on ensuring the cumulative count of censored events is calculated correctly for each treatment group and time point? Can this be extracted from proc lifetest directly? Should I be using LIFETEST CensoredSummary or Survival Plot instead? I've been trying a few things but the numbers just don't match the expected number of subjects per treatment group although it is a questionnaire dataset that is creating a high number of events (Time to First deterioration).
Thank you so much for any help or insights you can provide!
@sbxkoenk @WarrenKuhfeld Thank you both for your valuable inputs. I was able to finally get what I needed by running survival plot and limit estimates for both CNSR(0) and CNSR(1) extracting the numbers that I needed for my table. I will definitely take a closer look at this macro and documentation for future reference! Amazing work.
/* First run for cnsr(0) */
ods trace on;
ods graphics on;
ods listing close;
ods output CensoredSummary = CensoredSummary0
ProductLimitEstimates = SurvivalEstimates0
SurvivalPlot = SurvivalPlot0;
proc lifetest data=adtteqs plots=survival(atrisk= 0 to 48 by 3) timelist=(0,3,6,9,12,15,18,21,24,27,30,33,36,39,42,45,48) reduceout outsurv=survdata0;
time avalM * cnsr(0);
strata trtp;
by &Subgrp;
run;
/* Close the ODS destination to save the first set of outputs */
ods trace off;
ods graphics off;
ods listing;
/* Second run for cnsr(1) */
ods trace on;
ods graphics on;
ods listing close;
ods output CensoredSummary = CensoredSummary1
ProductLimitEstimates = SurvivalEstimates1
SurvivalPlot = SurvivalPlot1;
proc lifetest data=adtteqs plots=survival(atrisk= 0 to 48 by 3) timelist=(0,3,6,9,12,15,18,21,24,27,30,33,36,39,42,45,48) reduceout outsurv=survdata1;
time avalM * cnsr(1);
strata trtp;
by &Subgrp;
run;
/* Close the ODS destination to save the second set of outputs */
ods trace off;
ods graphics off;
ods listing;
data survivalplot (drop=Atrisk stratum_ survival event censored stratum tAtrisk);
set survivalplot0;
where cmiss(tAtRisk) = 0;
Number_At_Risk =atrisk;
stratum_=stratum;
trtp=Stratum_;
run;
proc sort data= survivalplot(rename=stratumnum=stratum);
by paramcd stratum time;
data SurvEst_nevents (drop=Timelist Failed Left survival censor Failure StdErr AvalM);
format TimeList avalM 8.;
set SurvivalEstimates1;
Time = Timelist;
Number_of_Events =Failed;
Overall_Survival= survival;
run;
proc sort data= SurvEst_nevents;
by paramcd stratum time;
data SurvEst_ncensevents (drop=Timelist Failed Left survival censor Failure StdErr AvalM);
format TimeList avalM 8.;
set SurvivalEstimates0;
Time = Timelist;
Number_Censored =Failed;
run;
proc sort data= SurvEst_ncensevents;
by paramcd stratum time;
data merged_data(drop=stratum);
retain paramcd time Overall_Survival Number_At_Risk Number_Censored Number_of_Events Trtp;
merge survivalplot (in=a)
SurvEst_nevents (in=b)
SurvEst_ncensevents (in=c);
by paramcd stratum time;
if a & b & c;
run;
Take the VALung dataset here -->
SAS® 9.4 and SAS® Viya® 3.5 Programming Documentation
SAS/STAT 15.3 User's Guide
The LIFETEST Procedure
Example 77.1 Product-Limit Estimates and Tests of Association
https://go.documentation.sas.com/doc/en/pgmsascdc/9.4_3.5/statug/statug_lifetest_examples01.htm
Is this what you are after?
ods output ProductLimitEstimates=SurvivalEstimates
CensoredSummary =CensoredSummary ;
proc lifetest data=VALung method=KM plots=survival
TIMELIST=(10, 50, 100, 150, 200, 250, 300, 400, 500)
outsurv=abc /*REDUCEOUT*/ ;
time SurvTime * censor(1);
strata Therapy;
run;
ods trace off;
proc format library=work;
value intvs
low - 10 = '010'
10 <- 50 = '050'
50 <- 100 = '100'
100 <- 150 = '150'
150 <- 200 = '200'
200 <- 250 = '250'
250 <- 300 = '300'
300 <- 400 = '400'
400 <- 500 = '500'
500 - high = 'high' ;
run;
proc means data=abc sum nway noprint;
CLASS Therapy SurvTime;
format SurvTime intvs.;
VAR _CENSOR_;
output out=xyz(drop=_TYPE_) sum= / autoname;
run;
/* Calculate the cumulative count of censored subjects */
data cumulative_censored;
set xyz;
retain CumulativeCensored 0;
by Therapy SurvTime;
if first.Therapy then CumulativeCensored=0;
CumulativeCensored = CumulativeCensored + _CENSOR__Sum;
run;
/* end of program */
Koen
Hi @sbxkoenk for the specific task I just want the information from the SurvivalPlot dataset after playing around with code. However I need time to be 0 to 48 by 3 to represent values at 3 month intervals. can I do that or do I need to do pre/post processing?. Currently it has individual records for every AVALM (AVAL in Months). Can timelist be applied to survivalplot dataset?
ods trace on;
ods graphics on;
ods listing close;
ods output SurvivalPlot = SurvivalPlot;
proc lifetest data = adtteqs plots=survival TIMELIST=(0,3, 6, 9, 12, 18, 24, 36, 48 );;
time avalM * cnsr(1);
strata trtp;
by &Subgrp;
run;
ods trace off;
ods graphics off;
This is close but I still need cumulative number of events
data adtteqs;
set adtteqs_ (drop=trtp);
avalM = aval / 365.25 * 12;
avalM = round(avalM, 3);
where &selpop. &filter.;
run;
proc sort; by &Subgrp usubjid; run;
ods trace on;
ods graphics on;
ods listing close;
ods output SurvivalPlot = SurvivalPlot;
proc lifetest data = adtteqs plots=survival TIMELIST=(0,3, 6, 9, 12, 18, 24, 36, 48 );;
time avalM * cnsr(1);
strata trtp;
by &Subgrp;
run;
ods trace off;
ods graphics off;
data SurvivalPlot;
set SurvivalPlot;
where cmiss(Survival) = 0;
run;
It might be easier
Like this example:
SAS Help Center: Product-Limit Estimates and Tests of Association
SAS Help Center: Example Code
Koen
I am rarely using PROC LIFETEST myself, so ... I am just trying ... and finding out.
Suppose you have code like this:
ods output SurvivalPlot = SurvivalPlot;
ods output CensoredSummary = CensoredSummary;
proc lifetest data = VALung plots=survival(atrisk=0 to 580 by 30)
TIMELIST=(0 to 580 by 30) REDUCEOUT
outsurv=aaaa;
time SurvTime * censor(1);
strata Therapy;
*by &Subgrp;
run;
It appears the SurvivalPlot dataset is not affected by:
, but the outsurv= dataset "aaaa" is !
This being said, there's plenty of literature ... with code examples ... around changing the (template of the) The Kaplan Meier survival plot. I guess you can also change the binning of survival time.
Just do some searches on Google like:
Good luck,
Koen
The KM plot modification work was a work in progress for many years. If you search all the sites mentioned in the reply, you will find old conference papers and other obsolete material. The last set of work is the SAS/STAT 15.1 chapter. https://support.sas.com/documentation/onlinedoc/stat/151/kaplan.pdf . I believe this is the final word on this topic (at least from SAS) as there is no one left there to work on it any more. In that final chapter, I wrote about every user-requested graph modification of which I was informed.
Also, independent of my work, Jefferey Meyers has done some great work on this topic.
Thanks @WarrenKuhfeld. That's a useful addition!
For the work of Jeffrey Meyers on this topic (customizing KM plot),
... see here :
Home > SAS Communities Library > Kaplan-Meier Survival Plotting Macro %NEWSURV
Kaplan-Meier Survival Plotting Macro %NEWSURV
Started 07-20-2018 | Modified 12-01-2022
by @JeffMeyers
https://communities.sas.com/t5/SAS-Communities-Library/Kaplan-Meier-Survival-Plotting-Macro-NEWSURV/...
Koen
@sbxkoenk @WarrenKuhfeld Thank you both for your valuable inputs. I was able to finally get what I needed by running survival plot and limit estimates for both CNSR(0) and CNSR(1) extracting the numbers that I needed for my table. I will definitely take a closer look at this macro and documentation for future reference! Amazing work.
/* First run for cnsr(0) */
ods trace on;
ods graphics on;
ods listing close;
ods output CensoredSummary = CensoredSummary0
ProductLimitEstimates = SurvivalEstimates0
SurvivalPlot = SurvivalPlot0;
proc lifetest data=adtteqs plots=survival(atrisk= 0 to 48 by 3) timelist=(0,3,6,9,12,15,18,21,24,27,30,33,36,39,42,45,48) reduceout outsurv=survdata0;
time avalM * cnsr(0);
strata trtp;
by &Subgrp;
run;
/* Close the ODS destination to save the first set of outputs */
ods trace off;
ods graphics off;
ods listing;
/* Second run for cnsr(1) */
ods trace on;
ods graphics on;
ods listing close;
ods output CensoredSummary = CensoredSummary1
ProductLimitEstimates = SurvivalEstimates1
SurvivalPlot = SurvivalPlot1;
proc lifetest data=adtteqs plots=survival(atrisk= 0 to 48 by 3) timelist=(0,3,6,9,12,15,18,21,24,27,30,33,36,39,42,45,48) reduceout outsurv=survdata1;
time avalM * cnsr(1);
strata trtp;
by &Subgrp;
run;
/* Close the ODS destination to save the second set of outputs */
ods trace off;
ods graphics off;
ods listing;
data survivalplot (drop=Atrisk stratum_ survival event censored stratum tAtrisk);
set survivalplot0;
where cmiss(tAtRisk) = 0;
Number_At_Risk =atrisk;
stratum_=stratum;
trtp=Stratum_;
run;
proc sort data= survivalplot(rename=stratumnum=stratum);
by paramcd stratum time;
data SurvEst_nevents (drop=Timelist Failed Left survival censor Failure StdErr AvalM);
format TimeList avalM 8.;
set SurvivalEstimates1;
Time = Timelist;
Number_of_Events =Failed;
Overall_Survival= survival;
run;
proc sort data= SurvEst_nevents;
by paramcd stratum time;
data SurvEst_ncensevents (drop=Timelist Failed Left survival censor Failure StdErr AvalM);
format TimeList avalM 8.;
set SurvivalEstimates0;
Time = Timelist;
Number_Censored =Failed;
run;
proc sort data= SurvEst_ncensevents;
by paramcd stratum time;
data merged_data(drop=stratum);
retain paramcd time Overall_Survival Number_At_Risk Number_Censored Number_of_Events Trtp;
merge survivalplot (in=a)
SurvEst_nevents (in=b)
SurvEst_ncensevents (in=c);
by paramcd stratum time;
if a & b & c;
run;
Registration is now open for SAS Innovate 2025 , our biggest and most exciting global event of the year! Join us in Orlando, FL, May 6-9.
Sign up by Dec. 31 to get the 2024 rate of just $495.
Register now!
Learn how use the CAT functions in SAS to join values from multiple variables into a single value.
Find more tutorials on the SAS Users YouTube channel.
Ready to level-up your skills? Choose your own adventure.