*************************************************************************************************************
Statitical Test
Mantel and Byar 1974 JASA
Evaluation of response-time data involving transient status: an illustration using heart-transplant data

Kaplan Meier Curve
Simon and Makuch 1984 Statistics in medicine                                                               
A non-parametric graphical representation of the relationship responder bias

Macro Parameters                                           
      data		: Analysis dataset                                 
      time		: Variable name of duration between enrollment and event (Death = 1 or Censor = 0)              
      wtime		: Variable name of duration between enrollment and status transition, like heart transportation or response to chemotherapy
				  If no transprtation is received, null should be imputted.
      landmark		: Numeric value to set a landmark timepoint (Default = 0)    
      outkm		: Output rtf for Kaplan-Meier curve
SOFTWARE/VERSION	: SAS 9.2
Auther			: Yuki Yoshida (yyoshida@shionogi.com)  
*************************************************************************************************************;
%macro SandM(data=, time=, wtime=, landmark=0, outKM=);
***Total N***;
data _NULL_;
 set  &data end=eof;
  if eof then call symputx('obs',_N_);
run;

data __trans;
 set &data;
 keep id &wtime;
run;

data __original;
 set &data;
 if &wtime = . then group = 0;
 else               group = 1;
 keep id &time Status group;
run;

data __dat001;
 set __trans(in = a rename =(&wtime = &time)) __original;
 if a = 1 then do; status =2; group = 1; end;
run; 

proc means data = __dat001 noprint nway;
 var id;
 class &time status group;
 output out= __dat002(drop = _:) n = n;
run;

proc transpose data = __dat002 out = __dat003 prefix = col;
 var n;
 by &time;
 id status group;
run;

data __dat004;
 set __dat003; retain Den1 &obs Den2 0;
  if col10 = . then col10 = 0;
  if col21 = . then col21 = 0;
  if col11 = . then col11 = 0;
  if col00 = . then col00 = 0;
  if col01 = . then col01 = 0;

  LagEN = lag(Col10);
  LagCN = lag(Col00);
  LagET = lag(Col11);
  LagCT = lag(Col01);

/*  Denominator No Tran*/
  if _n_ = 1 then Den1 = Den1 - col21;
  else            Den1 = Den1 - Col21 - LagEN - LagCN;
/*Denominator Tran*/
  if _n_ = 1 then Den2 = Den2 + col21;
  else            Den2 = Den2 + Col21 - LagET - LagCT;

/*Event ratio;*/
  if den1 ^in(. 0) then S_N_at_timeN = 1- col10/Den1;
  if den2 ^in(. 0) then S_T_at_timeN = 1- col11/Den2;
  if den1  in(. 0) then S_N_at_timeN = 1;
  if den2  in(. 0) then S_T_at_timeN = 1;

/*Label*/
  label col10 = 'N of Event in State1'
        col21 = 'N of State2'
		col11 = 'N of Event in State2'
		col00 = 'N of Cens in State1'
		col01 = 'N of Cens in State2';
run;

/*Mentel & Bayr*/
data __dat004_1;
 set __dat004;
  N_notE1 = den1 - col10;
  N_notE2 = den2 - col11;
  keep time Den1 Den2 col10 col11 N_notE1 N_notE2;
run;

data __dat004_2;
 set __dat004_1 end=eof; retain Event num den 0;  
  eq1 = Den2/(Den1+Den2)*(col11+col10);
  eq2 = (Den1*Den2)/(Den1+Den2)/(Den1+Den2);
   if (Den1 + Den2 -1) ^= 0 then do;
    factor = (Den1+Den2-col11-col10)*(col11+col10)/(Den1 + Den2 -1);
   end;
   else factor = 0;
  eq2F = eq2 * factor;

  Event = event + col11;
  Num   = Num + eq1;
  Den   = Den + eq2F;

  Chisq = (abs(event - num) -0.5)**2/den;
  p     = 1-probchi(chisq,1);
  if eof = 1 then output; 
  keep chisq p;
  label chisq = 'Chisq'
        p     = 'P-value';
  format p pvalue6.4;
run;


/*Landmark*/
data __dat005;
 set __dat004;
  if &time >=&landmark then output;
run;

/*Survival*/
data __dat005_1;
 set __dat005;retain s1 1 s2 1;
 s1 = s1 * S_N_at_timeN;
 s2 = s2 * S_T_at_timeN;
 label s1 = 'State1 Survive'
       s2 = 'State2 Survive'
       den1 = 'N in State1'
       den2 = 'N in State2';
 keep &time den1 col10 col00 col01 col21 den2 col11 s1 s2;
run;

/*figure*/
proc transpose data = __dat005_1 out = __dat006;
 var s1 s2;
 by &time;
run;

/*Start From 1*/
data __temp001;
 set __dat006;
 time = -0.00001;
 COL1 = 1;
 if _n_ <= 2;
run;

data __dat007;
 set __temp001 __dat006;
  label COL1 = 'Survival'
        time = 'Time'
		_LABEL_ = 'Group';
run;

title;
footnote;
option nodate;
ods rtf file = "&outKM" bodytitle nogtitle;

title1 "Mantel and Byar's test";
proc print data = __dat004_2 noobs label;
 var Chisq P;
run;

title1 "Simon and Makuch's representation";
proc sgplot data=__dat007;
 step x=&time y=COL1 / group=_LABEL_;
run;
ods rtf close;
title;
footnote;
%mend;