OPTIONS LINESIZE=80;
TITLE1 'Three-way Factorial Design  (all fixed)';
TITLE2 'Factor A (CONCENT) : 3 Levels of Concentration % of Hardwood';
TITLE3 'Factor B (CTIME)   : 2 LEVELS OF Cooking time (3 hrs, 4 hrs)';
TITLE4 'Factor C (PRESS)   : 3 LEVELS OF Pressure (400, 500, 650) ';
DATA PULPS;
/* @@ means a loop in reading input variable sequence */
INPUT CONCENT CTIME PRESS STRENGTH @@;
CARDS;
2 3.0 400 196.6    2 3.0 400 196.0
2 3.0 500 197.7    2 3.0 500 196.0
2 3.0 650 199.8    2 3.0 650 199.4
2 4.0 400 198.4    2 4.0 400 198.6
2 4.0 500 199.6    2 4.0 500 200.4
2 4.0 650 200.6    2 4.0 650 200.9
4 3.0 400 198.5    4 3.0 400 197.2
4 3.0 500 196.0    4 3.0 500 196.9
4 3.0 650 198.4    4 3.0 650 197.6
4 4.0 400 197.5    4 4.0 400 198.1
4 4.0 500 198.7    4 4.0 500 198.0
4 4.0 650 199.6    4 4.0 650 199.0
8 3.0 400 197.5    8 3.0 400 196.6
8 3.0 500 195.6    8 3.0 500 196.2
8 3.0 650 197.4    8 3.0 650 198.1
8 4.0 400 197.6    8 4.0 400 198.4
8 4.0 500 197.0    8 4.0 500 197.8
8 4.0 650 198.5    8 4.0 650 199.8
;
RUN;

/* Print the original data set                   */
PROC PRINT;
RUN;

/* Create Table & tabulate input data set        */
/* Nicety, but not absolutely necessary for      */
/* the analysis - you may skip                   */
PROC TABULATE;
     CLASS CONCENT CTIME PRESS;
     VAR STRENGTH;
     LABEL CONCENT = 'HARDWOOD CONCENTRATION %'
     CTIME = 'COOKING TIME'
     PRESS = 'PRESSURE'
     STRENGTH = 'STRENGTH OF PAPER (2 REPLICATES)';
     TABLE CONCENT, CTIME*PRESS*STRENGTH/CONDENSE;
     KEYLABEL SUM='Y(IJK)';
RUN;

/* General Linear Model (GLM)                    */
PROC GLM;
     CLASS CONCENT CTIME PRESS ;
     MODEL  STRENGTH = CONCENT | CTIME | PRESS ;
     OUTPUT OUT = A1  P = YHAT  R = RESID;
RUN;

PROC PLOT DATA=A1;
     PLOT RESID*YHAT;
RUN;

 
/* Compute Normal Scores of input data, 'A1'     */
/* using cum. normal function/Residual           */
/* BLOM -> yi = þ-1(ri -3/8)/(n+1/4)              */
/*         ri = rank of the ith obs              */
PROC RANK DATA = A1 NORMAL = BLOM OUT = NPLOT1;
     VAR RESID;
     RANKS NSCORE;
RUN;
 
PROC PLOT DATA = NPLOT1;
     PLOT NSCORE*RESID;
RUN;