/* SAS Example:  One observation per Treatment */

/* We use the insurance premium data (given in Chapter 20) */
/* presented in the example in class.                      */

DATA insur;
INPUT premium citysize region;
cards;
  140  1  1
  100  1  2
  210  2  1
  180  2  2
  220  3  1
  200  3  2
;
run;

/* The following is a macro to perform Tukey's test of additivity.  */
/* Just copy it into the SAS editor and invoke it with the name of  */
/* your data set, your factors, and your response, as shown in the  */
/* example below.                                                   */

/******************************************************************/
/******************* MACRO BEGINS HERE ****************************/
%MACRO tukeyaddit(dataset=, factor1=, factor2=, response=);
ods listing close;
proc glm data=&dataset;
  class &factor1 &factor2;
  model &response = &factor1 &factor2;
  ods output overallanova=overall modelanova=model;
run;
quit;
ods listing;
ods output close;
data _null_;
  set overall;
  if source='Corrected Total' then call symput('overall', ss);
run;
data _null_;
  set model ;
  if hypothesistype=1 and source="&factor2" then call symput('ssa', ss);
  if hypothesistype=1 and source="&factor1" then call symput('ssb', ss);
  if hypothesistype=1 and source="&factor2" then call symput('dfa', df);
  if hypothesistype=1 and source="&factor1" then call symput('dfb', df);
run;
%put here is &overall &ssa &ssb &dfa &dfb; /* the statement will appear in the log file so you can check the calculations */
proc sql;
  create table temp1 as
  select &response, &factor1, &factor2 , mean(&response) as yj
  from &dataset
  group by &factor1;
quit;
proc sql;
  create table temp2 as
  select *, mean(&response) as yi
  from temp1
  group by &factor2;
quit;
proc sql noprint; 
   select mean(&response) into :meanp from temp1;
quit;
%put here is  &meanp; /* check value in log file */
proc sql noprint;
  select sum( (yi - &meanp)*(yj - &meanp)*&response ) into :total from temp2;
quit;
%put here is &total; /*check value in log file*/
data final;
  msa = &ssa/(&dfb+1);
  msb = &ssb/(&dfa+1);
  ssab = (&total*&total) / ( msa*msb );
  ssrem = &overall - &ssa - &ssb - ssab;
  f = ssab/( ssrem/((&dfa+1)*(&dfb+1) - (&dfa+1) - (&dfb+1)) );
  p_value = 1- cdf('F',f,  1, (&dfa+1)*(&dfb+1) - (&dfa+1) - (&dfb+1) );
run;
proc print data=final;
run;
%MEND tukeyaddit;
/******************** MACRO ENDS HERE *****************************/
/******************************************************************/

/* Invoking the macro with our values for the dataset, the first factor, */
/* the second factor, and the response:                                  */


%tukeyaddit(dataset = insur, factor1 = citysize, factor2 = region, response = premium)

/* We see that the F* = 6.75 and the P-value for this test is 0.2339.  */
/* Therefore using alpha=0.05, we fail to reject the null.  It is      */
/* reasonable to assume there is no interaction between citysize and   */
/* region, and so we may safely use the no-interaction model.          */


/* Fitting the no-interaction model using PROC GLM: */

PROC GLM data=insur;
CLASS citysize region;
MODEL premium = citysize region;
run;

/*********************************************************************************/

/* Alternate approach to getting Tukey's test for additivity */

PROC GLM DATA=insur NOPRINT;
CLASS citysize region;
MODEL premium = citysize region;
OUTPUT OUT=diag PREDICTED=yhat;
run;

DATA diag; SET diag;
nonadd = yhat**2;
run;

PROC GLM DATA=diag;
CLASS citysize region;
MODEL premium = citysize region nonadd / SS1 SOLUTION;
RUN;

/* The F-test about the 'nonadd' term will be Tukey's test for additivity.*/


/* Fitting the no-interaction model using PROC GLM: */

PROC GLM data=insur;
CLASS citysize region;
MODEL premium = citysize region;
run;