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(1201.42 + ,NA + ,1201.52699174648 + ,1638.54139885537 + ,1156.65 + ,1157.2125 + ,1201.42 + ,1157.41646441400 + ,1497.05173429845 + ,1504.775 + ,1722.05 + ,1196.99925 + ,1721.96970488708 + ,1595.11391392574 + ,1509.55 + ,1918.2475 + ,1249.504325 + ,1918.17344273383 + ,1810.81733236927 + ,1635.05 + ,1930.4575 + ,1316.3786425 + ,1930.48017222676 + ,1804.12932358113 + ,1856.9 + ,1264.1775 + ,1377.78652825 + ,1263.99969486656 + ,1735.32747882431 + ,1329.875 + ,1456.3725 + ,1366.425625425 + ,1456.29994042014 + ,1507.49011430111 + ,1361.5 + ,2168.985 + ,1375.4203128825 + ,2168.79723376114 + ,1765.26576138598 + ,1038 + ,1983.765 + ,1454.77678159425 + ,1983.89302511237 + ,1934.23830420204 + ,1537 + ,1672.695 + ,1507.67560343483 + ,1672.92990971789 + ,1750.23627132221 + ,1834.525 + ,1938.575 + ,1524.17754309134 + ,1938.52541901855 + ,1727.68342690282 + ,2001.9 + ,1307.6425 + ,1565.61728878221 + ,1307.50607718124 + ,1736.51705207453 + ,1638 + ,1523.3425 + ,1539.81980990399 + ,1523.21360851071 + ,1352.55910680891 + ,1466.4 + ,1928.39 + ,1538.17207891359 + ,1928.25609237535 + ,1581.65625637482 + ,1775.4 + ,2208.435 + ,1577.19387102223 + ,2208.49199362838 + ,1889.40575429802 + ,2066.75 + ,2290.175 + ,1640.31798392001 + ,2290.220544191 + ,1958.87684384425 + ,2250.05 + ,2578.245 + ,1705.30368552801 + ,2578.19488211031 + ,1977.21467180681 + ,2564.65 + ,1152.84 + ,1792.59781697521 + ,1153.01584440369 + ,1782.30313239910 + ,1980 + ,1398.7575 + ,1728.62203527769 + ,1398.7598374143 + ,1351.50862585308 + ,1476.375 + ,1393.9175 + ,1695.63558174992 + ,1394.14187074320 + ,1919.78381416109 + ,1667.75 + ,1972.2525 + ,1665.46377357493 + ,1972.09241634781 + ,1622.49344537309 + ,914.15 + ,2410.4775 + ,1696.14264621743 + ,2410.25436087483 + ,1876.91928728301 + ,2102.025 + ,2363.27 + ,1767.57613159569 + ,2363.22252889700 + ,2058.90251344068 + ,2366.75 + ,1341.6075 + ,1827.14551843612 + ,1341.65173666235 + ,1699.68606008630 + ,1588 + ,1437.0425 + ,1778.59171659251 + ,1437.06301486057 + ,1590.45807134907 + ,1617) + ,dim=c(5 + ,25) + ,dimnames=list(c('Actuals' + ,'Croston' + ,'ETS' + ,'Arima' + ,'Kaneka') + ,1:25)) > y <- array(NA,dim=c(5,25),dimnames=list(c('Actuals','Croston','ETS','Arima','Kaneka'),1:25)) > 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.7871319 0 1.0000000 1.0462518 0.19446051 2 0.0100000 1 0.2128681 0.2440819 0.05280422 > postscript(file="/var/www/html/rcomp/tmp/1gu051274873909.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/29li81274873909.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 1201.420 1386.852 -185.432115 2 1157.213 1386.852 -229.639615 3 1722.050 1386.852 335.197885 4 1918.247 2140.940 -222.692083 5 1930.457 2140.940 -210.482083 6 1264.178 1386.852 -122.674615 7 1456.372 1386.852 69.520385 8 2168.985 2140.940 28.045417 9 1983.765 2140.940 -157.174583 10 1672.695 1386.852 285.842885 11 1938.575 2140.940 -202.364583 12 1307.642 1386.852 -79.209615 13 1523.342 1386.852 136.490385 14 1928.390 2140.940 -212.549583 15 2208.435 2140.940 67.495417 16 2290.175 2140.940 149.235417 17 2578.245 2140.940 437.305417 18 1152.840 1386.852 -234.012115 19 1398.757 1386.852 11.905385 20 1393.918 1386.852 7.065385 21 1972.253 2140.940 -168.687083 22 2410.477 2140.940 269.537917 23 2363.270 2140.940 222.330417 24 1341.608 1386.852 -45.244615 25 1437.043 1386.852 50.190385 > myr <- residuals(m) > myp <- predict(m) > postscript(file="/var/www/html/rcomp/tmp/39li81274873909.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/45dxh1274873909.tab") > > try(system("convert tmp/1gu051274873909.ps tmp/1gu051274873909.png",intern=TRUE)) character(0) > try(system("convert tmp/29li81274873909.ps tmp/29li81274873909.png",intern=TRUE)) character(0) > try(system("convert tmp/39li81274873909.ps tmp/39li81274873909.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 0.922 0.450 1.210