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(1206.771667,1513.683333,1894.13,1514.458333,1788.191667,1951.223333,1727.026667,1552.248333,1569.528333,2203.916667,2402.373333,1611.8,1324.156667,1852.461667,2412.903333,1664,958.0283333) > 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 1627.453 1083.011 2171.896 794.8002 2460.106 19 1627.453 1080.295 2174.611 790.6473 2464.259 20 1627.453 1077.593 2177.313 786.5149 2468.392 21 1627.453 1074.904 2180.002 782.4027 2472.504 22 1627.453 1072.229 2182.678 778.3104 2476.596 23 1627.453 1069.566 2185.341 774.2378 2480.669 24 1627.453 1066.915 2187.991 770.1845 2484.722 25 1627.453 1064.277 2190.629 766.1502 2488.756 26 1627.453 1061.652 2193.254 762.1348 2492.772 27 1627.453 1059.039 2195.868 758.1380 2496.768 $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 262.925966 488.191738 402.819108 10.237443 23.364007 1.011039 Forecasts: Point Forecast 18 1627.453 19 1627.453 20 1627.453 21 1627.453 22 1627.453 23 1627.453 24 1627.453 25 1627.453 26 1627.453 27 1627.453 > postscript(file="/var/www/html/rcomp/tmp/1g4sj1273760417.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/2g4sj1273760417.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/3uw8r1273760417.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/4gwof1273760417.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/58ooi1273760417.tab") > try(system("convert tmp/1g4sj1273760417.ps tmp/1g4sj1273760417.png",intern=TRUE)) character(0) > try(system("convert tmp/2g4sj1273760417.ps tmp/2g4sj1273760417.png",intern=TRUE)) character(0) > > #-SERVER-wessa.org > > > > proc.time() user system elapsed 1.386 0.317 1.582