" 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 "