data display;
input height $ width $ sale @@ ;
cards;
Bottom Regular 47
Bottom Regular 43
Bottom Wide 46
Bottom Wide 40
Middle Regular 62
Middle Regular 68
Middle Wide 67
Middle Wide 71
Top Regular 41
Top Regular 39
Top Wide 42
Top Wide 46
;

proc glm data=display;
   class height width;
   model sale = height width height*width;
  lsmeans height / stderr pdiff adjust=tukey;
   lsmeans width / stderr pdiff adjust=tukey;  
   contrast "regular vs wide" width -1 1;
   contrast "middle vs the others" height 0.5 -1 0.5; 
   run;

/* IF there is no significant interaction effect */

proc glm data=display;
   class height width;
   model sale = height width ;
   lsmeans height / stderr pdiff adjust=tukey;
   lsmeans width / stderr pdiff adjust=tukey;  
   contrast "regular vs wide" width -1 1;
   contrast "middle vs the others" height 0.5 -1 0.5; 
   run;

/* Data Detergent -- Two factors */
data detergent;
input brand $ temp $ dirt ;
cards;
Super Cold 4
Super Cold 5
Super Cold 6
Super Cold 5
Super Warm 7
Super Warm 9
Super Warm 8
Super Warm 12
Super Hot 10
Super Hot 12
Super Hot 11
Super Hot 9
Best Cold 6
Best Cold 6
Best Cold 4
Best Cold 4
Best Warm 13
Best Warm 15
Best Warm 12
Best Warm 12
Best Hot 12
Best Hot 13
Best Hot 10
Best Hot 13
;


/* run full model with interaction effect */

proc glm data = detergent;
   class brand temp;
   model dirt = brand temp brand*temp; 
   run;

/* IF the interaction effect is significant */
proc glm data=detergent;
   class brand temp;
   model dirt = brand temp brand*temp;
   lsmeans  brand*temp / pdiff stderr adjust=TUKEY;
   means brand*temp;
   run;
 
/* Interaction Plot - Method 2.*/
proc glm data=detergent;
   class brand temp;
   model dirt = brand | temp;         /*  full model including all interactions*/      
   lsmeans brand * temp ;
   contrast "Cold vs Hot&Warm by Best vs Super" brand*temp 1 -0.5 -0.5 
                                                          -1  0.5  0.5;      
   output out=diag p=py r=ry;         /* p=predicted value, r=residual */
   run;

/* residual plot */ 
proc plot data=diag;                         /* Diagnostic Plots */
   plot dirt*py;                            /* observed v. predicted */     
   plot py*ry;  
   run;

/* model checking*/
proc univariate data=diag normal;
   var ry;
   qqplot;
   run;