/******************************************************** THIS COMMAND FILE READS IN THE TECUMSEH DATA SET REARRANGES THE DATA FOR LONGITUDINAL ANALYSIS AND ILLUSTRATES HOW TO ANALYZE THE CIG DATA FOR ALL 3 TIMES. FILENAME: TECUMSEH_GEE.SAS **************************************************************/ /*Lab Ex 10: Read in SAS data set rearrange for GEE analysis*/ libname glmwkshp v9 "c:\temp\glmwkshp\"; data tecumseh; set glmwkshp.tecumseh; bmi1 = wtkg1/(htcm1/100)**2; bmi2 = wtkg2/(htcm2/100)**2; bmi3 = wtkg3/(htcm3/100)**2; if 0 <= bmi1 < 18.5 then bmicat1 = 1; if 18.5<= bmi1 < 25.0 then bmicat1 = 2; if 25.0<= bmi1 < 30.0 then bmicat1 = 3; if bmi1 >=30 then bmicat1 = 4; if 0 <= bmi2 < 18.5 then bmicat2 = 1; if 18.5<= bmi2 < 25.0 then bmicat2 = 2; if 25.0<= bmi2 < 30.0 then bmicat2 = 3; if bmi2 >=30 then bmicat2 = 4; if 0 <= bmi3 < 18.5 then bmicat3 = 1; if 18.5 <= bmi3 < 25.0 then bmicat3 = 2; if 25.0 <= bmi3 < 30.0 then bmicat3 = 3; if bmi3 >=30 then bmicat3 = 4; l_bmi1 =log(bmi1); l_bmi2 =log(bmi2); l_bmi3 =log(bmi3); if cigday1 = 0 then cigcat1 = 0; if cigday1 in (1,2) then cigcat1 = 1; if cigday1 = 3 then cigcat1 = 2; if cigday1 = 4 then cigcat1 = 3; if cigday1 in (5,6,7,8) then cigcat1 = 4; if cigday2 = 0 then cigcat2 = 0; if cigday2 in (1,2) then cigcat2 = 1; if cigday2 = 3 then cigcat2 = 2; if cigday2 = 4 then cigcat2 = 3; if cigday2 in (5,6,7,8) then cigcat2 = 4; if cigday3 = 0 then cigcat3 = 0; if cigday3 in (1,2) then cigcat3 = 1; if cigday3 = 3 then cigcat3 = 2; if cigday3 = 4 then cigcat3 = 3; if cigday3 in (5,6,7,8) then cigcat3 = 4; if marital1 = 1 then marit1 = 1; if marital1 in (2,4,5) then marit1 = 2; if marital1 = 3 then marit1 = 3; if marital2 = 1 then marit2 = 1; if marital2 in (2,4,5) then marit2 = 2; if marital2 = 3 then marit2 = 3; if marital3 = 1 then marit3 = 1; if marital3 in (2,4,5) then marit3 = 2; if marital3 = 3 then marit3 = 3; if sex=1 then male=1; if sex=2 then male=0; run; data teclong; set tecumseh; age=age1; agegrp=agegrp1; marital=marit1; ed=ed1; cig=cig1; cigcat=cigcat1; beer=beer1; bald=bald1; bmi=bmi1; l_bmi=l_bmi1; bmicat = bmicat1; round=1; output; /*Note: there is no variable for Beer at round 2*/ age=age2; agegrp=agegrp2; marital=marit2; ed=ed2; cig=cig2; cigcat=cigcat2; bald=bald2; bmi=bmi2; l_bmi=l_bmi2; bmicat = bmicat2; round=2; output; age=age3; agegrp=agegrp3; marital=marit3; ed=ed3; cig=cig3; cigcat=cigcat3; beer=beer3; bald=bald3; bmi=bmi3; l_bmi=l_bmi3; bmicat = bmicat3; round=3; output; keep id age age1 male age agegrp marital ed cig cigcat beer bald bmi1 bmi l_bmi bmicat round; run; /*Descriptive Statistics for Each Round*/ options nolabel; proc means data=teclong n mean std ; class round; title "Descriptive Statistics for Each Round"; run; proc freq data=teclong; tables ed round; run; /*Lab Ex 9: Run Proc Genmod for binomially dist dependent variable*/ proc genmod data=teclong descending; where age1 between 20 and 60; class id ed round; model cig = male age ed round / dist=binomial link=logit type3; repeated subject=id / type=un within=round corrw; title "GEE Model for Binomial Distribution: Logistic Link"; run; /*GEE Analysis for Poisson distributed dependent variable*/ proc genmod data=teclong ; where age1 between 20 and 60 and round in (1,3); class id ed round; model beer = male age ed round / dist=poisson type3; repeated subject=id / type=un within=round corrw; title "GEE Model for Poisson Distribution: Log Link"; run; /************************************************************* Special Topics: 1. GLIMMIX macro. This approach approximates the non-linear function. It gives us a subject-specific approach. ***************************************************************/ /*First, open the glimmix macro and submit it, or simply define to SAS where it is located using a filename statement*/ data tec_glmm; set teclong; if nmiss(beer,ed,male,age)=0 and 20 <= age1 <= 60; keep id round ed beer male age cig; run; proc means data=tec_glmm; run; filename glmm "c:\temp\glmwkshp\glimmix.sas"; %include glmm; %glimmix(data=tec_glmm, procopt=noclprint, stmts=%str( where round in (1,3); class ed round id; model beer = male age ed round / solution; random int / subject=id; lsmeans ed round;), error=poisson, link=log); %glimmix(data=tec_glmm, procopt=noclprint, stmts=%str( where round in (1,2,3); class ed round id; model cig = male age ed round / solution; random int / subject=id; lsmeans ed round;), error=poisson, link=log); /************************************************************* Special Topics: 2. Proc NLMIXED. This approach uses a numerical integration method to approximate the random effects problem for a non-linear function. ***************************************************************/ data teclong2; set teclong; /*Create Dummy variables for Categorical Predictors*/ if ed not=. then do; lt_hs = (ed=1); hs = (ed=2); end; round1 = (round=1); round2 = (round=2); round3 = (round=3); run; proc freq data=teclong2; tables lt_hs hs round1 round2 round3; run; proc nlmixed data=teclong noad qpoints=3; where age1 between 20 and 60 and nmiss(cig,ed,male,age,lt_hs,hs)=0; parms b0 -.1672 b1 .8737 b2 -.0166 b3 .2340 b4 .0637 b5 .3092 b6 .2062 ; eta = b0 + b1*male + b2*age + b3*lt_hs + b4*hs + b5*round1 +b6*round2 +u; exp_eta=exp(eta); p = 1/(1+exp(-eta)); model cig ~ binary(p); random u ~ normal(0,sigma**2) subject=id; predict p out=p; title "NLMIXED for Binomial Regression"; run; proc nlmixed data=teclong2 noad qpoints=3; where age1 between 20 and 60 and nmiss(beer,ed,marital,male,age)=0 and round in (1,3); parms b0 .0075 b1 0.7534 b2 -.0092 b3 .4848 b4 .2734 b5 -.8410; eta = b0 + b1*male + b2*age + b3*lt_hs + b4*hs + b5*round1 + u; lambda=exp(eta); model beer ~ poisson(lambda); random u ~ normal(0,sigma**2) subject=id; title "NLMIXED for Poisson Regression"; run;