R version 2.8.0 (2008-10-20) Copyright (C) 2008 The R Foundation for Statistical Computing ISBN 3-900051-07-0 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. Natural language support but running in an English locale R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > x <- c(361.58,363.19,363.61,364.14,365.51,365.51,365.5,365.5,364.59,364.63,364.54,363.67,365.22,369.05,370.45,370.46,370.46,370.58,370.58,370.22,370.21,370.29,370.29,370.2,370.2,372.55,374.51,375.58,375.75,375.75,375.75,375.69,375.76,377.5,377.51,377.74,369.82,373.1,374.55,375.01,374.81,375.31,375.31,375.39,375.59,376.26,377.18,377.26,377.26,381.87,387.09,387.14,388.78,389.16,389.16,389.42,389.49,388.97,388.97,389.09,389.09,391.76,390.96,391.76,392.8,393.06,393.06,393.26,393.87,394.47,394.57,394.57,394.57,399.57,406.13,407.03,409.46,409.9,409.9,410.14,410.54,410.69,410.79,410.97) > par3 = 'multiplicative' > par2 = 'Triple' > par1 = '12' > #'GNU S' R Code compiled by R2WASP v. 1.0.44 () > #Author: Prof. Dr. P. Wessa > #To cite this work: Wessa P., (2010), Exponential Smoothing (v1.0.4) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_exponentialsmoothing.wasp/ > #Source of accompanying publication: > #Technical description: > par1 <- as.numeric(par1) > if (par2 == 'Single') K <- 1 > if (par2 == 'Double') K <- 2 > if (par2 == 'Triple') K <- par1 > nx <- length(x) > nxmK <- nx - K > x <- ts(x, frequency = par1) > if (par2 == 'Single') fit <- HoltWinters(x, gamma=F, beta=F) > if (par2 == 'Double') fit <- HoltWinters(x, gamma=F) > if (par2 == 'Triple') fit <- HoltWinters(x, seasonal=par3) > fit Holt-Winters exponential smoothing with trend and multiplicative seasonal component. Call: HoltWinters(x = x, seasonal = par3) Smoothing parameters: alpha: 0.8830462 beta : 0.02329537 gamma: 1 Coefficients: [,1] a 413.9777769 b 0.7537160 s1 0.9914268 s2 1.0021214 s3 1.0075816 s4 1.0060531 s5 1.0053703 s6 1.0032777 s7 1.0021041 s8 1.0013627 s9 0.9992616 s10 0.9973972 s11 0.9955049 s12 0.9927344 > myresid <- x - fit$fitted[,'xhat'] > postscript(file="/var/www/html/freestat/rcomp/tmp/1uuht1294340136.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow=c(2,1)) > plot(fit,ylab='Observed (black) / Fitted (red)',main='Interpolation Fit of Exponential Smoothing') > plot(myresid,ylab='Residuals',main='Interpolation Prediction Errors') > par(op) > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/2ctpb1294340136.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > p <- predict(fit, par1, prediction.interval=TRUE) > np <- length(p[,1]) > plot(fit,p,ylab='Observed (black) / Fitted (red)',main='Extrapolation Fit of Exponential Smoothing') > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/3s8rh1294340136.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow = c(2,2)) > acf(as.numeric(myresid),lag.max = nx/2,main='Residual ACF') > spectrum(myresid,main='Residals Periodogram') > cpgram(myresid,main='Residal Cumulative Periodogram') > qqnorm(myresid,main='Residual Normal QQ Plot') > qqline(myresid) > par(op) > dev.off() null device 1 > > #Note: the /var/www/html/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/html/freestat/rcomp/createtable") > > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Estimated Parameters of Exponential Smoothing',2,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Parameter',header=TRUE) > a<-table.element(a,'Value',header=TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'alpha',header=TRUE) > a<-table.element(a,fit$alpha) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'beta',header=TRUE) > a<-table.element(a,fit$beta) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'gamma',header=TRUE) > a<-table.element(a,fit$gamma) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/www/html/freestat/rcomp/tmp/4a8k71294340136.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Interpolation Forecasts of Exponential Smoothing',4,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'t',header=TRUE) > a<-table.element(a,'Observed',header=TRUE) > a<-table.element(a,'Fitted',header=TRUE) > a<-table.element(a,'Residuals',header=TRUE) > a<-table.row.end(a) > for (i in 1:nxmK) { + a<-table.row.start(a) + a<-table.element(a,i+K,header=TRUE) + a<-table.element(a,x[i+K]) + a<-table.element(a,fit$fitted[i,'xhat']) + a<-table.element(a,myresid[i]) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/freestat/rcomp/tmp/5faq51294340136.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Extrapolation Forecasts of Exponential Smoothing',4,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'t',header=TRUE) > a<-table.element(a,'Forecast',header=TRUE) > a<-table.element(a,'95% Lower Bound',header=TRUE) > a<-table.element(a,'95% Upper Bound',header=TRUE) > a<-table.row.end(a) > for (i in 1:np) { + a<-table.row.start(a) + a<-table.element(a,nx+i,header=TRUE) + a<-table.element(a,p[i,'fit']) + a<-table.element(a,p[i,'lwr']) + a<-table.element(a,p[i,'upr']) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/freestat/rcomp/tmp/6n6er1294340136.tab") > > try(system("convert tmp/1uuht1294340136.ps tmp/1uuht1294340136.png",intern=TRUE)) character(0) > try(system("convert tmp/2ctpb1294340136.ps tmp/2ctpb1294340136.png",intern=TRUE)) character(0) > try(system("convert tmp/3s8rh1294340136.ps tmp/3s8rh1294340136.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 1.232 0.671 1.386