%let location=C:\Documents and Settings\woodph\My Documents\glmwood\kutnersolutions\Chapter7;
PROC IMPORT OUT= WORK.bodyfat 
            DATAFILE= "&location\chapter7.xls" 
            DBMS=EXCEL REPLACE;
     SHEET="CH07TA01$"; 
     GETNAMES=YES;
     MIXED=NO;
     SCANTEXT=YES;
     USEDATE=YES;
     SCANTIME=YES;
RUN;
proc print data=bodyfat;
title "data listing refer to Table 7.1";run;
%include "&location\scatter.sas";
symbol i=rlci95 v=dot;
%scatter(data = bodyfat, var = triceps thigh midarm bodyfat);run;
proc glm;
model bodyfat=triceps;
title "Refer 7.2a, anova summary table";run;
proc glm;
model bodyfat=thigh;
title "Refer 7.2b, anova summary table";run;
proc glm;
model bodyfat=triceps thigh;
title "Refer 7.2c, anova summary table";run;
proc glm;
model bodyfat=triceps thigh midarm;
title "Refer 7.2d, anova summary table";run;
proc glm data=bodyfat outstat=ssvaluesr;
model bodyfat=triceps thigh;
title "proc glm - two predictors used to calculate Type I SS";run;
proc print data=ssvaluesr;run;
data ssvaluesr;set ssvaluesr;if _type_="SS1" or _type_="ERROR";
group="bodyfat=triceps thigh";
proc glm data=bodyfat outstat=ssvaluesl;
model bodyfat=thigh triceps;
title "proc glm - two predictors used to calculate Type I SS";run;
data ssvaluesl;set ssvaluesl;if _type_="SS1" or _type_="ERROR";
group="bodyfat=thigh triceps";
proc print;
title "ss left";run;
data ssvalues;set ssvaluesr ssvaluesl;
proc sort data=ssvalues out=ssvalues;by group;
proc print;by group;var ss _source_;
title "Type I SS for two models Refer to Figure 7.1, page 261";run;
proc gchart;
vbar group/ discrete sumvar=ss subgroup=_source_;
title "Type I SS for two competing models, Refer to Figure 7.1, page 261";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh midarm;
bmidarm: test midarm=0;
title "test of single regression weight via t-test and F- Refer to page 263";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh midarm;
b2b3: test thigh=midarm=0;
title "test of two regression weights via F- Refer to page 263";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh midarm;
b1b2b3: test triceps=thigh=midarm=0;
title "test of all regression weights via F- Refer to page 266 & note correspondence to model F";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh midarm;
eqb2b3: test thigh=midarm;
title "test of equality of two regression weights- refer to equation 7.32 page 268";
proc reg data=bodyfat;
model bodyfat=triceps thigh midarm;
b1_5b2_m3: test triceps=5, thigh=-3;
title "test of multiple regression weights against 'canonical' values Refer to eq. 7.34";
title2 "except that I used different numbers in the bodyfat example";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh midarm/pcorr1 partial stb;
title "coefficients of partial determination. Refer to page 270 & comments section";run;
proc reg data=bodyfat;
model bodyfat=thigh triceps midarm/pcorr1;
title "coefficients of partial determination. Refer to page 270";run;
proc reg data=bodyfat;
model bodyfat=triceps;
output out=y21 r=rbody;
title "computation of partial correlations Refer to p. 271";
proc reg data=y21;
model thigh=triceps;
output out=y21 r=rthigh;
proc corr data=y21;var rbody rthigh;
title "partial correlation first body fat example p. 271";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh;
output out=y3_12 r=rbodyfat;
proc reg data=y3_12;
model midarm=triceps thigh;
output out=y3_12 r=rmidarm;
proc corr data=y3_12;var rbodyfat rmidarm;
title "partial correlation second body fat example p. 271";run;
proc reg data=bodyfat;
model bodyfat=thigh;
output out=y1_2 r=rbodyfat;
proc reg data=y1_2;
model triceps=thigh;
output out=y3_12 r=rtriceps;
proc corr data=y3_12;var rbodyfat rtriceps;
title "partial correlation third body fat example p. 271";run;
PROC IMPORT OUT= WORK.Dwaine 
            DATAFILE= "&location\chapter7.xls" 
            DBMS=EXCEL REPLACE;
     SHEET="CH07TA05$"; 
     GETNAMES=YES;
     MIXED=NO;
     SCANTEXT=YES;
     USEDATE=YES;
     SCANTIME=YES;
