* SAS Lab 9 Exercise;

PROC IMPORT out=S3A3.respiratory
   datafile="Data\Respiratory.txt"
       DBMS=DLM
       REPLACE;
   Getnames=yes;
   Datarow=2;
run;

* Part (a);
PROC GLM Data=S3A3.respiratory;
    Class Race Gender;
	Model FEV1=Race Gender Race*Gender Height 
               Height*Race Height*Gender Height*Race*Gender /solution;
run;

* Part (b);
* I will do this entirely in SAS but it is often easier to calculate;
* the test statistic using the printed output on your calculator;

* Full Model;
PROC GLM Data=S3A3.respiratory plots=none;
    Class Race Gender;
	Model FEV1=Race Gender Race*Gender Height 
               Height*Race Height*Gender Height*Race*Gender /solution;
	ODS output OverallAnova=full_anova;
run;

Data full_anova;
   set full_anova;
   If source='Error' then call symput ('full_sse', ss);
   If source='Error' then call symput ('full_df', df); 
run;

* Reduced Model;
PROC GLM Data=S3A3.respiratory;
    Class Race Gender;
	Model FEV1=Race Gender Race*Gender Height /solution;
	ODS output OverallAnova=Reduced_anova;
run;
quit;

Data Reduced_anova;
   set Reduced_anova;
   If source='Error' then call symput ('red_sse', ss);
   If source='Error' then call symput ('red_df', df); 
run;

* Consctruct an output dataset with the statistic and p-value;
Data Ftest;
   num=(&red_sse-&full_sse)/(&red_df-&full_df);
   num_df=&red_df-&full_df;
   den=&full_sse/&full_df;
   den_df=&full_df;
   F=num/den;
   pvalue=1-cdf('F', F, num_df, den_df);
run;

PROC PRINT data=Ftest;
run;

* We see this gives a test statistic F=4.75471 and p-value p=0.00473;

* Part (c);
* Reduced Model which only includes Gender, Height and their interaction;
PROC GLM Data=S3A3.respiratory;
    Class Race Gender;
	Model FEV1=Gender Height Height*Gender /solution;
	ODS output OverallAnova=Reduced_anova1;
run;
quit;

Data Reduced_anova1;
   set Reduced_anova1;
   If source='Error' then call symput ('red_sse', ss);
   If source='Error' then call symput ('red_df', df); 
run;

* Consctruct an output dataset with the statistic and p-value;
Data Ftest1;
   num=(&red_sse-&full_sse)/(&red_df-&full_df);
   num_df=&red_df-&full_df;
   den=&full_sse/&full_df;
   den_df=&full_df;
   F=num/den;
   pvalue=1-cdf('F', F, num_df, den_df);
run;

PROC PRINT data=Ftest1;
run;

* This gives F=5.39924 (df=4,63) and p=0.00084;