proc format;
label drink 1='0 drinks'
            2='Less than 1 drink'
            3='1 to 2 drinks'
            4='3 to 5 drinks'
            5='6 or more drinks';
data malformation; *See page 89;
input consumption score present absent @@;
total=present+absent;
format consumption drink.;
datalines;
1 0 48 17066 2 0.5 38 14464 3 1.5 5 788 4 4.0 1 126 5 7.0 1 37
;
run;
proc logistic data=malformation;
class consumption/param=ref;
model present/total=consumption;
run;
*Changing the reference level with a keyword;
proc logistic data=malformation;
class consumption/param=ref ref=first;
model present/total=consumption;
run;
*Changing the reference level to option other than keyword choices;
proc logistic data=malformation;
class consumption (ref='3')/param=ref;
model present/total=consumption;
run;
*Testing contrasts;
proc logistic data=malformation; 
class consumption / param=ref ref=first; 
model present/total = consumption; 
contrast "beta3-beta1" consumption -1 0 1 0;
contrast "exp(beta3-beta1)" consumption -1 0 1 0 / estimate=exp;
contrast "beta3-beta1" consumption -1 0 1 0 / estimate;
run;
*Tests of trend;
proc genmod data=malformation; 
model present/total = consumption / dist=bin link=logit; run;
proc genmod data=malformation; 
model present/total = consumption / dist=bin link=identity; run;
proc genmod data=malformation; 
model present/total = score / dist=bin link=logit; run;
proc genmod data=malformation; 
model present/total = score / dist=bin link=identity; run;

*Cochran-Armitage trend test based on scores--an alternate approach;
proc sort data=malformation; by consumption score;
proc transpose data=malformation out=maltable (rename=(col1=Count _name_=Malformation)); by consumption score;
var present absent;
run;
proc print data=maltable; run;
proc freq data=maltable;
weight count;
tables malformation*score/chisq trend;
*chisq will provide us with the M^2 test statistic;
*trend provides the Cochran-Armitage trend test;
run;