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 <- array(list(62 + ,NA + ,NA + ,NA + ,32 + ,30 + ,62 + ,62 + ,62 + ,30 + ,31 + ,58.8 + ,58.8 + ,58.8 + ,70 + ,50 + ,56.02 + ,56.02 + ,56.02 + ,30 + ,33 + ,55.418 + ,55.418 + ,55.418 + ,30 + ,12 + ,53.1762 + ,53.1762 + ,53.1762 + ,10 + ,20 + ,49.05858 + ,49.05858 + ,49.05858 + ,30 + ,30 + ,46.152722 + ,46.152722 + ,46.152722 + ,30 + ,21.5 + ,44.5374498 + ,44.5374498 + ,44.5374498 + ,38 + ,23 + ,42.23370482 + ,42.23370482 + ,42.23370482 + ,20 + ,13.5 + ,40.310334338 + ,40.310334338 + ,40.310334338 + ,10 + ,0.5 + ,37.6293009042 + ,37.6293009042 + ,37.6293009042 + ,11 + ,12 + ,33.91637081378 + ,33.91637081378 + ,33.91637081378 + ,12 + ,10 + ,31.724733732402 + ,31.724733732402 + ,31.724733732402 + ,10 + ,70.5 + ,29.5522603591618 + ,29.5522603591618 + ,29.5522603591618 + ,30 + ,33.6470343232456 + ,33.6470343232456 + ,33.6470343232456 + ,31 + ,20.5 + ,33.2823308909211 + ,33.2823308909211 + ,33.2823308909211 + ,30 + ,12 + ,32.0040978018290 + ,32.0040978018290 + ,32.0040978018290 + ,12 + ,20 + ,30.0036880216461 + ,30.0036880216461 + ,30.0036880216461 + ,20 + ,45 + ,29.0033192194815 + ,29.0033192194815 + ,29.0033192194815 + ,10 + ,11.505 + ,30.6029872975333 + ,30.6029872975333 + ,30.6029872975333 + ,50 + ,0 + ,28.69318856778 + ,28.69318856778 + ,28.69318856778 + ,10 + ,10 + ,28.69318856778 + ,28.69318856778 + ,28.69318856778 + ,11 + ,5.5 + ,24.3853361009109 + ,24.3853361009109 + ,24.3853361009109 + ,1 + ,27.5 + ,22.6527364586255 + ,22.6527364586255 + ,22.6527364586255 + ,31 + ,0.5 + ,23.1011419666157 + ,23.1011419666157 + ,23.1011419666157 + ,1 + ,7 + ,20.994595040843 + ,20.994595040843 + ,20.994595040843 + ,20 + ,0 + ,19.6813007736305 + ,19.6813007736305 + ,19.6813007736305 + ,0 + ,2.5 + ,19.6813007736305 + ,19.6813007736305 + ,19.6813007736305 + ,0 + ,0 + ,16.5008829011108 + ,16.5008829011108 + ,16.5008829011108 + ,0 + ,0 + ,16.5008829011108 + ,16.5008829011108 + ,16.5008829011108 + ,0 + ,6.025 + ,16.5008829011108 + ,16.5008829011108 + ,16.5008829011108 + ,1 + ,1 + ,13.2638773777770 + ,13.2638773777770 + ,13.2638773777770 + ,0 + ,0 + ,12.3268664380437 + ,12.3268664380437 + ,12.3268664380437 + ,0 + ,0 + ,12.3268664380437 + ,12.3268664380437 + ,12.3268664380437 + ,0 + ,0 + ,12.3268664380437 + ,12.3268664380437 + ,12.3268664380437 + ,0 + ,0 + ,12.3268664380437 + ,12.3268664380437 + ,12.3268664380437 + ,0 + ,2 + ,12.3268664380437 + ,12.3268664380437 + ,12.3268664380437 + ,2 + ,0 + ,8.77286013713574 + ,8.77286013713574 + ,8.77286013713574 + ,2 + ,6 + ,8.77286013713574 + ,8.77286013713574 + ,8.77286013713574 + ,0 + ,20 + ,8.09773184771271 + ,8.09773184771271 + ,8.09773184771271 + ,0 + ,0 + ,8.82386207874147 + ,8.82386207874147 + ,8.82386207874147 + ,0 + ,0 + ,8.82386207874147 + ,8.82386207874147 + ,8.82386207874147 + ,0 + ,0 + ,8.82386207874147 + ,8.82386207874147 + ,8.82386207874147 + ,0 + ,7 + ,8.82386207874147 + ,8.82386207874147 + ,8.82386207874147 + ,5 + ,35 + ,7.3149510902026 + ,7.3149510902026 + ,7.3149510902026 + ,0 + ,0 + ,8.86359100991655 + ,8.86359100991655 + ,8.86359100991655 + ,0 + ,0 + ,8.86359100991655 + ,8.86359100991655 + ,8.86359100991655 + ,0 + ,0 + ,8.86359100991655 + ,8.86359100991655 + ,8.86359100991655 + ,0 + ,1 + ,8.86359100991655 + ,8.86359100991655 + ,8.86359100991655 + ,0) + ,dim=c(5 + ,50) + ,dimnames=list(c('Actuals' + ,'Croston' + ,'ETS' + ,'Arima' + ,'Kaneka') + ,1:50)) > y <- array(NA,dim=c(5,50),dimnames=list(c('Actuals','Croston','ETS','Arima','Kaneka'),1:50)) > for (i in 1:dim(x)[1]) + { + for (j in 1:dim(x)[2]) + { + y[i,j] <- as.numeric(x[i,j]) + } + } > par2 = 'No' > par1 = '1' > #'GNU S' R Code compiled by R2WASP v. 1.0.44 () > #Author: Dr. Ian E. Holliday > #To cite this work: Ian E. Holliday, 2009, YOUR SOFTWARE TITLE (vNUMBER) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_YOURPAGE.wasp/ > #Source of accompanying publication: > #Technical description: > library(rpart) > library(partykit) Loading required package: grid Loading required package: mvtnorm > par1 <- as.numeric(par1) > autoprune <- function ( tree, method='Minimum CV'){ + xerr <- tree$cptable[,'xerror'] + cpmin.id <- which.min(xerr) + if (method == 'Minimum CV Error plus 1 SD'){ + xstd <- tree$cptable[,'xstd'] + errt <- xerr[cpmin.id] + xstd[cpmin.id] + cpSE1.min <- which.min( errt < xerr ) + mycp <- (tree$cptable[,'CP'])[cpSE1.min] + } + if (method == 'Minimum CV') { + mycp <- (tree$cptable[,'CP'])[cpmin.id] + } + return (mycp) + } > conf.multi.mat <- function(true, new) + { + if ( all( is.na(match( levels(true),levels(new) ) )) ) + stop ( 'conflict of vector levels') + multi.t <- list() + for (mylev in levels(true) ) { + true.tmp <- true + new.tmp <- new + left.lev <- levels (true.tmp)[- match(mylev,levels(true) ) ] + levels(true.tmp) <- list ( mylev = mylev, all = left.lev ) + levels(new.tmp) <- list ( mylev = mylev, all = left.lev ) + curr.t <- conf.mat ( true.tmp , new.tmp ) + multi.t[[mylev]] <- curr.t + multi.t[[mylev]]$precision <- + round( curr.t$conf[1,1] / sum( curr.t$conf[1,] ), 2 ) + } + return (multi.t) + } > x <- t(y) > k <- length(x[1,]) > n <- length(x[,1]) > x1 <- cbind(x[,par1], x[,1:k!=par1]) > mycolnames <- c(colnames(x)[par1], colnames(x)[1:k!=par1]) > colnames(x1) <- mycolnames #colnames(x)[par1] > m <- rpart(as.data.frame(x1)) > par2 [1] "No" > if (par2 != 'No') { + mincp <- autoprune(m,method=par2) + print(mincp) + m <- prune(m,cp=mincp) + } > m$cptable CP nsplit rel error xerror xstd 1 0.36057508 0 1.0000000 1.0776705 0.3429171 2 0.17711860 1 0.6394249 0.9656173 0.3252052 3 0.01488550 2 0.4623063 0.8074581 0.2233549 4 0.01000000 3 0.4474208 0.7791990 0.2232878 > postscript(file="/var/www/html/rcomp/tmp/1obdj1275298823.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(as.party(m),tp_args=list(id=FALSE)) > dev.off() null device 1 > postscript(file="/var/www/html/rcomp/tmp/2obdj1275298823.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plotcp(m) > dev.off() null device 1 > cbind(y=m$y,pred=predict(m),res=residuals(m)) y pred res 1 62.000000 37.800369 24.1996312 2 30.000000 37.800369 -7.8003688 3 31.000000 37.800369 -6.8003688 4 50.000000 37.800369 12.1996312 5 33.000000 37.800369 -4.8003688 6 12.000000 21.969608 -9.9696083 7 20.000000 37.800369 -17.8003688 8 30.000000 37.800369 -7.8003688 9 21.500000 37.800369 -16.3003688 10 23.000000 21.969608 1.0303917 11 13.500000 21.969608 -8.4696083 12 0.500000 21.969608 -21.4696083 13 12.000000 21.969608 -9.9696083 14 10.000000 21.969608 -11.9696083 15 70.500000 37.800369 32.6996312 16 33.647034 21.969608 11.6774260 17 33.282331 21.969608 11.3127226 18 32.004098 21.969608 10.0344895 19 30.003688 37.800369 -7.7966808 20 29.003319 21.969608 7.0337109 21 30.602987 21.969608 8.6333790 22 28.693189 21.969608 6.7235803 23 28.693189 21.969608 6.7235803 24 24.385336 21.969608 2.4157278 25 22.652736 21.969608 0.6831282 26 23.101142 21.969608 1.1315337 27 20.994595 21.969608 -0.9750133 28 19.681301 21.969608 -2.2883075 29 19.681301 21.969608 -2.2883075 30 16.500883 13.822318 2.6785653 31 16.500883 13.822318 2.6785653 32 16.500883 13.822318 2.6785653 33 13.263877 13.822318 -0.5584402 34 12.326866 13.822318 -1.4954511 35 12.326866 13.822318 -1.4954511 36 12.326866 13.822318 -1.4954511 37 12.326866 13.822318 -1.4954511 38 12.326866 13.822318 -1.4954511 39 8.772860 8.642351 0.1305088 40 8.772860 8.642351 0.1305088 41 8.097732 8.642351 -0.5446194 42 8.823862 8.642351 0.1815108 43 8.823862 8.642351 0.1815108 44 8.823862 8.642351 0.1815108 45 8.823862 8.642351 0.1815108 46 7.314951 8.642351 -1.3274002 47 8.863591 8.642351 0.2212397 48 8.863591 8.642351 0.2212397 49 8.863591 8.642351 0.2212397 50 8.863591 8.642351 0.2212397 > myr <- residuals(m) > myp <- predict(m) > postscript(file="/var/www/html/rcomp/tmp/3y2c41275298823.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow=c(2,2)) > plot(myr,ylab='residuals') > plot(density(myr),main='Residual Kernel Density') > plot(myp,myr,xlab='predicted',ylab='residuals',main='Predicted vs Residuals') > plot(density(myp),main='Prediction Kernel Density') > 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,'Model Performance',6,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'#',header=TRUE) > a<-table.element(a,'Complexity',header=TRUE) > a<-table.element(a,'split',header=TRUE) > a<-table.element(a,'relative error',header=TRUE) > a<-table.element(a,'CV error',header=TRUE) > a<-table.element(a,'CV S.D.',header=TRUE) > a<-table.row.end(a) > for (i in 1:length(m$cptable[,1])) { + a<-table.row.start(a) + a<-table.element(a,i,header=TRUE) + a<-table.element(a,round(m$cptable[i,'CP'],3)) + a<-table.element(a,m$cptable[i,'nsplit']) + a<-table.element(a,round(m$cptable[i,'rel error'],3)) + a<-table.element(a,round(m$cptable[i,'xerror'],3)) + a<-table.element(a,round(m$cptable[i,'xstd'],3)) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/rcomp/tmp/4yv8i1275298823.tab") > > try(system("convert tmp/1obdj1275298823.ps tmp/1obdj1275298823.png",intern=TRUE)) character(0) > try(system("convert tmp/2obdj1275298823.ps tmp/2obdj1275298823.png",intern=TRUE)) character(0) > try(system("convert tmp/3y2c41275298823.ps tmp/3y2c41275298823.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 1.098 0.495 4.565