statgeek/surv_analysis.sas

## surv_analysis.sas
/*
Author F. Khurshed
Date: June 15, 2011
Purpose: To automatically create a summary table with info and KM Graph

Output created is an RTF file with table and graph below. THe first table is for the overall curve and then for each categorical
variable specified.

*/


%macro surv_results(dsetin=, varlist=, surv_time=, censor=, censor_val=, maxt=, atrisk=atrisk, time_list=, fnout=, dsetout=);
/*
dsetin=dataset to be analyzed
varlist=list of variables for the analysis to be done
surv_time=variable that contains the length of survival
censor=censoring variable
censor_val=values that indicate censoring, currently handles only 1, but may need to handle more at some point
atrisk=atrisk or blank. If atrisk then show at risk numbers in plot, if not then don't show.
fnout=name of output file.rtf
dsetout=output dataset with summary survival values at each timepoint
*/


*Delete old datasets just in case theyre hanging around;
proc datasets nodetails nolist;
delete &dsetout;
run; quit;

ods graphics on;


%let i=1;
%do %while (%scan(&varlist, &i, " ")  ne %str());
%let var=%scan(&varlist, &i, " ");

*Put which iteration for log;
%put &i. &var.;


ods select none;
proc lifetest data=&dsetin timelist=(&time_list) outsurv=ple reduceout;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
	strata &var;
ods output quartiles=quart_v censoredsummary=csumm homtests=pvalue1 (where=(Test='Log-Rank'));
run;


*Format data in num (%) format;
data ple;
	length p1 $30.;
	set ple;
	p1=cat(compress(put(survival, percent7.1))," (", compress(put(sdf_lcl, percent7.1)), ", ", compress(put(sdf_ucl, percent7.1)), ")");
	format timelist 8.;
run;

proc sort data=ple; by &var; run;

*Flip data to structure I want;
proc transpose data=ple out=ple2 (drop= _name_) prefix=time_;
	by &var ;
	id timelist;
	var p1;
run;

*Add in median data;
data median;
	set quart_v;
	if percent ne 50 then delete;
	median=cats(compress(put(estimate, 5.1))," (", compress(put(lowerlimit, 5.1)), ", ", compress(put(upperlimit, 5.1)), ")");
	if estimate=. then median='';
run;

proc sort data=median;
	by &var;
run;

proc transpose data=median out=median2 (drop=_name_) prefix=median;
	by &var;
	var median;
run;

proc sort data=csumm;
where not missing(stratum);
by &var;
run;


*Merge data together;

data result_&var;
	length Variable $50. Category $50.;
	merge csumm (drop=control_var) ple2 median2;
	by &var;
	if _n_=1 then set pvalue1 (keep=ProbChiSq);
	drop stratum;
	Variable="&var.";
	Category=&var.;

	drop &var;
run;

*Append to dsetout for output table;

proc append base=&dsetout data=result_&var;
run;


*Clean up files so no mistakes in next loop;
proc datasets;
	delete ple quart_v csumm sig_test ple2 median median2;
run;
quit;

*Increment counter;
%let i=%eval(&i+1);
%end;

*Get results for Overall Study;
proc lifetest data=&dsetin plots=(survival) timelist=(&time_list) maxtime=&maxt outsurv=ple reduceout;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
	ods output quartiles=quart_v censoredsummary=csumm;
run;


*Format data in num (%) format;
data ple;
	length p1 $30.;
	set ple;
	p1=cat(put(survival, percent7.1)," (", put(sdf_lcl, percent7.1), ",", put(sdf_ucl, percent7.1), ")");
	format timelist 8.;
run;

*Flip data to structure I want;
proc transpose data=ple out=ple2 (drop= _name_) prefix=time_;
	id timelist;
	var p1;
run;

*Add in median data;
data median;
	set quart_v;
	if percent ne 50 then delete;
	median=cats(put(estimate, 5.1)," (", put(lowerlimit, 5.1), ",", put(upperlimit, 5.1), ")");
	if estimate=. then median='';
run;


proc transpose data=median out=median2 (drop=_name_) prefix=median;
	var median;
run;

*Merge data together;

data result_overall;
	merge csumm ple2 median2;
	drop stratum;
	Variable='Overall';
	Category='Overall';
run;

