x <- array(list(272545 ,1747 ,69 ,483 ,32033 ,3 ,116 ,144 ,179444 ,1209 ,64 ,429 ,20654 ,4 ,127 ,133 ,222373 ,1844 ,69 ,673 ,16346 ,16 ,106 ,162 ,218443 ,2683 ,104 ,1137 ,35926 ,2 ,133 ,148 ,167843 ,1228 ,51 ,374 ,10621 ,1 ,64 ,88 ,70849 ,631 ,28 ,179 ,10024 ,3 ,89 ,129 ,506574 ,4627 ,123 ,2251 ,43068 ,0 ,122 ,128 ,33186 ,381 ,19 ,111 ,1271 ,0 ,22 ,67 ,216660 ,2063 ,59 ,740 ,34416 ,7 ,117 ,132 ,213274 ,1758 ,44 ,595 ,20318 ,0 ,82 ,120 ,307153 ,2132 ,109 ,800 ,24409 ,0 ,147 ,169 ,237633 ,2128 ,114 ,660 ,20648 ,7 ,184 ,210 ,164292 ,1667 ,68 ,635 ,12347 ,8 ,113 ,122 ,364402 ,2965 ,79 ,1172 ,21857 ,4 ,171 ,191 ,244103 ,2098 ,84 ,674 ,11034 ,10 ,87 ,162 ,384448 ,4904 ,178 ,1692 ,33433 ,0 ,199 ,223 ,325587 ,2242 ,68 ,811 ,35902 ,6 ,139 ,140 ,323652 ,2977 ,157 ,1168 ,22355 ,4 ,92 ,144 ,176082 ,1438 ,55 ,507 ,31219 ,3 ,85 ,111 ,266736 ,2347 ,87 ,689 ,21983 ,8 ,193 ,199 ,278265 ,2522 ,70 ,837 ,40085 ,0 ,160 ,187 ,442703 ,2889 ,103 ,1270 ,18507 ,1 ,144 ,144 ,180393 ,1447 ,41 ,462 ,16278 ,5 ,84 ,89 ,189897 ,1717 ,54 ,601 ,24662 ,9 ,208 ,208 ,234247 ,3362 ,121 ,1242 ,31452 ,1 ,154 ,165 ,237452 ,2898 ,125 ,1025 ,32580 ,0 ,139 ,146 ,267268 ,2828 ,127 ,1062 ,22883 ,5 ,127 ,158 ,270787 ,1972 ,86 ,618 ,27652 ,0 ,148 ,154 ,155915 ,1495 ,51 ,559 ,9845 ,0 ,99 ,117 ,342564 ,2840 ,69 ,1062 ,20190 ,0 ,135 ,158 ,282172 ,2299 ,76 ,913 ,46201 ,3 ,171 ,183 ,216584 ,1909 ,76 ,643 ,10971 ,6 ,149 ,186 ,318563 ,2091 ,84 ,779 ,34811 ,1 ,178 ,185 ,98672 ,971 ,37 ,322 ,3029 ,4 ,137 ,141 ,386258 ,3293 ,95 ,1243 ,38941 ,4 ,151 ,156 ,273950 ,2764 ,56 ,1186 ,4958 ,0 ,127 ,159 ,425120 ,3682 ,120 ,1324 ,32344 ,0 ,151 ,161 ,227636 ,1918 ,83 ,640 ,19433 ,2 ,89 ,139 ,115658 ,947 ,33 ,284 ,12558 ,1 ,46 ,55 ,349863 ,3433 ,194 ,1210 ,36524 ,2 ,153 ,163 ,324178 ,3246 ,79 ,1490 ,26041 ,10 ,122 ,145 ,178083 ,1692 ,67 ,667 ,16637 ,9 ,111 ,148 ,195153 ,1735 ,73 ,635 ,28395 ,5 ,108 ,115 ,177694 ,1771 ,61 ,479 ,16747 ,6 ,142 ,174 ,153778 ,2496 ,82 ,1022 ,9105 ,1 ,45 ,73 ,455168 ,5501 ,151 ,2068 ,11941 ,2 ,131 ,147 ,78800 ,918 ,42 ,330 ,7935 ,2 ,66 ,82 ,208051 ,2228 ,76 ,648 ,19499 ,0 ,180 ,201 ,348077 ,4051 ,118 ,1367 ,22938 ,10 ,165 ,181 ,175523 ,2081 ,54 ,868 ,25314 ,3 ,146 ,164 ,224591 ,1875 ,74 ,588 ,28524 ,0 ,137 ,158 ,24188 ,496 ,24 ,218 ,2694 ,0 ,7 ,12 ,372238 ,2537 ,314 ,833 ,20867 ,8 ,157 ,163 ,65029 ,744 ,17 ,255 ,3597 ,5 ,61 ,67 ,101097 ,1161 ,64 ,454 ,5296 ,3 ,41 ,52 ,279012 ,3027 ,58 ,1108 ,32982 ,1 ,120 ,134 ,317644 ,2526 ,84 ,662 ,38975 ,5 ,228 ,230 ,340471 ,3705 ,185 ,1119 ,42721 ,5 ,137 ,145 ,358958 ,2667 ,141 ,1058 ,41455 ,0 ,150 ,153 ,252529 ,2175 ,83 ,822 ,23923 ,12 ,127 ,139 ,370628 ,3949 ,140 ,1302 ,26719 ,10 ,161 ,178 ,304468 ,3165 ,117 ,1145 ,53405 ,12 ,73 ,101 ,265870 ,2939 ,113 ,1185 ,12526 ,11 ,97 ,169 ,264889 ,2610 ,88 ,931 ,26584 ,8 ,142 ,163 ,228595 ,1426 ,66 ,557 ,37062 ,2 ,125 ,139 ,216027 ,1646 ,65 ,436 ,25696 ,0 ,87 ,116 ,198798 ,1971 ,132 ,596 ,24634 ,6 ,128 ,137 ,238146 ,2746 ,145 ,837 ,27269 ,9 ,148 ,167 ,234891 ,2308 ,81 ,848 ,25270 ,2 ,116 ,135 ,175816 ,1684 ,69 ,625 ,24634 ,5 ,89 ,102 ,239314 ,2537 ,68 ,865 ,17828 ,13 ,154 ,173 ,73566 ,893 ,32 ,385 ,3007 ,6 ,67 ,88 ,242622 ,2195 ,84 ,718 ,20065 ,7 ,171 ,175 ,187167 ,1695 ,53 ,705 ,24648 ,2 ,90 ,133 ,209049 ,2061 ,63 ,732 ,21588 ,2 ,133 ,148 ,360592 ,2329 ,86 ,988 ,25217 ,4 ,144 ,169 ,342846 ,2695 ,92 ,1077 ,30927 ,3 ,133 ,140 ,207650 ,1809 ,107 ,524 ,18487 ,6 ,125 ,154 ,206500 ,2290 ,62 ,697 ,18050 ,2 ,134 ,148 ,182357 ,1791 ,64 ,644 ,17696 ,0 ,110 ,134 ,153613 ,1678 ,46 ,622 ,17326 ,1 ,89 ,109 ,456979 ,4023 ,124 ,1227 ,39361 ,0 ,138 ,175 ,145943 ,1369 ,69 ,653 ,9648 ,5 ,99 ,99 ,280366 ,2308 ,104 ,656 ,26759 ,2 ,92 ,122 ,80953 ,870 ,25 ,437 ,7905 ,0 ,27 ,28 ,150216 ,1966 ,54 ,822 ,4527 ,0 ,77 ,101 ,167878 ,1459 ,59 ,423 ,41517 ,6 ,137 ,139 ,369718 ,3795 ,205 ,1489 ,21261 ,1 ,137 ,143 ,322454 ,2673 ,116 ,929 ,36099 ,0 ,122 ,206 ,179797 ,3085 ,104 ,1044 ,39039 ,1 ,159 ,171 ,262883 ,2367 ,91 ,792 ,13841 ,1 ,85 ,138 ,262793 ,2209 ,77 ,678 ,23841 ,3 ,138 ,148 ,189142 ,1829 ,63 ,597 ,8589 ,9 ,90 ,114 ,275997 ,3087 ,74 ,1099 ,15049 ,1 ,135 ,140 ,328875 ,2559 ,82 ,966 ,39038 ,4 ,147 ,156 ,189252 ,1624 ,36 ,555 ,30391 ,3 ,139 ,140 ,222504 ,1607 ,51 ,552 ,39932 ,5 ,127 ,127 ,287386 ,2109 ,79 ,778 ,43840 ,0 ,104 ,141 ,389104 ,4015 ,151 ,1322 ,43146 ,12 ,248 ,251 ,397681 ,3705 ,108 ,1415 ,50099 ,13 ,106 ,114 ,287748 ,2714 ,136 ,853 ,40312 ,8 ,176 ,198 ,294320 ,2325 ,65 ,848 ,32616 ,0 ,130 ,155 ,186856 ,1999 ,179 ,640 ,11338 ,0 ,59 ,138 ,43287 ,602 ,14 ,214 ,7409 ,4 ,64 ,71 ,185468 ,2146 ,80 ,716 ,18213 ,4 ,36 ,84 ,235352 ,2325 ,146 ,795 ,45873 ,0 ,98 ,167 ,268077 ,2617 ,48 ,1170 ,39844 ,0 ,125 ,155 ,305195 ,2688 ,90 ,1048 ,28317 ,0 ,124 ,150 ,143356 ,1207 ,72 ,399 ,24797 ,0 ,83 ,112 ,154165 ,3102 ,88 ,906 ,7471 ,0 ,127 ,161 ,307000 ,1869 ,68 ,609 ,27259 ,4 ,143 ,149 ,298039 ,2304 ,88 ,688 ,23201 ,0 ,115 ,164 ,23623 ,398 ,11 ,156 ,238 ,0 ,0 ,0 ,195817 ,2205 ,73 ,779 ,28830 ,0 ,103 ,155 ,61857 ,530 ,25 ,192 ,3913 ,4 ,30 ,32 ,163766 ,1596 ,48 ,457 ,9935 ,0 ,119 ,169 ,414506 ,3083 ,117 ,1195 ,27738 ,1 ,102 ,140 ,21054 ,387 ,16 ,146 ,338 ,0 ,0 ,0 ,252805 ,2137 ,52 ,866 ,13326 ,5 ,77 ,111 ,31961 ,492 ,22 ,200 ,3988 ,0 ,9 ,25 ,317367 ,3450 ,115 ,1230 ,24347 ,3 ,137 ,146 ,240153 ,2089 ,65 ,696 ,27111 ,7 ,163 ,183 ,175083 ,1658 ,88 ,491 ,3938 ,13 ,146 ,181 ,152043 ,1685 ,53 ,670 ,17416 ,3 ,84 ,107 ,38214 ,568 ,34 ,276 ,1888 ,0 ,21 ,27 ,216299 ,2059 ,42 ,716 ,18700 ,2 ,151 ,163 ,357602 ,2792 ,82 ,1021 ,36809 ,0 ,187 ,198 ,198104 ,1395 ,61 ,481 ,24959 ,0 ,171 ,205 ,410803 ,3590 ,80 ,1582 ,37343 ,4 ,167 ,187 ,316105 ,2387 ,97 ,820 ,21849 ,0 ,145 ,187 ,397297 ,3334 ,124 ,1153 ,49809 ,3 ,175 ,186 ,187992 ,1250 ,35 ,473 ,21654 ,0 ,137 ,151 ,102424 ,1121 ,42 ,401 ,8728 ,0 ,100 ,131 ,286327 ,2880 ,335 ,954 ,20920 ,4 ,150 ,155 ,407378 ,4104 ,170 ,1447 ,27195 ,4 ,163 ,172 ,143860 ,1759 ,54 ,546 ,1037 ,15 ,149 ,160 ,391854 ,4138 ,132 ,1728 ,42570 ,2 ,161 ,172 ,157429 ,1831 ,77 ,689 ,17672 ,4 ,112 ,143 ,258751 ,1787 ,48 ,590 ,34245 ,2 ,135 ,151 ,282399 ,2535 ,94 ,897 ,16786 ,1 ,124 ,158 ,217665 ,1816 ,113 ,613 ,20954 ,0 ,45 ,125 ,366774 ,3873 ,116 ,1548 ,16378 ,9 ,120 ,145 ,236660 ,2181 ,88 ,759 ,31852 ,1 ,126 ,145 ,173260 ,2035 ,63 ,716 ,2805 ,3 ,78 ,79 ,323545 ,2960 ,99 ,955 ,38086 ,11 ,136 ,174 ,168994 ,1915 ,57 ,720 ,21166 ,5 ,179 ,192 ,253330 ,2648 ,86 ,1023 ,34672 ,2 ,118 ,132 ,301703 ,2633 ,105 ,818 ,36171 ,1 ,147 ,159 ,1 ,2 ,0 ,0 ,0 ,9 ,0 ,0 ,14688 ,207 ,10 ,85 ,2065 ,0 ,0 ,0 ,98 ,5 ,1 ,0 ,0 ,0 ,0 ,0 ,455 ,8 ,2 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,246435 ,2116 ,84 ,737 ,19354 ,2 ,88 ,133 ,382374 ,3286 ,154 ,1080 ,22124 ,3 ,129 ,204 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,203 ,4 ,4 ,0 ,0 ,0 ,0 ,0 ,7199 ,151 ,5 ,74 ,556 ,0 ,0 ,0 ,46660 ,474 ,20 ,259 ,2089 ,0 ,13 ,15 ,17547 ,141 ,5 ,69 ,2658 ,0 ,4 ,4 ,116678 ,1047 ,42 ,285 ,1813 ,0 ,76 ,152 ,969 ,29 ,2 ,0 ,0 ,0 ,0 ,0 ,206501 ,1822 ,68 ,591 ,17372 ,2 ,71 ,125) ,dim=c(8 ,164) ,dimnames=list(c('Y' ,'X1' ,'X2' ,'X3' ,'X4' ,'X5' ,'X6' ,'X7') ,1:164)) y <- array(NA,dim=c(8,164),dimnames=list(c('Y','X1','X2','X3','X4','X5','X6','X7'),1:164)) for (i in 1:dim(x)[1]) { for (j in 1:dim(x)[2]) { y[i,j] <- as.numeric(x[i,j]) } } par20 = '' par19 = '' par18 = '' par17 = '' par16 = '' par15 = '' par14 = '' par13 = '' par12 = '' par11 = '' par10 = '' par9 = '' par8 = '' par7 = '' par6 = '' par5 = '' par4 = 'no' par3 = '2' par2 = 'quantiles' par1 = '1' library(party) library(Hmisc) par1 <- as.numeric(par1) par3 <- as.numeric(par3) x <- data.frame(t(y)) is.data.frame(x) x <- x[!is.na(x[,par1]),] k <- length(x[1,]) n <- length(x[,1]) colnames(x)[par1] x[,par1] if (par2 == 'kmeans') { cl <- kmeans(x[,par1], par3) print(cl) clm <- matrix(cbind(cl$centers,1:par3),ncol=2) clm <- clm[sort.list(clm[,1]),] for (i in 1:par3) { cl$cluster[cl$cluster==clm[i,2]] <- paste('C',i,sep='') } cl$cluster <- as.factor(cl$cluster) print(cl$cluster) x[,par1] <- cl$cluster } if (par2 == 'quantiles') { x[,par1] <- cut2(x[,par1],g=par3) } if (par2 == 'hclust') { hc <- hclust(dist(x[,par1])^2, 'cen') print(hc) memb <- cutree(hc, k = par3) dum <- c(mean(x[memb==1,par1])) for (i in 2:par3) { dum <- c(dum, mean(x[memb==i,par1])) } hcm <- matrix(cbind(dum,1:par3),ncol=2) hcm <- hcm[sort.list(hcm[,1]),] for (i in 1:par3) { memb[memb==hcm[i,2]] <- paste('C',i,sep='') } memb <- as.factor(memb) print(memb) x[,par1] <- memb } if (par2=='equal') { ed <- cut(as.numeric(x[,par1]),par3,labels=paste('C',1:par3,sep='')) x[,par1] <- as.factor(ed) } table(x[,par1]) colnames(x) colnames(x)[par1] x[,par1] if (par2 == 'none') { m <- ctree(as.formula(paste(colnames(x)[par1],' ~ .',sep='')),data = x) } #Note: the /var/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab load(file="/var/wessaorg/rcomp/createtable") if (par2 != 'none') { m <- ctree(as.formula(paste('as.factor(',colnames(x)[par1],') ~ .',sep='')),data = x) if (par4=='yes') { a<-table.start() a<-table.row.start(a) a<-table.element(a,'10-Fold Cross Validation',3+2*par3,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'',1,TRUE) a<-table.element(a,'Prediction (training)',par3+1,TRUE) a<-table.element(a,'Prediction (testing)',par3+1,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Actual',1,TRUE) for (jjj in 1:par3) a<-table.element(a,paste('C',jjj,sep=''),1,TRUE) a<-table.element(a,'CV',1,TRUE) for (jjj in 1:par3) a<-table.element(a,paste('C',jjj,sep=''),1,TRUE) a<-table.element(a,'CV',1,TRUE) a<-table.row.end(a) for (i in 1:10) { ind <- sample(2, nrow(x), replace=T, prob=c(0.9,0.1)) m.ct <- ctree(as.formula(paste('as.factor(',colnames(x)[par1],') ~ .',sep='')),data =x[ind==1,]) if (i==1) { m.ct.i.pred <- predict(m.ct, newdata=x[ind==1,]) m.ct.i.actu <- x[ind==1,par1] m.ct.x.pred <- predict(m.ct, newdata=x[ind==2,]) m.ct.x.actu <- x[ind==2,par1] } else { m.ct.i.pred <- c(m.ct.i.pred,predict(m.ct, newdata=x[ind==1,])) m.ct.i.actu <- c(m.ct.i.actu,x[ind==1,par1]) m.ct.x.pred <- c(m.ct.x.pred,predict(m.ct, newdata=x[ind==2,])) m.ct.x.actu <- c(m.ct.x.actu,x[ind==2,par1]) } } print(m.ct.i.tab <- table(m.ct.i.actu,m.ct.i.pred)) numer <- 0 for (i in 1:par3) { print(m.ct.i.tab[i,i] / sum(m.ct.i.tab[i,])) numer <- numer + m.ct.i.tab[i,i] } print(m.ct.i.cp <- numer / sum(m.ct.i.tab)) print(m.ct.x.tab <- table(m.ct.x.actu,m.ct.x.pred)) numer <- 0 for (i in 1:par3) { print(m.ct.x.tab[i,i] / sum(m.ct.x.tab[i,])) numer <- numer + m.ct.x.tab[i,i] } print(m.ct.x.cp <- numer / sum(m.ct.x.tab)) for (i in 1:par3) { a<-table.row.start(a) a<-table.element(a,paste('C',i,sep=''),1,TRUE) for (jjj in 1:par3) a<-table.element(a,m.ct.i.tab[i,jjj]) a<-table.element(a,round(m.ct.i.tab[i,i]/sum(m.ct.i.tab[i,]),4)) for (jjj in 1:par3) a<-table.element(a,m.ct.x.tab[i,jjj]) a<-table.element(a,round(m.ct.x.tab[i,i]/sum(m.ct.x.tab[i,]),4)) a<-table.row.end(a) } a<-table.row.start(a) a<-table.element(a,'Overall',1,TRUE) for (jjj in 1:par3) a<-table.element(a,'-') a<-table.element(a,round(m.ct.i.cp,4)) for (jjj in 1:par3) a<-table.element(a,'-') a<-table.element(a,round(m.ct.x.cp,4)) a<-table.row.end(a) a<-table.end(a) table.save(a,file="/var/wessaorg/rcomp/tmp/1ggr81324406163.tab") } } m postscript(file="/var/wessaorg/rcomp/tmp/2pafv1324406163.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) plot(m) dev.off() postscript(file="/var/wessaorg/rcomp/tmp/3kvx81324406163.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) plot(x[,par1] ~ as.factor(where(m)),main='Response by Terminal Node',xlab='Terminal Node',ylab='Response') dev.off() if (par2 == 'none') { forec <- predict(m) result <- as.data.frame(cbind(x[,par1],forec,x[,par1]-forec)) colnames(result) <- c('Actuals','Forecasts','Residuals') print(result) } if (par2 != 'none') { print(cbind(as.factor(x[,par1]),predict(m))) myt <- table(as.factor(x[,par1]),predict(m)) print(myt) } postscript(file="/var/wessaorg/rcomp/tmp/4ig241324406163.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) if(par2=='none') { op <- par(mfrow=c(2,2)) plot(density(result$Actuals),main='Kernel Density Plot of Actuals') plot(density(result$Residuals),main='Kernel Density Plot of Residuals') plot(result$Forecasts,result$Actuals,main='Actuals versus Predictions',xlab='Predictions',ylab='Actuals') plot(density(result$Forecasts),main='Kernel Density Plot of Predictions') par(op) } if(par2!='none') { plot(myt,main='Confusion Matrix',xlab='Actual',ylab='Predicted') } dev.off() if (par2 == 'none') { detcoef <- cor(result$Forecasts,result$Actuals) a<-table.start() a<-table.row.start(a) a<-table.element(a,'Goodness of Fit',2,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Correlation',1,TRUE) a<-table.element(a,round(detcoef,4)) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'R-squared',1,TRUE) a<-table.element(a,round(detcoef*detcoef,4)) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'RMSE',1,TRUE) a<-table.element(a,round(sqrt(mean((result$Residuals)^2)),4)) a<-table.row.end(a) a<-table.end(a) table.save(a,file="/var/wessaorg/rcomp/tmp/5hj4x1324406163.tab") a<-table.start() a<-table.row.start(a) a<-table.element(a,'Actuals, Predictions, and Residuals',4,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'#',header=TRUE) a<-table.element(a,'Actuals',header=TRUE) a<-table.element(a,'Forecasts',header=TRUE) a<-table.element(a,'Residuals',header=TRUE) a<-table.row.end(a) for (i in 1:length(result$Actuals)) { a<-table.row.start(a) a<-table.element(a,i,header=TRUE) a<-table.element(a,result$Actuals[i]) a<-table.element(a,result$Forecasts[i]) a<-table.element(a,result$Residuals[i]) a<-table.row.end(a) } a<-table.end(a) table.save(a,file="/var/wessaorg/rcomp/tmp/69odg1324406163.tab") } if (par2 != 'none') { a<-table.start() a<-table.row.start(a) a<-table.element(a,'Confusion Matrix (predicted in columns / actuals in rows)',par3+1,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'',1,TRUE) for (i in 1:par3) { a<-table.element(a,paste('C',i,sep=''),1,TRUE) } a<-table.row.end(a) for (i in 1:par3) { a<-table.row.start(a) a<-table.element(a,paste('C',i,sep=''),1,TRUE) for (j in 1:par3) { a<-table.element(a,myt[i,j]) } a<-table.row.end(a) } a<-table.end(a) table.save(a,file="/var/wessaorg/rcomp/tmp/7ggik1324406163.tab") } try(system("convert tmp/2pafv1324406163.ps tmp/2pafv1324406163.png",intern=TRUE)) try(system("convert tmp/3kvx81324406163.ps tmp/3kvx81324406163.png",intern=TRUE)) try(system("convert tmp/4ig241324406163.ps tmp/4ig241324406163.png",intern=TRUE))