/* SAS example:  Repeated Measures Analysis */
/* The data set is the wine judging data from chapter 27 */


/* Each of 6 subjects (judges) rates 4 wines in a random sequence. */

/* The rating is the response.         */
/* The wine types are the treatments.  */

/* Main question:  Is there a significant difference */
/* in the mean ratings for the 4 wines?              */

DATA winedata;
INPUT rating judge wine;
cards;
  20  1  1
  24  1  2
  28  1  3
  28  1  4
  15  2  1
  18  2  2
  23  2  3
  24  2  4
  18  3  1
  19  3  2
  24  3  3
  23  3  4
  26  4  1
  26  4  2
  30  4  3
  30  4  4
  22  5  1
  24  5  2
  28  5  3
  26  5  4
  19  6  1
  21  6  2
  27  6  3
  25  6  4
;
run;

/* Informally checking the "compound symmetry" assumption: */

/* Getting the within-subjects variance-covariance matrix: */

PROC SORT data=winedata;
BY judge;
PROC TRANSPOSE data=winedata out = newdata PREFIX = winenum;
BY judge;
ID wine;
VAR rating;
run;

PROC CORR DATA=newdata COV NOPROB NOCORR NOSIMPLE;
VAR winenum1-winenum4;  /* These identify the four wines in SAS */
run;

/* The covariance matrix is printed. */

/* We check whether the variances (the diagonal elements) */
/* are roughly of similar magnitude. */

/* We also check whether all the covariances (the off-diagonal elements) */
/* are roughly of similar magnitude. */

/* If so, the compound symmetry assumption is reasonable. */

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

/* The analysis is identical to that of the RCBD: */

PROC GLM DATA=winedata;
  CLASS judge wine;
  MODEL rating = judge wine;
  RANDOM judge;
  MEANS wine / TUKEY CLDIFF alpha=0.05;
  OUTPUT out=diagnost p=yhat r=resid;
run;

/* All the same assumption checks that were used in the RCBD can be used here. */

/* Normal Q-Q plot: */
PROC UNIVARIATE noprint data=diagnost;
  QQPLOT resid / normal;
run;

/* Interaction plot (like Figure 27.2 in book): */


PROC SGPLOT data=diagnost;
  SERIES X=wine Y=rating / GROUP=judge;
RUN;

/******************* 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 = winedata, factor1 = judge, factor2 = wine, response = rating)

quit;