proc append data=result_overall base=&dsetout;
run;


*Clean up files so no mistakes in next loop;
proc datasets;
	delete ple quart_v csumm sig_test ple2 median median2;
run;
quit;


/*Start output of results*/
options nonumber nodate;
ods listing close;
ods rtf file="&fnout" startpage=no;

title 'Overall Survival';
ods select all;
proc print data=result_overall noobs label;
	label median1='Median (95% CI)';
run;

ods rtf startpage=no;

ods select survivalplot;
ods noproctitle;

proc lifetest data=&dsetin plots=(survival(&atrisk)) timelist=(&time_list) maxtime=&maxt;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
run;

title;
ods rtf startpage=now;

%let i=1;
%do %while (%scan(&varlist, &i, " ") ^=%str());
%let var=%scan(&varlist, &i, " ");

*Put which iteration for log;
%put &i. &var.;

*Need to get label to output in legend;
data _null_;
	set &dsetin (obs=1);
	call symput('legend_label', vlabel(&var));
run;


*Display results;
proc print data=result_&var noobs label;
	label median1='Median (95% CI)';
run;

ods rtf startpage=no;

ods select survivalplot;
ods noproctitle;
proc lifetest data=&dsetin plots=(survival(&atrisk test)) maxtime=&maxt;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
	strata &var;
run;

ods rtf startpage=now;

*Increment counter;
%let i=%eval(&i+1);
%end;

ods rtf close;
ods listing;


ods graphics off;

%mend surv_results;


data VALung;
   drop check m;
   retain Therapy Cell;
   infile cards column=column;
   length Check $ 1;
   label SurvTime='Failure or Censoring Time'
      Kps='Karnofsky Index'
      DiagTime='Months till Randomization'
      Age='Age in Years'
      Prior='Prior Treatment?'
      Cell='Cell Type'
      Therapy='Type of Treatment'
      Treatment='Treatment Indicator';
   M=Column;
   input Check $ @@;
   if M>Column then M=1;
   if Check='s'|Check='t' then input @M Therapy $ Cell $ ;
   else input @M SurvTime Kps DiagTime Age Prior @@;
   if SurvTime > .;
   censor=(SurvTime<0);
   SurvTime=abs(SurvTime);
   Treatment=(Therapy='test');
   datalines;
standard squamous
 72 60  7 69  0   411 70  5 64 10   228 60  3 38  0   126 60  9 63 10
118 70 11 65 10    10 20  5 49  0    82 40 10 69 10   110 80 29 68  0
314 50 18 43  0  -100 70  6 70  0    42 60  4 81  0     8 40 58 63 10
144 30  4 63  0   -25 80  9 52 10    11 70 11 48 10
standard small
 30 60  3 61  0   384 60  9 42  0     4 40  2 35  0    54 80  4 63 10
 13 60  4 56  0  -123 40  3 55  0   -97 60  5 67  0   153 60 14 63 10
 59 30  2 65  0   117 80  3 46  0    16 30  4 53 10   151 50 12 69  0
 22 60  4 68  0    56 80 12 43 10    21 40  2 55 10    18 20 15 42  0
139 80  2 64  0    20 30  5 65  0    31 75  3 65  0    52 70  2 55  0
287 60 25 66 10    18 30  4 60  0    51 60  1 67  0   122 80 28 53  0
 27 60  8 62  0    54 70  1 67  0     7 50  7 72  0    63 50 11 48  0
392 40  4 68  0    10 40 23 67 10
standard adeno
  8 20 19 61 10    92 70 10 60  0    35 40  6 62  0   117 80  2 38  0
132 80  5 50  0    12 50  4 63 10   162 80  5 64  0     3 30  3 43  0
 95 80  4 34  0
standard large
177 50 16 66 10   162 80  5 62  0   216 50 15 52  0   553 70  2 47  0
278 60 12 63  0    12 40 12 68 10   260 80  5 45  0   200 80 12 41 10
156 70  2 66  0  -182 90  2 62  0   143 90  8 60  0   105 80 11 66  0
103 80  5 38  0   250 70  8 53 10   100 60 13 37 10
test squamous
999 90 12 54 10   112 80  6 60  0   -87 80  3 48  0  -231 50  8 52 10
242 50  1 70  0   991 70  7 50 10   111 70  3 62  0     1 20 21 65 10
587 60  3 58  0   389 90  2 62  0    33 30  6 64  0    25 20 36 63  0
357 70 13 58  0   467 90  2 64  0   201 80 28 52 10     1 50  7 35  0
 30 70 11 63  0    44 60 13 70 10   283 90  2 51  0    15 50 13 40 10
