* Data from Lindsey (1995) from Andersen (1977); * This includes both sas GRAPH and SGPLOT for graphics: For SGPLOT it is same basic commands that are in the document on graphing in R and SAS; data lcancer; input age $ 1-5 age_midpt city $ cases population; lpop = log(population); rate = cases/population; lograte = log(rate); age_sq = age_midpt*age_midpt; frederic=0; if city='Frederic' then frederic=1; datalines; 40-54 47 Fredericia 11 3059 55-59 57 Fredericia 11 800 60-64 62 Fredericia 11 710 65-69 67 Fredericia 10 581 70-74 72 Fredericia 11 509 >75 75 Fredericia 10 605 40-54 47 Horsens 13 2879 55-59 57 Horsens 6 1083 60-64 62 Horsens 15 923 65-69 67 Horsens 10 834 70-74 72 Horsens 12 634 >75 75 Horsens 2 782 40-54 47 Kolding 4 3142 55-59 57 Kolding 8 1050 60-64 62 Kolding 7 895 65-69 67 Kolding 11 702 70-74 72 Kolding 9 535 >75 75 Kolding 12 659 40-54 47 Vejle 5 2520 55-59 57 Vejle 7 878 60-64 62 Vejle 10 839 65-69 67 Vejle 14 631 70-74 72 Vejle 8 539 >75 75 Vejle 7 619 ; run; goptions reset=(axis, legend, pattern, symbol, title, footnote) norotate hpos=0 vpos=0 htext=2.5 ftext=swiss ctext= target= gaccess= gsfmode= ; goptions device=WIN ctext=black graphrc interpol=join; axis1 color=black width=2.0 label=('Age (mid-point of ranges)') order=45 to 80 by 5 ; axis2 color=black width=2.0 label=(angle=90 'Rate = cases/population') order=0 to .025 by .005; axis3 color=black width=2.0 label=(angle=90 'Log(Rate)'); legend1 position=(top outside center) frame cshadow=pink ; symbol1 color=blue value=dot i=join; symbol2 color=green value=dot i=join; symbol3 color=red value=dot i=join; symbol4 color=orange value=dot i=join; proc gplot data=WORK.LCANCER; plot rate * age_midpt = city / legend=legend1 haxis=axis1 vaxis=axis2 frame ; run; proc gplot data=WORK.LCANCER; plot lograte * age_midpt = city / legend=legend1 haxis=axis1 vaxis=axis3 frame ; run; ods graphics/imagefmt=jpeg; proc sgplot; title 'How to do this with sg'; scatter x=age_midpt y=rate / group=city; series x=age_midpt y=rate / group=city; run; proc genmod data=lcancer order=data; class city age; model cases = city age / link=log dist=poisson offset=lpop type3; title1 'Poission loglinear Model for Rates'; title2 'cases = city age'; run; /*** New code ***/ ods graphics/imagefmt=jpeg; proc sgplot; title 'How to do this with sg'; scatter x=age_midpt y=rate / group=city; series x=age_midpt y=rate / group=city; run; /*** End new code ***/ proc genmod data=lcancer order=data; class city ; model cases = city age_midpt / link=log dist=poisson offset=lpop type3; output out=model2 pred=fitted upper=up lower=lo; title1 'Poission loglinear Model for Rates'; title2 'cases = city age_midpt'; run; /*** New Code ****/ data fix; set model2; prate = fitted/population; run; ods graphics/imagefmt=jpeg; proc sgplot data=fix; title 'How to do this with sg'; scatter x=age_midpt y=rate / group=city; series x=age_midpt y=prate / group=city; run; /*** End new code ***/ data tmp; set model2; FittedRate = fitted/population; run; proc print; run; data plotmod2; input age Frate Ffit Hrate Hfit Krate Kfit Vrate Vfit; datalines; 47 0.003596 0.004919 0.004515 0.003625 0.001273 0.003466 0.001984 0.003847 57 0.013750 0.008678 0.005540 0.006395 0.007619 0.006113 0.007973 0.006785 62 0.015493 0.011526 0.016251 0.008493 0.007821 0.008120 0.011919 0.009012 67 0.017212 0.015308 0.011990 0.011281 0.015670 0.010784 0.022187 0.011969 72 0.021611 0.020331 0.018927 0.014982 0.016822 0.014323 0.014842 0.015897 75 0.016529 0.024105 0.002558 0.017764 0.018209 0.016982 0.011309 0.018848 run; goptions reset=(axis, legend, pattern, symbol, title, footnote) norotate hpos=0 vpos=0 htext=2.5 ftext=swiss ctext= target= gaccess= gsfmode= ; goptions device=WIN ctext=black graphrc interpol=join; axis1 color=black width=2.0 label=('Age (mid-point of ranges)') order=45 to 80 by 5 ; axis2 color=black width=2.0 label=(angle=90 'Rate = cases/population') order=0 to .025 by .005; axis3 color=black width=2.0 label=(angle=90 'Log(Rate)'); legend1 position=(top outside center) frame cshadow=pink ; legend1 position=(top outside center) frame cshadow=pink label=none value=('Frederica' 'Horsens' 'Kolding' 'Vejle' ); symbol1 color=blue value=dot i=none; symbol2 color=red value=dot i=none; symbol3 color=green value=dot i=none; symbol4 color=orange value=dot i=none; symbol5 color=blue value=none i=join; symbol6 color=red value=none i=join; symbol7 color=green value=none i=join; symbol8 color=orange value=none i=join; proc gplot data=plotmod2; plot (Frate Hrate Krate Vrate Ffit Hfit Kfit Vfit)* age / overlay legend=legend1 haxis=axis1 vaxis=axis2 frame ; run; proc genmod data=lcancer order=data; class city ; model cases = city age_midpt age_sq / link=log dist=poisson offset=lpop type3; output out=model3 pred=fitted upper=up lower=lo; title1 'Poission loglinear Model for Rates (quadratic) '; title2 'cases = city age_midpt age_midpt^2'; run; data tmp; set model3; FittedRate = fitted/population; run; proc print; var age_midpt rate fittedRate; run; data model3plot; input age Frate Ffit Hrate Hfit Krate Kfit Vrate Vfit; datalines; 47 0.003596 0.003433 0.004515 0.002461 0.001273 0.002358 0.001984 0.002610 57 0.013750 0.011949 0.005540 0.008567 0.007619 0.008208 0.007973 0.009085 62 0.015493 0.016889 0.016251 0.012109 0.007821 0.011601 0.011919 0.012841 67 0.017212 0.019838 0.011990 0.014224 0.015670 0.013627 0.022187 0.015083 72 0.021611 0.019365 0.018927 0.013885 0.016822 0.013303 0.014842 0.014724 75 0.016529 0.017465 0.002558 0.012522 0.018209 0.011997 0.011309 0.013279 run; proc gplot data=model3plot; plot (Frate Hrate Krate Vrate Ffit Hfit Kfit Vfit)* age / overlay legend=legend1 haxis=axis1 vaxis=axis2 frame ; title ; run; proc genmod data=lcancer order=data; class frederic ; model cases = frederic age_midpt age_sq / link=log dist=poisson offset=lpop type3 obstats; title1 'Poission loglinear Model for Rates (quadratic) '; title2 'cases = frederic age_midpt age_sq '; output out=finalmodel pred=fitted upper=up lower=lo; run; proc print data=finalmodel; run; data toplot; input age_midpt Frate Fpop Ffit Hrate Hpop Hfit Krate Kpop Kfit Vrate Vpop Vfit; Ffitr=Ffit/Fpop; Hfitr=Hfit/Hpop; Kfitr=Kfit/Kpop; Vfitr=Vfit/Vpop; datalines; 47 0.003596 3059 10.4877 0.004515 2879 7.1068 0.001273 3142 7.7561 0.001984 2520 6.2206 57 0.013750 800 9.5555 0.005540 1083 9.3138 0.007619 1050 9.0300 0.007973 878 7.5508 62 0.015493 710 11.9910 0.016251 923 11.2236 0.007821 895 10.8832 0.011919 839 10.2022 67 0.017212 581 11.5295 0.011990 834 11.9162 0.015670 702 10.0302 0.022187 631 9.0157 72 0.021611 509 9.8627 0.018927 634 8.8451 0.016822 535 7.4639 0.014842 539 7.5197 75 0.016529 605 10.5737 0.002558 782 9.8404 0.018209 659 8.2926 0.011309 619 7.7893 run; goptions reset=(axis, legend, pattern, symbol, title, footnote) norotate hpos=0 vpos=0 htext=2.5 ftext=swiss ctext= target= gaccess= gsfmode= ; goptions device=WIN ctext=black graphrc interpol=join; axis1 color=black width=2.0 label=('Age (mid-point of ranges)') ; axis2 color=black width=2.0 label=(angle=90 'Rate = cases/population') ; axis3 color=black width=2.0 label=(angle=90 'Log(Rate)'); legend1 position=(top outside center) frame cshadow=pink label=none value=('Frederica' 'Horsens' 'Kolding' 'Vejle' 'Frederica Fitted' 'Other cities fitted'); symbol1 color=blue value=dot i=none; symbol2 color=red value=dot i=none; symbol3 color=green value=dot i=none; symbol4 color=orange value=dot i=none; symbol5 color=blue value=none i=join; symbol6 color=black value=none i=join; proc gplot data=toplot; plot (Frate Hrate Krate Vrate Ffitr Vfitr)* age_midpt / overlay legend=legend1 haxis=axis1 vaxis=axis2 frame ; run;