options nocenter ls = 96 ps = 54 mprint mlogic symbolgen;

ods html body = "example4.htm";

title1 "Program for Chapter Four Example";

title2 "UIS Data";

data uisData;
infile "C:\Research\HosmerAndLemeshow\Data\UISDat.txt";
input id age beck ivhx ndrugtx race treat site dfree;
*** Create indicator variables for IV drug use history at admission.;
*** 1 = never, 2 = previous, and 3 = recent.;
if ivhx = 2 then
   ivhx_2 = 1;
else
   ivhx_2 = 0;
if ivhx = 3 then
   ivhx_3 = 1;
else
   ivhx_3 = 0;
*** Create fractional polynomial transformations.;
ndrgfp1 = 1 / (( ndrugtx + 1 ) / 10 );
ndrgfp2 = ndrgfp1 * log (( ndrugtx + 1 ) / 10 );
*** Create some interaction terms.;
age_ndrgfp1 = age * ndrgfp1;
race_site = race * site;
run;

*proc print data = uisData;
*run;

title2 "Check Constructed Variables";

proc freq data = uisData;
table ivhx * ivhx_2 * ivhx_3 / list missing;
run;

title2 "Table 4.1 on Page 105";

%let depVars = age beck ndrugtx ivhx_2 ivhx_3 race treat site;

%put depVars = &depVars.;

proc logistic data = uisData descending;
model dfree = age / plrl waldrl;
units age = 10;
run;

proc logistic data = uisData descending;
model dfree = beck / plrl waldrl;
units beck = 5;
run;

%macro univariate;
%do depVar = 3 %to 8;
   proc logistic data = uisData descending;
   model dfree = %scan ( &depVars., %eval ( &depVar. ));
   run;
%end;
%mend univariate;

%univariate;

title2 "Table 4.2 on Page 106";

proc logistic data = uisData descending;
model dfree = age ndrugtx ivhx_2 ivhx_3 race treat site;
run;

title2 "Table 4.7 on Page 113";

proc logistic data = uisData descending;
model dfree = age ndrgfp1 ndrgfp2 ivhx_2 ivhx_3 race treat site;
run;

title2 "Table 4.9 on Page 115";

proc logistic data = uisData descending;
model dfree = age ndrgfp1 ndrgfp2 ivhx_2 ivhx_3 race treat site age_ndrgfp1 race_site;
output out = probs predicted = phat;
run;

title2 "Table 5.1 on Page 150";

proc sort data = probs;
by phat;
run;

data decile;
set probs;
if _n_ <= 58 then
   decile = 1;
else if _n_ <= 115 then
   decile = 2;
else if _n_ <= 173 then
   decile = 3;
else if _n_ <= 230 then
   decile = 4;
else if _n_ <= 288 then
   decile = 5;
else if _n_ <= 345 then
   decile = 6;
else if _n_ <= 403 then
   decile = 7;
else if _n_ <= 460 then
   decile = 8;
else if _n_ <= 518 then
   decile = 9;
else if _n_ <= 575 then
   decile = 10;
else
   decile = .;
run;

*proc print data = decile;
*var decile phat dfree;
*run;

proc means data = decile max sum;
class decile;
var phat;
run;

proc freq data = decile;
table decile * dfree / nocol norow nopct;
run;

ods html close;