/* code to plot observations against block means */
/*  assumes that the block variable is called BLOCK */
/*    the treatments are identified by TRT */
/*    and the response is called Y */

/* in general stuff in CAPITOL letters may need updating */

options formdlim = '-' nodate;

data raw;
  infile 'multblock.txt';
  input block trt y;
  logy =  log(y);
run;
        
proc sort;
  by BLOCK;

proc print;
  title 'Original data';
run;
             
proc means noprint;
  by BLOCK;
  var Y;
  output out = means mean = blockmean;
run;
     
proc print;
  title 'Data set with block averages';
run;
     
data new;
  merge raw means;
  by BLOCK;

proc print;
  title 'Merged data set';
run;
            
proc plot;
  plot Y*blockmean = TRT;
  title 'plot of observations against block  means';
run;
            
proc glm data = raw;
  class block trt;
  model y = block trt;
  output out = resids r = resid p = yhat;

proc plot;
  plot resid*yhat;
  title 'Residual plot';
run;

proc plot data =  new;
  plot logY*blockmean = TRT;
  title 'plot of log observations against block means';
run;
            
proc glm data = new;
  class block trt;
  model logy = block trt;
  output out = resids r = resid p = yhat;

proc plot;
  plot resid*yhat;
  title 'Residual plot for log transformed obs';
run;