/*************************************************
 SAS EXAMPLE   -- CONTINGENCY TABLES (CROSS-TABS)
                  LOGISTIC REGRESSION
 PROCS USED:
      PROC LOGISTIC
      PROC FREQ
      PROC GENMOD
 FILENAME: logistic.sas
*************************************************/
options yearcutoff=1900;
options pageno=1 title formdlim=" ";

data bcancer;
 infile "d:\510\2006\data\brca.dat" lrecl=300;
 input idnum 1-4 stopmens 5 agestop1 6-7 numpreg1 8-9 agebirth 10-11
       mamfreq4 12 @13 dob mmddyy8. educ 21-22
	   totincom 23 smoker 24 weight1 25-27;

 format dob mmddyy10.;
 if dob = "09SEP99"D then dob=.;
 if stopmens=9 then stopmens=.;
 if agestop1 = 88 or agestop1=99 then agestop1=.;
 if agebirth =99 then agebirth=.;
 if numpreg1=99 then numpreg1=.;
 if mamfreq4=9 then mamfreq4=.;
 if educ=99 then educ=.;
 if totincom=8 or totincom=9 then totincom=.;
 if smoker=9 then smoker=.;
 if weight1=999 then weight1=.;

 if stopmens = 1 then menopause=1;
 if stopmens = 2 then menopause=0;

 yearbirth = year(dob);
 age = int(("01JAN1997"d - dob)/365.25);


 if educ not=. then do;
   if educ in (1,2,3,4) then edcat = 1;
   if educ in (5,6)     then edcat = 2;
   if educ in (7,8)     then edcat = 3;

   highed = (educ in (6,7,8));
 end;

 if age not=. then do;
   if age <50 then agecat=1;
   if age >=50 and age < 60 then agecat=2;
   if age >=60 and age < 70 then agecat=3;
   if age >=70 then agecat=4;

   if age < 50 then over50 = 0;
   if age >=50 then over50 = 1;

   if age >= 50 then highage = 1;
   if age <  50 then highage = 2;
 end;
run;

title "Descriptive Statistics";
proc means data=bcancer n nmiss min max mean std;
run;

title "Logistic Regression with a Continuous Predictor";
proc logistic data=bcancer descending;
  model menopause = age / rsquare;
  units age = 1 5 10;
run;
/************************************************************************
ods rtf file = "d:\510\2006\handouts\SAS\logistic_regression\logistic.rtf";
ods graphics on;
title "Logistic Regression with a Continuous Predictor";
proc logistic data=bcancer descending;
  model menopause = age / rsquare;
  units age = 1 5 10;
  output out=pdat dfbetas= _all_ 
                  difchisq = d_chisq
                  difdev = d_dev
                  reschi = res_chisq
                  resdev = res_dev; 
  graphics estprob;
run;
ods graphics off;
ods rtf close;
title "Check the outlier";
proc print data=pdat;
  where res_chisq not=. and res_chisq < -5;
run;

*************************************************************************/
title "Oneway Frequencies";
proc freq data=bcancer;
  tables dob;
  tables stopmens menopause;
  tables educ edcat;
  tables age agecat over50 highage;
run;

/*Crosstabs of HIGHAGE by STOPMENS*/
title "2 x 2 Table";
title2 "HIGHAGE Coded as 1, 2";
proc freq data=bcancer2;
  tables highage*stopmens / relrisk chisq;
run;

title "Logistic Regression with Dummy Variable Predictor";
title2 "Use Dummy Variable, Coded as 0, 1";
proc logistic data=bcancer2 descending;
  model menopause =  over50/ rsquare;
run;

title "Relationship of Education Categories to Menopause";
proc freq data=bcancer;
  tables edcat*menopause / chisq;
run;

title "Logistic Regression to Predict Menopause From Education";
proc logistic data=bcancer descending;
  class edcat(ref="1") / param = ref;
  model menopause =  edcat/ rsquare;
run;

title "Relationship of AGECAT to MENOPAUSE";
proc freq data=bcancer;
  tables agecat*menopause/ chisq nocol nopercent;
run;

title "Logistic Regression with AGECAT";
title2 "This Analysis Does not Work";
title3 "Check out the Parameter Estimates and Standard Errors";
proc logistic data=bcancer descending;
  class agecat(ref="1") / param = ref;
  model menopause =  agecat/ rsquare;
run;

/*Recode Agecat into AGECAT3 with 3 categories*/
data bcancer2;
  set bcancer;
  if age not=. then do;
    if age < 50 then agecat3 = 1;
	if age >=50 and age < 60 then agecat3 = 2;
	if age >=60 then agecat3 = 3;
  end;
run;

title "Logistic Regression with Ordinal Categorical Predictor";
title2 "This Analysis Works";
proc logistic data=bcancer2 descending;
  class agecat3(ref="1") / param = ref;
  model menopause =  agecat3/ rsquare;
run;

title "Logistic Regression with Several Predictors";
proc logistic data=bcancer descending;
  class edcat(ref="1") / param = ref;
  model menopause = age edcat smoker totincom numpreg1 
                    / rsquare;
run;

title "Logistic Regression Using Proc Genmod";
proc genmod data=bcancer descending;
  class edcat(ref="1") / param = ref;
  model menopause = age edcat smoker totincom numpreg1 
                    / dist=bin type3;
run;