test small
 25 30  2 69  0  -103 70 22 36 10    21 20  4 71  0    13 30  2 62  0
 87 60  2 60  0     2 40 36 44 10    20 30  9 54 10     7 20 11 66  0
 24 60  8 49  0    99 70  3 72  0     8 80  2 68  0    99 85  4 62  0
 61 70  2 71  0    25 70  2 70  0    95 70  1 61  0    80 50 17 71  0
 51 30 87 59 10    29 40  8 67  0
test adeno
 24 40  2 60  0    18 40  5 69 10   -83 99  3 57  0    31 80  3 39  0
 51 60  5 62  0    90 60 22 50 10    52 60  3 43  0    73 60  3 70  0
  8 50  5 66  0    36 70  8 61  0    48 10  4 81  0     7 40  4 58  0
140 70  3 63  0   186 90  3 60  0    84 80  4 62 10    19 50 10 42  0
 45 40  3 69  0    80 40  4 63  0
test large
 52 60  4 45  0   164 70 15 68 10    19 30  4 39 10    53 60 12 66  0
 15 30  5 63  0    43 60 11 49 10   340 80 10 64 10   133 75  1 65  0
111 60  5 64  0   231 70 18 67 10   378 80  4 65  0    49 30  3 37  0
;
run;


%let temp=cell prior;

