R version 2.9.0 (2009-04-17) Copyright (C) 2009 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. 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(275.0916667,401.3423333,260.6333333,286.1,368.0833333,385.1916667,181.5333333,145.1,203.1833333,227.5833333,239.0833333,109.175,231.0833333,216.9166667,229.8583333,272.7916667,155.0833333) > par10 = '0.1' > par9 = '3' > par8 = 'dumresult' > par7 = 'dum' > par6 = '12' > par5 = 'ZZZ' > par4 = 'NA' > par3 = 'NA' > par2 = 'Croston' > par1 = 'Input box' > if(par3!='NA') par3 <- as.numeric(par3) else par3 <- NA > if(par4!='NA') par4 <- as.numeric(par4) else par4 <- NA > par6 <- as.numeric(par6) #Seasonal Period > par9 <- as.numeric(par9) #Forecast Horizon > par10 <- as.numeric(par10) #Alpha > library(forecast) Loading required package: tseries Loading required package: quadprog Loading required package: zoo Attaching package: 'zoo' The following object(s) are masked from package:base : as.Date.numeric This is forecast 2.03 > if (par1 == 'CSV') { + xarr <- read.csv(file=paste('tmp/',par7,'.csv',sep=''),header=T) + numseries <- length(xarr[1,])-1 + n <- length(xarr[,1]) + nmh <- n - par9 + nmhp1 <- nmh + 1 + rarr <- array(NA,dim=c(n,numseries)) + farr <- array(NA,dim=c(n,numseries)) + parr <- array(NA,dim=c(numseries,8)) + colnames(parr) = list('ME','RMSE','MAE','MPE','MAPE','MASE','ACF1','TheilU') + for(i in 1:numseries) { + sindex <- i+1 + x <- xarr[,sindex] + if(par2=='Croston') { + if (i==1) m <- croston(x,alpha=par10) + if (i==1) mydemand <- m$model$demand[] + fit <- croston(x[1:nmh],h=par9,alpha=par10) + } + if(par2=='ARIMA') { + m <- auto.arima(ts(x,freq=par6),d=par3,D=par4) + mydemand <- forecast(m) + fit <- auto.arima(ts(x[1:nmh],freq=par6),d=par3,D=par4) + } + if(par2=='ETS') { + m <- ets(ts(x,freq=par6),model=par5) + mydemand <- forecast(m) + fit <- ets(ts(x[1:nmh],freq=par6),model=par5) + } + try(rarr[,i] <- mydemand$resid,silent=T) + try(farr[,i] <- mydemand$mean,silent=T) + if (par2!='Croston') parr[i,] <- accuracy(forecast(fit,par9),x[nmhp1:n]) + if (par2=='Croston') parr[i,] <- accuracy(fit,x[nmhp1:n]) + } + write.csv(farr,file=paste('tmp/',par8,'_f.csv',sep='')) + write.csv(rarr,file=paste('tmp/',par8,'_r.csv',sep='')) + write.csv(parr,file=paste('tmp/',par8,'_p.csv',sep='')) + } > if (par1 == 'Input box') { + numseries <- 1 + n <- length(x) + if(par2=='Croston') { + m <- croston(x) + mydemand <- m$model$demand[] + } + if(par2=='ARIMA') { + m <- auto.arima(ts(x,freq=par6),d=par3,D=par4) + mydemand <- forecast(m) + } + if(par2=='ETS') { + m <- ets(ts(x,freq=par6),model=par5) + mydemand <- forecast(m) + } + summary(m) + } Forecast method: Croston's method Model Information: $demand Point Forecast Lo 80 Hi 80 Lo 95 Hi 95 18 235.3446 131.2061 339.4830 76.07858 394.6105 19 235.3446 130.6867 340.0024 75.28423 395.4049 20 235.3446 130.1699 340.5192 74.49381 396.1953 21 235.3446 129.6556 341.0335 73.70725 396.9819 22 235.3446 129.1438 341.5453 72.92449 397.7646 23 235.3446 128.6344 342.0547 72.14550 398.5436 24 235.3446 128.1275 342.5616 71.37020 399.3189 25 235.3446 127.6229 343.0662 70.59855 400.0906 26 235.3446 127.1207 343.5684 69.83050 400.8586 27 235.3446 126.6209 344.0683 69.06600 401.6231 $period Point Forecast Lo 80 Hi 80 Lo 95 Hi 95 18 1 1 1 1 1 19 1 1 1 1 1 20 1 1 1 1 1 21 1 1 1 1 1 22 1 1 1 1 1 23 1 1 1 1 1 24 1 1 1 1 1 25 1 1 1 1 1 26 1 1 1 1 1 27 1 1 1 1 1 In-sample error measures: ME RMSE MAE MPE MAPE MASE -24.8419431 82.5079163 66.5650298 -23.4185170 34.6466813 0.9140726 Forecasts: Point Forecast 18 235.3446 19 235.3446 20 235.3446 21 235.3446 22 235.3446 23 235.3446 24 235.3446 25 235.3446 26 235.3446 27 235.3446 > postscript(file="/var/www/html/rcomp/tmp/1zbnn1273759001.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow=c(2,1)) > if (par2=='Croston') plot(m) > if ((par2=='ARIMA') | par2=='ETS') plot(forecast(m)) > plot(mydemand$resid,type='l',main='Residuals', ylab='residual value', xlab='time') > par(op) > dev.off() null device 1 > postscript(file="/var/www/html/rcomp/tmp/2x5te1273759001.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow=c(2,2)) > acf(mydemand$resid, lag.max=n/3, main='Residual ACF', ylab='autocorrelation', xlab='time lag') > pacf(mydemand$resid,lag.max=n/3, main='Residual PACF', ylab='partial autocorrelation', xlab='time lag') > cpgram(mydemand$resid, main='Cumulative Periodogram of Residuals') > qqnorm(mydemand$resid); qqline(mydemand$resid, col=2) > par(op) > dev.off() null device 1 > > #Note: the /var/www/html/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/html/rcomp/createtable") > > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Demand Forecast',6,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Point',header=TRUE) > a<-table.element(a,'Forecast',header=TRUE) > a<-table.element(a,'95% LB',header=TRUE) > a<-table.element(a,'80% LB',header=TRUE) > a<-table.element(a,'80% UB',header=TRUE) > a<-table.element(a,'95% UB',header=TRUE) > a<-table.row.end(a) > for (i in 1:length(mydemand$mean)) { + a<-table.row.start(a) + a<-table.element(a,i+n,header=TRUE) + a<-table.element(a,as.numeric(mydemand$mean[i])) + a<-table.element(a,as.numeric(mydemand$lower[i,2])) + a<-table.element(a,as.numeric(mydemand$lower[i,1])) + a<-table.element(a,as.numeric(mydemand$upper[i,1])) + a<-table.element(a,as.numeric(mydemand$upper[i,2])) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/rcomp/tmp/3el3e1273759001.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Actuals and Interpolation',3,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Time',header=TRUE) > a<-table.element(a,'Actual',header=TRUE) > a<-table.element(a,'Forecast',header=TRUE) > a<-table.row.end(a) > for (i in 1:n) { + a<-table.row.start(a) + a<-table.element(a,i,header=TRUE) + a<-table.element(a,x[i]) + a<-table.element(a,x[i] - as.numeric(m$resid[i])) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/rcomp/tmp/4z3j21273759001.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'What is next?',1,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,hyperlink(paste('http://www.wessa.net/Patrick.Wessa/rwasp_demand_forecasting_simulate.wasp',sep=''),'Simulate Time Series','',target='')) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,hyperlink(paste('http://www.wessa.net/Patrick.Wessa/rwasp_demand_forecasting_croston.wasp',sep=''),'Generate Forecasts','',target='')) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,hyperlink(paste('http://www.wessa.net/Patrick.Wessa/rwasp_demand_forecasting_analysis.wasp',sep=''),'Forecast Analysis','',target='')) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/www/html/rcomp/tmp/5sui41273759001.tab") > try(system("convert tmp/1zbnn1273759001.ps tmp/1zbnn1273759001.png",intern=TRUE)) character(0) > try(system("convert tmp/2x5te1273759001.ps tmp/2x5te1273759001.png",intern=TRUE)) character(0) > > #-SERVER-wessa.org > > > > proc.time() user system elapsed 1.384 0.334 1.555