/* Plot -ln(-ln) survival plots to assess proportional hazard assumption  */
/* 12/17/96  J.M.Garrett                            STB-35 gr2            */
/* Form:  stphplot, by(x) adjust(covlist) options                         */
/* Must have stset data before using this command                         */ 
/* Options required: by(x) or strata(x) with adjust(covlist)              */
/* Options allowed:  adjust, symbol, connect, titles, labels,  noshow     */   

program define stphplot
  version 5.0
  st_is
  #delimit ;
    local options "BY(string) STRata(string) ADJust(string) NOSHow
      L1title(string) L2title(string) *" ;
  #delimit cr 
  local if "opt"
  parse "`*'"
  preserve
  local t: char _dta[st_t]
  local d: char _dta[st_d]
  if "`by'"=="" & "`strata'"=="" {
     disp in red "Either by() or strata() are required"
     exit
     }
  if "`by'"~="" {local x="`by'"}
  if "`strata'"~="" {local x="`strata'"}
  if "`by'"~="" & "`adjust'"=="" {local gensurv=1}
  if "`by'"~="" & "`adjust'"~="" {local gensurv=2}
  if "`strata'"~="" & "`adjust'"~="" {local gensurv=3}
  if "`strata'"~="" & "`adjust'"=="" {
     disp in red "strata() requires adjust(); perhaps you mean by()"
     disp in blue "r(198);"
     exit
     }
  capture keep `if'
  keep `x' `d' `t' `adjust' 
  quietly drop if `x'==.
  quietly tab `x', gen(x)
  local numcat=_result(2)
  local xlbl : variable label `x'
  if "`xlbl'"=="" {local xlbl="`x'"}

* If there are covariates, drop obs. with missing values, center
  parse "`adjust'", parse(" ")
  local i=1
  while "`1'"~=""  {
     local cov`i'="`1'"
     quietly drop if `cov`i''==.
     quietly recast float `cov`i''
     quietly sum `cov`i''
     quietly replace `cov`i''=`cov`i''- _result(3)
     local i=`i'+1
     macro shift
     }

* generate loglog line for each level of class variable
  if `gensurv'==1 {sts gen surv=s, by(`x')}
  if `gensurv'==2 {sts gen surv=s, by(`x') adjust(`adjust')}
  if `gensurv'==3 {sts gen surv=s, strata(`x') adjust(`adjust')}
  tempvar logt
  quietly gen `logt'=log(`t')
  label var `logt' "ln(`t')"
  sort `logt'
  local i=1
  while `i'<=`numcat'  {
    local varlbl`i' : variable label x`i'
    quietly gen surv`i'=-log(-log(surv)) if x`i'==1
    label var surv`i' "`varlbl`i''"
    local i=`i'+1
    }
  drop surv
 
* graph the results
  if "`noshow'"=="" {st_show}
  * if "`title'"=="" {local title "      Test of PH Assumption -- by `x'"}
  if "`l1title'"=="" {local l1title "-Ln[-Ln(Survival Probabilities)]"}
  if "`l2title'"=="" {local l2title "By Categories of `xlbl'"}
  graph surv* `logt', l1("`l1title'") l2("`l2title'") `options'
end