R version 2.12.2 (2011-02-25) Copyright (C) 2011 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: x86_64-redhat-linux-gnu (64-bit) 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(101397,97994,94044,91159,87239,89235,118647,125620,125154,117529,109459,108483,107137,104699,100804,96066,91971,93228,120144,127233,127166,118194,109940,106683,102834,99882,96666,92540,88744,89321,115870,122401,122030,113802,105791,103076,98658,96945,92497,90687,88796,90015,113228,118711,117460,106556,97347,92657,93118,89037,83570,81693,75956,73993,97088,102394,96549,89727,82336,82653,82303,79596,74472,73562,66618,69029,89899,93774,90305,83799,80320,82497,84420,84646,84186,83269,77793,81145,101691,107357,104253,95963,91432,94324,93855,92183,87600,83641,78195,79604,100846,105293,102518,93132,87479,85476) > 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.9004081 beta : 0.0913872 gamma: 1 Coefficients: [,1] a 85093.8510931 b -452.0784540 s1 0.9896519 s2 0.9629889 s3 0.9222114 s4 0.8888677 s5 0.8365510 s6 0.8537640 s7 1.1038227 s8 1.1718693 s9 1.1560452 s10 1.0812993 s11 1.0162815 s12 1.0044909 > myresid <- x - fit$fitted[,'xhat'] > postscript(file="/var/www/wessaorg/rcomp/tmp/12tsm1305809432.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/wessaorg/rcomp/tmp/2gnat1305809432.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/wessaorg/rcomp/tmp/32zlg1305809432.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/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/wessaorg/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/wessaorg/rcomp/tmp/4kf8e1305809432.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/wessaorg/rcomp/tmp/5srip1305809432.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/wessaorg/rcomp/tmp/6zau21305809432.tab") > > try(system("convert tmp/12tsm1305809432.ps tmp/12tsm1305809432.png",intern=TRUE)) character(0) > try(system("convert tmp/2gnat1305809432.ps tmp/2gnat1305809432.png",intern=TRUE)) character(0) > try(system("convert tmp/32zlg1305809432.ps tmp/32zlg1305809432.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 2.430 0.250 2.766