%surv_results(dsetin=valung, varlist=&temp, surv_time=survtime, censor=censor, censor_val=0, maxt=1000, atrisk=, time_list=200 400,  fnout=C:\_localData\survival.rtf, dsetout=want);
	/*
	Author F. Khurshed
	Date: June 15, 2011
	Purpose: To automatically create a summary table with info and KM Graph

	Output created is an RTF file with table and graph below. THe first table is for the overall curve and then for each categorical
	variable specified.

	*/



	%macro surv_results(dsetin=, varlist=, surv_time=, censor=, censor_val=, maxt=, atrisk=atrisk, time_list=, fnout=, dsetout=);
	/*
	dsetin=dataset to be analyzed
	varlist=list of variables for the analysis to be done
	surv_time=variable that contains the length of survival
	censor=censoring variable
	censor_val=values that indicate censoring, currently handles only 1, but may need to handle more at some point
	atrisk=atrisk or blank. If atrisk then show at risk numbers in plot, if not then don't show.
	fnout=name of output file.rtf
	dsetout=output dataset with summary survival values at each timepoint
	*/


	*Delete old datasets just in case theyre hanging around;
	proc datasets nodetails nolist;
	delete &dsetout;
	run; quit;

	ods graphics on;




	%let i=1;
	%do %while (%scan(&varlist, &i, " ") ne %str());
	%let var=%scan(&varlist, &i, " ");

	*Put which iteration for log;
	%put &i. &var.;







	ods select none;
	proc lifetest data=&dsetin timelist=(&time_list) outsurv=ple reduceout;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
	strata &var;
	ods output quartiles=quart_v censoredsummary=csumm homtests=pvalue1 (where=(Test='Log-Rank'));
	run;




	*Format data in num (%) format;
	data ple;
	length p1 $30.;
	set ple;
	p1=cat(compress(put(survival, percent7.1))," (", compress(put(sdf_lcl, percent7.1)), ", ", compress(put(sdf_ucl, percent7.1)), ")");
	format timelist 8.;
	run;

	proc sort data=ple; by &var; run;

	*Flip data to structure I want;
	proc transpose data=ple out=ple2 (drop= _name_) prefix=time_;
	by &var ;
	id timelist;
	var p1;
	run;

	*Add in median data;
	data median;
	set quart_v;
	if percent ne 50 then delete;
	median=cats(compress(put(estimate, 5.1))," (", compress(put(lowerlimit, 5.1)), ", ", compress(put(upperlimit, 5.1)), ")");
	if estimate=. then median='';
	run;

	proc sort data=median;
	by &var;
	run;

	proc transpose data=median out=median2 (drop=_name_) prefix=median;
	by &var;
	var median;
	run;

	proc sort data=csumm;
	where not missing(stratum);
	by &var;
	run;



	*Merge data together;

	data result_&var;
	length Variable $50. Category $50.;
	merge csumm (drop=control_var) ple2 median2;
	by &var;
	if _n_=1 then set pvalue1 (keep=ProbChiSq);
	drop stratum;
	Variable="&var.";
	Category=&var.;

	drop &var;
	run;

	*Append to dsetout for output table;

	proc append base=&dsetout data=result_&var;
	run;



	*Clean up files so no mistakes in next loop;
	proc datasets;
	delete ple quart_v csumm sig_test ple2 median median2;
	run;
	quit;

	*Increment counter;
	%let i=%eval(&i+1);
	%end;

	*Get results for Overall Study;
	proc lifetest data=&dsetin plots=(survival) timelist=(&time_list) maxtime=&maxt outsurv=ple reduceout;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
	ods output quartiles=quart_v censoredsummary=csumm;
	run;


	*Format data in num (%) format;
	data ple;
	length p1 $30.;
	set ple;
	p1=cat(put(survival, percent7.1)," (", put(sdf_lcl, percent7.1), ",", put(sdf_ucl, percent7.1), ")");
	format timelist 8.;
	run;

	*Flip data to structure I want;
	proc transpose data=ple out=ple2 (drop= _name_) prefix=time_;
	id timelist;
	var p1;
	run;

	*Add in median data;
	data median;
	set quart_v;
	if percent ne 50 then delete;
	median=cats(put(estimate, 5.1)," (", put(lowerlimit, 5.1), ",", put(upperlimit, 5.1), ")");
	if estimate=. then median='';
	run;


	proc transpose data=median out=median2 (drop=_name_) prefix=median;
	var median;
	run;

	*Merge data together;

	data result_overall;
	merge csumm ple2 median2;
	drop stratum;
	Variable='Overall';
	Category='Overall';
	run;

	proc append data=result_overall base=&dsetout;
	run;


	*Clean up files so no mistakes in next loop;
	proc datasets;
	delete ple quart_v csumm sig_test ple2 median median2;
	run;
	quit;


	/Start output of results/
	options nonumber nodate;
	ods listing close;
	ods rtf file="&fnout" startpage=no;

	title 'Overall Survival';
	ods select all;
	proc print data=result_overall noobs label;
	label median1='Median (95% CI)';
	run;

	ods rtf startpage=no;

	ods select survivalplot;
	ods noproctitle;

	proc lifetest data=&dsetin plots=(survival(&atrisk)) timelist=(&time_list) maxtime=&maxt;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
	run;

	title;
	ods rtf startpage=now;

	%let i=1;
	%do %while (%scan(&varlist, &i, " ") ^=%str());
	%let var=%scan(&varlist, &i, " ");

	*Put which iteration for log;
	%put &i. &var.;

	*Need to get label to output in legend;
	data _null_;
	set &dsetin (obs=1);
	call symput('legend_label', vlabel(&var));
	run;


	*Display results;
	proc print data=result_&var noobs label;
	label median1='Median (95% CI)';
	run;

	ods rtf startpage=no;

	ods select survivalplot;
	ods noproctitle;
	proc lifetest data=&dsetin plots=(survival(&atrisk test)) maxtime=&maxt;
	time &surv_time*&censor(&censor_val) ;
	label &surv_time='Time (Months)';
	strata &var;
	run;

	ods rtf startpage=now;

	*Increment counter;
	%let i=%eval(&i+1);
	%end;

	ods rtf close;
	ods listing;


	ods graphics off;

	%mend surv_results;


	data VALung;
	drop check m;
	retain Therapy Cell;
	infile cards column=column;
	length Check $ 1;
	label SurvTime='Failure or Censoring Time'
	Kps='Karnofsky Index'
	DiagTime='Months till Randomization'
	Age='Age in Years'
	Prior='Prior Treatment?'
	Cell='Cell Type'
	Therapy='Type of Treatment'
	Treatment='Treatment Indicator';
	M=Column;
	input Check $ @@;
	if M>Column then M=1;
	if Check='s'\|Check='t' then input @M Therapy $ Cell $ ;
	else input @M SurvTime Kps DiagTime Age Prior @@;
	if SurvTime > .;
	censor=(SurvTime<0);
	SurvTime=abs(SurvTime);
	Treatment=(Therapy='test');
	datalines;
	standard squamous
	72 60 7 69 0 411 70 5 64 10 228 60 3 38 0 126 60 9 63 10
	118 70 11 65 10 10 20 5 49 0 82 40 10 69 10 110 80 29 68 0
	314 50 18 43 0 -100 70 6 70 0 42 60 4 81 0 8 40 58 63 10
	144 30 4 63 0 -25 80 9 52 10 11 70 11 48 10
	standard small
	30 60 3 61 0 384 60 9 42 0 4 40 2 35 0 54 80 4 63 10
	13 60 4 56 0 -123 40 3 55 0 -97 60 5 67 0 153 60 14 63 10
	59 30 2 65 0 117 80 3 46 0 16 30 4 53 10 151 50 12 69 0
	22 60 4 68 0 56 80 12 43 10 21 40 2 55 10 18 20 15 42 0
	139 80 2 64 0 20 30 5 65 0 31 75 3 65 0 52 70 2 55 0
	287 60 25 66 10 18 30 4 60 0 51 60 1 67 0 122 80 28 53 0
	27 60 8 62 0 54 70 1 67 0 7 50 7 72 0 63 50 11 48 0
	392 40 4 68 0 10 40 23 67 10
	standard adeno
	8 20 19 61 10 92 70 10 60 0 35 40 6 62 0 117 80 2 38 0
	132 80 5 50 0 12 50 4 63 10 162 80 5 64 0 3 30 3 43 0
	95 80 4 34 0
	standard large
	177 50 16 66 10 162 80 5 62 0 216 50 15 52 0 553 70 2 47 0
	278 60 12 63 0 12 40 12 68 10 260 80 5 45 0 200 80 12 41 10
	156 70 2 66 0 -182 90 2 62 0 143 90 8 60 0 105 80 11 66 0
	103 80 5 38 0 250 70 8 53 10 100 60 13 37 10
	test squamous
	999 90 12 54 10 112 80 6 60 0 -87 80 3 48 0 -231 50 8 52 10
	242 50 1 70 0 991 70 7 50 10 111 70 3 62 0 1 20 21 65 10
	587 60 3 58 0 389 90 2 62 0 33 30 6 64 0 25 20 36 63 0
	357 70 13 58 0 467 90 2 64 0 201 80 28 52 10 1 50 7 35 0
	30 70 11 63 0 44 60 13 70 10 283 90 2 51 0 15 50 13 40 10
	test small
	25 30 2 69 0 -103 70 22 36 10 21 20 4 71 0 13 30 2 62 0
	87 60 2 60 0 2 40 36 44 10 20 30 9 54 10 7 20 11 66 0
	24 60 8 49 0 99 70 3 72 0 8 80 2 68 0 99 85 4 62 0
	61 70 2 71 0 25 70 2 70 0 95 70 1 61 0 80 50 17 71 0
	51 30 87 59 10 29 40 8 67 0
	test adeno
	24 40 2 60 0 18 40 5 69 10 -83 99 3 57 0 31 80 3 39 0
	51 60 5 62 0 90 60 22 50 10 52 60 3 43 0 73 60 3 70 0
	8 50 5 66 0 36 70 8 61 0 48 10 4 81 0 7 40 4 58 0
	140 70 3 63 0 186 90 3 60 0 84 80 4 62 10 19 50 10 42 0
	45 40 3 69 0 80 40 4 63 0
	test large
	52 60 4 45 0 164 70 15 68 10 19 30 4 39 10 53 60 12 66 0
	15 30 5 63 0 43 60 11 49 10 340 80 10 64 10 133 75 1 65 0
	111 60 5 64 0 231 70 18 67 10 378 80 4 65 0 49 30 3 37 0
	;
	run;



	%let temp=cell prior;

	%surv_results(dsetin=valung, varlist=&temp, surv_time=survtime, censor=censor, censor_val=0, maxt=1000, atrisk=, time_list=200 400, fnout=C:\_localData\survival.rtf, dsetout=want);