use "lecture11.dta", clear

quietly xi: logistic dead i.agecat male hrsdm hrsca hrslung hrschf bmi smoke bath money walkjog push if valid==1

gen points = push+2*walkjog + 2*money+2*bath + 2*smoke+bmi+2*hrschf+2*hrslung+2*hrsca+hrsdm+2*male + 7*_Iagecat_6 + 5*_Iagecat_5 + 4*_Iagecat_4 + 3*_Iagecat_3+2*_Iagecat_2 + _Iagecat_1

* slide 2
quietly logistic dead points if valid==0
lroc if valid==0, scheme(s1color)

* slide 3
quietly logistic dead points if valid==1
lroc if valid==1, scheme(s1color)

* slide 4
twoway (kdensity points if dead==1 & valid==1, bw(2)) (kdensity points if dead==1 & valid==0, bw(2)) (kdensity points if dead==0  & valid==1,bw(1.3)) (kdensity points if dead==0  & valid==0,bw(1.3)), xtitle("Linear Predictor") ytitle(" ") legend(order(1 "Died, Valid" 2 "Died, Training" 3 "Alive Valid" 4 "Alive Training")) scheme(s1color)

* slide 6
xi: logistic dead i.agecat male hrsdm hrsca hrslung hrschf bmi smoke bath money walkjog push if valid==1

quietly xi: logistic dead i.agecat male hrsdm hrsca hrslung hrschf bmi smoke bath money walkjog push if valid==0
quietly predict prob

* slide 8
table points if valid==0, c(n dead mean dead mean prob)

* slide 9
table points if valid==1, c(n dead mean dead mean prob)

* slide 10
twoway (lowess dead prob if valid==0) (lowess dead prob if valid==1) (lfit prob prob),  ytitle(Observed Probability of Dying) xtitle(Predicted Probability by Model) scheme(s1color) legend(order(1 "Learning Set" 2 "Validation Set" ))