" Exercise 17.6

 Statistical Methods in Biology: Design and Analysis of Experiments and Regression
 by S.J. Welham, S.A. Gezan, S.J. Clark & A. Mead (2014)
 Chapman & Hall/CRC Press, Boca Raton, Florida. ISBN: 978-1-4398-0878-8
 Version 1, 11/10/2015 "

" Set working directory - change to location of your data file "
\SET [WORKINGDIRECTORY='d:/stats4biol/data/']

" Read data from working directory "
FILEREAD [NAME='FERTILIZER.DAT'; IMETHOD=read] FGROUPS=no,2(yes,no)

" Plot data " 
DGRAPH Y=Yield; X=N; PEN=Block

" Fit exponential model "
MODEL Y=Yield
FITCURVE [CURVE=exponential; FPROBABILITY=yes] Block+N

" Save standardized residuals and fitted values and produce residual plots "
RKEEP [RMETHOD=DEVIANCE] RESIDUALS=Resexp; FITTEDVALUES=Fitexp
DRESIDUALS [RESIDUALS=Resexp; FITTEDVALUES=Fitexp] METHOD=fittedvalues,normal,histogram,absresidual

" Plot residuals against N "
DGRAPH Y=Resexp; X=N

" Plot fitted curves "
RGRAPH

" Save residual sum of squares and df "
RKEEP DEVIANCE=ResSSexp; DF=ResDFexp

" Create factor version of N "
GROUPS VECTOR=N; FACTOR=fN

" Analyse as designed experiment with 5 groups to get pure error "
MODEL Y=Yield
FIT [PRINT=#,accumulated; FPROBABILITY=yes] Block+fN

" Save pure error SS and df "
RKEEP DEVIANCE=PESS; DF=PEDF

" Check lack of fit (LOF) for exponential model "
CALCULATE LOFSS=ResSSexp-PESS
CALCULATE LOFDF=ResDFexp-PEDF
CALCULATE PEMS=PESS/PEDF
CALCULATE F=(LOFSS/LOFDF)/PEMS
CALCULATE P=CUF(F;LOFDF;PEDF)
PRINT STRUCTURE=F,P,LOFDF,PEDF; DECIMALS=2(*,0)


" Try exponential plus linear "
MODEL Y=Yield
FITCURVE [CURVE=lexponential; FPROBABILITY=yes] Block+N

" Save standardized residuals and fitted values and produce residual plots "
RKEEP [RMETHOD=DEVIANCE] RESIDUALS=Reslpe; FITTEDVALUES=Fitlpe
DRESIDUALS [RESIDUALS=Reslpe; FITTEDVALUES=Fitlpe] METHOD=fittedvalues,normal,histogram,absresidual

" Plot residuals against N "
DGRAPH Y=Reslpe; X=N

" Plot fitted curves "
RGRAPH

" Save residual sum of squares and df "
RKEEP DEVIANCE=ResSSlpe; DF=ResDFlpe

" Check for improvement in fit against exponential model "
CALCULATE ChangeSS=ResSSexp-ResSSlpe
CALCULATE ChangeDF=ResDFexp-ResDFlpe
CALCULATE LpeMS=ResSSlpe/ResDFlpe
CALCULATE F=(ChangeSS/ChangeDF)/LpeMS
CALCULATE P=CUF(F;ResDFlpe;ResDFlpe)
PRINT STRUCTURE=F,P,ChangeDF,ResDFlpe; DECIMALS=2(*,0)

" Re-check lack-of-fit "
CALCULATE F=(ResSSlpe/ResDFlpe)/LpeMS
CALCULATE P=CUF(F;ResDFlpe;ResDFlpe)
PRINT STRUCTURE=F,P,ChangeDF,ResDFlpe; DECIMALS=2(*,0)

" End of File "