options pageno=1;
data threefct;
input rep a b c yield;
datalines;
1 0 0 0 25.7
2 0 0 0 25.4
3 0 0 0 23.8
4 0 0 0 22.0
1 0 0 1 31.8
2 0 0 1 29.5
3 0 0 1 28.7
4 0 0 1 26.4
1 0 0 2 34.6
2 0 0 2 37.2
3 0 0 2 29.1
4 0 0 2 23.7
1 0 1 0 27.7
2 0 1 0 30.3
3 0 1 0 30.2
4 0 1 0 33.2
1 0 1 1 38.0 
2 0 1 1 40.6
3 0 1 1 34.6
4 0 1 1 31.0
1 0 1 2 42.1
2 0 1 2 43.6
3 0 1 2 44.6
4 0 1 2 42.7
1 1 0 0 28.9
2 1 0 0 24.7
3 1 0 0 27.8
4 1 0 0 23.4
1 1 0 1 37.5
2 1 0 1 31.5
3 1 0 1 31.0
4 1 0 1 27.8
1 1 0 2 38.4
2 1 0 2 32.5
3 1 0 2 31.2
4 1 0 2 29.8
1 1 1 0 38.0
2 1 1 0 31.0
3 1 1 0 29.5
4 1 1 0 30.7
1 1 1 1 36.9
2 1 1 1 31.9
3 1 1 1 31.5
4 1 1 1 35.9
1 1 1 2 44.2
2 1 1 2 41.6
3 1 1 2 38.9
4 1 1 2 37.6
1 2 0 0 23.4
2 2 0 0 24.2
3 2 0 0 21.2
4 2 0 0 20.9
1 2 0 1 25.3
2 2 0 1 27.7
3 2 0 1 23.7
4 2 0 1 24.3
1 2 0 2 29.8
2 2 0 2 29.9
3 2 0 2 24.3
4 2 0 2 23.8
1 2 1 0 20.8
2 2 1 0 23.0
3 2 1 0 25.2
4 2 1 0 23.1
1 2 1 1 29.0
2 2 1 1 32.0
3 2 1 1 26.5
4 2 1 1 31.2
1 2 1 2 36.6
2 2 1 2 37.8
3 2 1 2 34.8
4 2 1 2 40.2
;;
proc print;
title 'Printout of 3-factor Data';
run;
*COMMENT:  ANOVA for A whole plot, BXC factorial arrangement for subplot.  
All factors are fixed effects;
proc anova;
classes rep a b c;
model yield=rep a a(rep) b c a*b a*c b*c a*b*c;
*COMMENT:  a(rep) is Error(a);
test h=a e=a(rep);
means a/lsd e=a(rep);
means b c/lsd;
means a*b a*c b*c a*b*c;
title1 'ANOVA for A whole plot, BxC factorial arrangement for
subplot';
title2 'All factors fixed';
run;
*COMMENT:  ANOVA for AxB factorial arrangement for whole plot, C subplot.
  All factors are fixed effects;
proc anova;
classes rep a b c;
model yield=rep a b a*b a*b(rep) c a*c b*c a*b*c;
*COMMENT:  a*b(rep) is Error(a);
test h=a b a*b e=a*b(rep);
means a b/lsd e=a*b(rep);
means c/lsd;
means a*b a*c b*c a*b*c;
title1 'ANOVA for AxB factorial arrangement for whole plot, C
suplot.';
title2 'All factors fixed';
run;