/* Read the bone density data from a text file */ data set1; infile 'c:\courses\st511\snew\bone2.density.txt'; input id group bmd weight age caffine calcium growmilk teenmilk twenmilk; run; /* Examine the contents of the file */ proc contents data=set1; run; /* Fit a model */ proc glm data=set1; class group; model bmd = age group / solution ss1 ss3; output out=set1 residual=r predicted=p; lsmeans group /stderr pdiff; run; /* Compute plotting positions for a normal probability plot */ proc rank data=set1 normal=blom out=set1; var r; ranks q; run; /* Windows users can use the following */ goptions device=WIN target=WINPRTC rotate=landscape; axis1 label=(f=swiss h=1.8 a=90 r=0 "Predicted BMD (g/cm^2)") value=(f=swiss h=1.8) w=3.0 length= 6.5 in; axis2 label=(f=swiss h=2.0 "Residuals") value=(f=swiss h=1.8) w= 3.0 length = 5.0 in; axis1 label=(f=swiss h=1.8 a=90 r=0 "q") order -3 to 3 by 1 value=(f=swiss h=1.8) w=3.0 length= 6.5 in; SYMBOL1 V=CIRCLE H=1.7 w=3 l=1 i=none ; SYMBOL2 V=DIAMOND H=1.7 w=3 l=3 i=join ; SYMBOL3 V=square H=1.7 w=3 l=9 i=join ; PROC GPLOT DATA=set1; PLOT r*p / vaxis=axis2 haxis=axis1; TITLE1 ls=0.01in H=2.0 F=swiss "Residuals vs. Prdeicted Values"; footnote ls=0.01in; RUN; PROC GPLOT DATA=set1; PLOT r*q / vaxis=axis2 haxis=axis1; TITLE1 ls=0.01in H=2.0 F=swiss "Normal Probability Plot"; footnote ls=0.01in; RUN;