" Exercise 13.2

 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, 06/12/2015 "

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

" Read data from working directory "
FILEREAD [NAME='CABBAGE.DAT'; IMETHOD=read] FGROUPS=no

" Plot data: number of leaves vs sample date "
DGRAPH [KEY=0] Y=NLeaves; X=Days

" Fit SLR to untransformed data "
MODEL Y=NLeaves
FIT [FPROBABILITY=yes; TPROBABILITY=yes] Days

" Plot fitted model "
RGRAPH

" Save fitted values and standardized residuals "
RKEEP FITTEDVALUES=fval; RESIDUALS=res

" Plot residuals against explanatory variate "
DRESIDUALS [RESIDUALS=res; FITTEDVALUES=fval; INDEX=Days] METHOD=index

" Residual plots based on standardized residuals "
DRESIDUALS [RESIDUALS=res; FITTEDVALUES=fval] METHOD=fittedvalues,normal,histogram,absresidual

" Transform response "
CALCULATE logN=LOG(NLeaves)

" Plot data with transformed response "
DGRAPH [KEY=0] Y=logN; X=Days

" Fit SLR on transformed scale "
MODEL Y=logN
FIT [FPROBABILITY=yes; TPROBABILITY=yes] Days

" Plot fitted model "
RGRAPH

" Save fitted values and standardized residuals "
RKEEP FITTEDVALUES=fvalt; RESIDUALS=rest

" Plot residuals against explanatory variate "
DRESIDUALS [RESIDUALS=rest; FITTEDVALUES=fvalt; INDEX=Days] METHOD=index

" Residual plots based on standardized residuals "
DRESIDUALS [RESIDUALS=rest; FITTEDVALUES=fvalt] METHOD=fittedvalues,normal,histogram,absresidual

" Create factor to check for LOF "
GROUPS VECTOR=Days; FACTOR=fDays

" Check for LOF "
MODEL Y=logN
FIT [PRINT=#,accumulated; FPROBABILITY=yes; TPROBABILITY=yes] Days+fDays

" Final model = SLR on transformed scale "
MODEL Y=logN
FIT [FPROBABILITY=yes; TPROBABILITY=yes] Days

" Predict across range of days with SE "
VARIATE [VALUES=0...46] xD
PREDICT [PREDICTIONS=Tpt; SE=Tset] Days; LEVELS=xD
VTABLE Tpt,Tset; VARIATE=Vpt,Vset

" Save residual df "
RKEEP DF=resdf

" Calculate 95% CI "
CALC Lowert = Vpt - EDT(0.975;resdf)*Vset
CALC Uppert = Vpt + EDT(0.975;resdf)*Vset

" Back-transform to original scale "
CALC Bpt,Blowert,Buppert = EXP(Vpt,Lowert,Uppert)

" Plot back-transformed SLR with 95% CI "
PEN 11,12; METHOD=line; SYMBOL=0
DGRAPH Bpt,Blowert,Buppert; xD; PEN=11,12,12; \
   DESC='Back-transformed SLR','Back-transformed 95% CI',''

" Compare to model fitted to mean numbers: intercept = 7.299, slope = 0.261 "
CALC fitmean = 7.299 + 0.261*xD
DGRAPH Bpt,fitmean,NLeaves; xD,xD,Days; PEN=11,12,1; \
   DESC='Back-transformed SLR','SLR fit to mean numbers','Number of leaves'


" Estimate rate of increase on original scale = slope*exp(fitted) "
" Write expression for rate "
EXPRESSION [VALUE=ERate='Days'*EXP('Constant'+'Days'*x)] ratecalc
" Create value for dummy variable x "
CALC x=0
" Calculate rate at x=0 "
RFUNCTION [PRINT=e,se; CALC=ratecalc] ERate

" Now calculate rates for x = 0-46 as a loop with SE using delta method "
VARIATE [NVAL=47] EstRate,SERate
FOR [INDEX=i] xx=0...46
  CALC x=xx
  RFUNCTION [PRINT=*; CALC=ratecalc; SE=se] ERate
  CALC EstRate$[i] = ERate
  CALC SERate$[i] = se
ENDFOR
CALC LRate,URate = EstRate + (-1,1)*EDT(0.975;resdf)*SERate
PRINT xD,LRate,EstRate,URate

" End of File "