RUN;
proc print;
title "listing of data- Dwaine Studios data";run;
data dwaine;set dwaine;
zsales=sales;
zchild=children;
zincome=income;
proc standard mean=1 std=1 data=dwaine out=dwaine;var zsales zchild zincome;
proc print;var zsales zchild zincome;
title "transformed data= refer to Table 7.5, page 277";
proc reg;
model zsales=zchild zincome;
title "Standardized regression via transformed variables Refer to table 7.5, p. 277";run;
proc reg;
model sales=children income/stb;
title "Same analysis based on raw data showing standardized scores as well. Refer to Table 7.5, p. 277";run;
PROC IMPORT OUT= WORK.WorkCrew 
            DATAFILE= "&location\chapter7.xls" 
            DBMS=EXCEL REPLACE;
     SHEET="CH07TA06$"; 
     GETNAMES=YES;
     MIXED=NO;
     SCANTEXT=YES;
     USEDATE=YES;
     SCANTIME=YES;
RUN;
proc print;
title "Work Crew Data- Refer to Table 7.6, page 279";run;
%scatter(data = workcrew, var = product crewsize bonus);run;
proc corr data=workcrew;
title "correlations between predictors and criterion Refer to discussion page 279-280";run;
proc reg;
model product=crewsize bonus/selection=rsquare;
title "Work Crew Data- All possible regressions Refer to Table 7.7, page 280";
title2 "(dividing regression SS by total SS in Table gives R-square)";run;
proc reg;
model product=crewsize bonus/selection=maxr;
title "Work Crew Data- Regressions Refer to Table 7.7, page 280";run;
proc reg;
model product=bonus;
title "Work Crew Data- Regressions Refer to Table 7.7, page 280";run;
PROC IMPORT OUT= WORK.preddata 
            DATAFILE= "&location\chapter7.xls" 
            DBMS=EXCEL REPLACE;
     SHEET="CH07TA08$"; 
     GETNAMES=YES;
     MIXED=NO;
     SCANTEXT=YES;
     USEDATE=YES;
     SCANTIME=YES;
RUN;
proc print;
title "Perfect Corr- Refer to Table 7.8, page 281";run;

proc iml;
use preddata;
read all into x;
a=x`*x;
deta=det(a);
echelon=echelon(a);
*Notice that the following statement gives an error message in the log:
Singular (no relation to the phone company) is another way of saying the
determinant is zero;
inverse=inv(a);
print "a=" a "determinant of a is " deta;
print "echelon of a is " echelon;
title "Note rank of score matrix as determined by echelon is only 2";
proc g3d data=preddata;
scatter x1*x2=y/grid shape="balloon";
title "perfect predictor correlation Refer to Figure 7.2, p. 282";run;
data plane;
do x=0 to 20 by 0.25;
      do y=0 to 20 by 0.25;
         z=-87+x+18*y;
		 z2=-7+9*x+2*y;
         output;
      end;
   end;
proc g3d;
plot y*x=z;
title "first regression surface- refer to Figure 7.2, page 282";run;
proc g3d;
plot y*x=z2;
title "second regression surface- refer to Figure 7.2, page 282";run;
proc reg data=bodyfat;
model bodyfat=triceps;
title "Refer page 285 r-square for triceps alone";run;
proc reg data=bodyfat;
model bodyfat=thigh triceps/pcorr1;
title "Refer page 285 r-square for triceps adjusted for thigh";run;
data example;
input y x1 x2;
cards;
20 5 25
20 10 30
0 5 5
1 10 10
;
proc corr;
title "example collinearity data- refer to page 285";run;
proc reg data=example;
model y=x2;
title "SS regression based on X2- refer to page 286";run;
proc glm data=example;
model y=x1 x2;
title "SS regression of X2 after controllingn for X1 refer to page 286";run;
proc reg data=bodyfat;
model bodyfat=triceps/seb;
title "Standard Error for b1 when only triceps is in model: Refer page 286";run;
proc reg data=bodyfat;
model bodyfat= thigh /seb;
title "Standard Error for b2 when only thigh is in model: Refer page 286";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh /seb;
title "Standard Error for b1 and b2 when triceps and thigh are in model: Refer page 286";run;
proc reg data=bodyfat;
model bodyfat=triceps thigh midarm/seb;
title "Standard Error for b1 and b2 when triceps, thigh, and midarm are in model: Refer page 286";run;
data newobs;
triceps=25;thigh=50;midarm=29;
data bodyfat2;set bodyfat newobs;
proc reg data=bodyfat2;
model bodyfat=triceps;
output out=values stdp=stderr p=pbodyfat;
title "standard error of predicted value: Refer page 286";run;
data values;set values;
if triceps=25;
proc print;
title "standard error of predicted value- Refer page 286";run;
proc reg data=bodyfat2;
model bodyfat=triceps thigh;
output out=values stdp=stderr p=pbodyfat;
title "standard error of predicted value: Refer page 286";run;
data values;set values;
if triceps=25;
proc print;
title "standard error of predicted value: Refer page 287";run;
proc reg data=bodyfat2;
model bodyfat=triceps thigh midarm;
output out=values stdp=stderr p=pbodyfat;
title "standard error of predicted value: Refer page 287";run;
data values;set values;
if triceps=25;
proc print;
title "standard error of predicted value: Refer page 287 (does not exactly match published value of .621)";run;