/* Polynomial Regression, Test for Lack of Fit */
options ls=75 nodate;

data coffee;
input disp sales @@;
disp2 = disp*disp; lof=disp;

cards;
0 508.1 0 498.4 1 568.2 1 577.3 2 651.7 2 657.0 4 755.3 4 758.9
5 787.6 5 792.1 6 841.4 6 831.8 7 854.7 7 871.4
;

proc plot; plot sales*disp;

/*fit a quadratic model.  Note that PROC REG requires that all covariates
be defined in the data step.  PROC GLM, on the other hand, allows us to
have powers in the model statement.*/

proc reg; 
   title1 "note SSR from linear model on df=1 (type I SS for disp)";
   title2 "note SSR from quadratic model on df=2";
   model sales = disp disp2/ss1 ss2;
run;



*lack of fit test via proc glm;

proc glm; 
   title "full model, note SSTrt on df=6";
   class disp; 
   model sales = disp;
run;


proc glm; 
model sales = disp disp*disp disp*disp*disp
disp*disp*disp*disp disp*disp*disp*disp*disp
disp*disp*disp*disp*disp*disp;
run;

/* Note that PROC GLM permits higher order terms  in the polynomial
using the * operator and by not declaring disp as a class variable.
PROC REG insists that higher order terms be available as variables,
for example as disp2 disp3 etc */

proc glm; 
   title "note ss[lof] of quadratic model on df=6-2=4";
class lof; model sales = disp disp*disp lof/ss1 solution;
proc glm; 
   title "note ss[lof] of linear model on df=6-1=5";
   class lof; 
   model sales = disp lof/ss1 solution;
run;

proc gplot;
  plot sales*disp/haxis=-1 to 8;
  symbol1 v=dot c=black i=rq;  *i=rq sets method of interpolation to quadratic regression;
run;
proc gplot;
  plot sales*disp/haxis=-1 to 8;
  symbol1 v=dot c=black i=rl;  *i=rl sets method of interpolation to linear regression;
run;

proc insight; run;  

/* Use SAS/INSIGHT to fit polynomial regressions that are
1) linear
2) quadratic
3) of degree 6 (equivalent to one-way ANOVA)
Compare the regression SS with that achievvd in PROC REG and PROC GLM 
*/