/*  SAS code for analyzing data from the battery
    life experiment on assignment 7 in 2002.  */

data set1;
  infile 'c:\courses\st511\sas\batteries.dat';
  input material temp y;
  run;


/*  Print the data  */

proc print data=set1;
run;



/* Compute sample means for all factor combinations
   with data.  Make a profile plot. */


proc sort data=set1; by material temp;
proc means data=set1 noprint; by material temp;
  var Y;
  output out=means mean=my;
  run;



goptions device=WIN  target=WINPRTC rotate=landscape;

/* UNIX users can use the following options:

goptions cback=white colors=(black)
    targetdevice=ps300 rotate=landscape;
*/

axis1 label=(f=swiss h=2.0)
      value=(f=swiss h=2.0)  w=3.0  length= 6.0 in;

axis2 label=(f=swiss h=1.8 a=90 r=0)
      value=(f=swiss h=2.0) w= 3.0  length = 4.5 in;

       SYMBOL1 V=circle H=2.0 w=3 l=1 i=join ;
       SYMBOL2 V=diamond H=2.0 w=3 l=3 i=join ;
       SYMBOL3 V=square H=2.0 w=3 l=9 i=join ;

 proc gplot data=means;
   plot my*temp=material / vaxis=axis2 haxis=axis1;
   title  H=3.0 F=swiss "Sample means";
       label my='Mean Response ';
       label material = 'material ';
 run;


/*  Perform analysis of variance where the temp
    factor is entered into the model before the material
    factor.  Use the LSMEANS statement to compare means
    for different combinations of material and temp factors  */

  proc glm data=set1;
    class material temp;
    model y =  temp material temp*material / solution
               ss1 ss2 ss3 e e1 e2 e3 p;
    lsmeans temp*material / pdiff tdiff stderr;
    run;

/*  Perform analysis of variance where the material factor
    is entered into the model before the temp
    factor.  Use the LSMEANS statement to compare means
    for different combinations of material and temp factors  */

  proc glm data=set1;
    class material temp;
    model y = material temp material*temp / solution
               ss1 ss2 ss3 e e1 e2 e3;
    output out=set2 predicted=yhat residual=resid;
    run;


 /* Compute a normal probability plot for
        the resdiuals  */

 proc rank data=set2 normal=blom out=set2;
   var resid; ranks q;  run;


axis1 label=(f=swiss h=2.0)
      value=(f=swiss h=2.0)  w=3.0  length= 6.0 in;

axis2 label=(f=swiss h=1.8 a=90 r=0)
      value=(f=swiss h=2.0) w= 3.0  length = 4.5 in;

       SYMBOL1 V=circle H=2.0 w=3 l=1 i=none ;


 proc gplot data=set2;
   plot resid*q / vaxis=axis2 haxis=axis1;
   title  H=3.0 F=swiss "Normal Probability Plot";
       label resid='Ordered Residuals ';
       label q = 'Standard Normal Quantiles ';
 run;


 /* Check variance homogeneity  */


axis1 label=(f=swiss h=2.0)
      value=(f=swiss h=2.0)  w=3.0  length= 6.0 in;

axis2 label=(f=swiss h=1.8 a=90 r=0)
      value=(f=swiss h=2.0) w= 3.0  length = 4.5 in;

       SYMBOL1 V=circle H=2.0 w=3 l=1 i=none ;


 proc gplot data=set2;
   plot resid*yhat / vaxis=axis2 haxis=axis1;
   title  H=3.0 F=swiss "Residual Plot";
       label resid='Residuals ';
       label q = 'Estimated Means ';
 run;