* use afb.dta

generate event = 1 - censored

global time afb
global event event
global xlist edu sibs white immig lowinc16 city16


describe $time $event $xlist
summarize $time $event $xlist

* Set data as survival time
stset $time, failure($event)
stdescribe
stsum

* List of hazard and survival functions - reproducing Table 8.2
* in W&B textbook
ltable afb, hazard survival

* Nonparametric estimation

* Graph of survival function (Kaplan-Meier survival curve)
sts graph, survival 

* Graph of Kaplan-Meier Failure Function
sts graph, failure 

* Graph of smoothed hazard function
sts graph, hazard

* Graph of cumulative hazard ratio (Nelson-Aalen cumulative hazard curve)
sts graph, cumhaz

* Kaplan-Meier survival curves for the white groups
sts graph, by(white) 

* Test for equality of survival functions between white groups
* using logrank and Wilcoxon tests
sts test white, logrank
sts test white, wilcoxon

* Kaplan-Meier survival curves for the immig groups
sts graph, by(immig) 

* Test for equality of survival functions between immigrant groups
* using logrank and Wilcoxon tests
sts test immig, logrank
sts test immig, wilcoxon

* Kaplan-Meier survival curves for the loinc16 groups
sts graph, by(lowinc16) 

* Test for equality of survival functions between loinc16 groups
* using logrank and Wilcoxon tests
sts test lowinc16, logrank
sts test lowinc16, wilcoxon

* Kaplan-Meier survival curves for the city16 groups
sts graph, by(city16)

* Test for equality of survival functions between city16 groups
* using logrank and Wilcoxon tests
sts test city16, logrank
sts test city16, wilcoxon

* Parametric models (accelerared failure time (ACF) parametrization)
* Here duplicating the results in Table 8.6 in W&B textbook

* Exponential regression coefficients and fitted hazard curve 
streg $xlist, dist(exponential) time nolog
stcurve, hazard

* Weibull regression coefficients and fitted hazard curve 
streg $xlist, dist(weibull) time nolog
stcurve, hazard

* Lognormal regression coefficients and fitted hazard curve 
streg $xlist, dist(lognormal) nolog
stcurve, hazard

* Lognormal regression (Gamma Fraility) coefficients and fitted hazard curve
* Note: In order to run this streg command you must have run the previous
* lognormal streg command  
streg $xlist, frailty(gamma) time nolog
stcurve, hazard 

* Cox proportional hazard model coefficients and hazard rates
* Notice the nonparametrically fitted hazard curve is possibly a better
* match of the smoothed hazard function that was fitted above by the 
* the nonparametric "sts graph, hazard" command.  
stcox $xlist, nohr
stcurve, hazard
stcox $xlist
* Test of Cox Proportional Hazards Assumption in Cox Regession
estat phtest, log detail
* The PH assumption of the Cox Regression is not supported by the data