R version 2.8.0 (2008-10-20)
Copyright (C) 2008 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.

  Natural language support but running in an English locale

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(105.31
+ ,1576.23
+ ,29.29
+ ,710.45
+ ,105.63
+ ,1546.37
+ ,28.99
+ ,720
+ ,106.02
+ ,1545.05
+ ,28.91
+ ,720
+ ,105.85
+ ,1552.34
+ ,29.29
+ ,720
+ ,106.57
+ ,1594.3
+ ,30.96
+ ,754.78
+ ,106.48
+ ,1605.78
+ ,30.57
+ ,802.73
+ ,106.60
+ ,1673.21
+ ,30.59
+ ,845.24
+ ,106.75
+ ,1612.94
+ ,31.39
+ ,893.91
+ ,106.69
+ ,1566.34
+ ,31.28
+ ,931.43
+ ,106.69
+ ,1530.17
+ ,31.1
+ ,940
+ ,106.93
+ ,1582.54
+ ,31.7
+ ,947.73
+ ,107.21
+ ,1702.16
+ ,32.57
+ ,960
+ ,107.88
+ ,1701.93
+ ,32.49
+ ,996.96
+ ,108.84
+ ,1811.15
+ ,32.46
+ ,1000
+ ,108.96
+ ,1924.2
+ ,32.3
+ ,1000
+ ,109.52
+ ,2034.25
+ ,32.97
+ ,1000
+ ,108.45
+ ,2011.13
+ ,32.9
+ ,1013.04
+ ,108.67
+ ,2013.04
+ ,32.93
+ ,1095.24
+ ,108.96
+ ,2151.67
+ ,33.72
+ ,1159.09
+ ,108.76
+ ,1902.09
+ ,33.33
+ ,1200
+ ,107.85
+ ,1944.01
+ ,33.44
+ ,1200
+ ,108.78
+ ,1916.67
+ ,33.89
+ ,1282.61
+ ,107.51
+ ,1967.31
+ ,34.34
+ ,1513.64
+ ,108.83
+ ,2119.88
+ ,33.56
+ ,1669.05
+ ,111.54
+ ,2216.38
+ ,32.67
+ ,1700
+ ,111.74
+ ,2522.83
+ ,32.57
+ ,1700
+ ,112.04
+ ,2647.64
+ ,33.23
+ ,1700
+ ,111.74
+ ,2631.23
+ ,32.85
+ ,1665.91
+ ,111.81
+ ,2693.41
+ ,32.61
+ ,1650
+ ,111.86
+ ,3021.76
+ ,32.57
+ ,1650
+ ,114.23
+ ,2953.67
+ ,32.98
+ ,1619.57
+ ,114.80
+ ,2796.8
+ ,31.33
+ ,1599.05
+ ,115.17
+ ,2672.05
+ ,29.8
+ ,1572.73
+ ,115.11
+ ,2251.23
+ ,28.06
+ ,1470
+ ,114.43
+ ,2046.08
+ ,25.47
+ ,1268
+ ,114.66
+ ,2420.04
+ ,24.65
+ ,1217.39
+ ,115.11
+ ,2608.89
+ ,23.94
+ ,1154.09
+ ,117.74
+ ,2660.47
+ ,23.89
+ ,984
+ ,118.18
+ ,2493.98
+ ,23.54
+ ,900
+ ,118.56
+ ,2541.7
+ ,24.28
+ ,900
+ ,117.63
+ ,2554.6
+ ,25.51
+ ,916.67
+ ,117.71
+ ,2699.61
+ ,27.03
+ ,957.73
+ ,117.46
+ ,2805.48
+ ,27.09
+ ,966.09
+ ,117.37
+ ,2956.66
+ ,27.3
+ ,980
+ ,117.34
+ ,3149.51
+ ,27.11
+ ,990.91
+ ,117.09
+ ,3372.5
+ ,26.39
+ ,1000.91
+ ,116.65
+ ,3379.33
+ ,27.54
+ ,1042.38
+ ,116.71
+ ,3517.54
+ ,26.85
+ ,1142.61
+ ,116.82
+ ,3527.34
+ ,26.82
+ ,1214.29
+ ,117.33
+ ,3281.06
+ ,25.9
+ ,1218
+ ,117.95
+ ,3089.65
+ ,24.96
+ ,1202.61
+ ,123.53
+ ,3222.76
+ ,25.4
+ ,1200
+ ,124.91
+ ,3165.76
+ ,24.38
+ ,1228.57
+ ,125.99
+ ,3232.43
+ ,24.73
+ ,1195.91
+ ,126.29
+ ,3229.54
+ ,25.43
+ ,1180
+ ,125.68
+ ,3071.74
+ ,26.04
+ ,1210.91
+ ,125.52
+ ,2850.17
+ ,25.59
+ ,1272.27)
+ ,dim=c(4
+ ,57)
+ ,dimnames=list(c('PC&S'
+ ,'PCacao'
+ ,'PSuiker'
+ ,'Pnoten')
+ ,1:57))
>  y <- array(NA,dim=c(4,57),dimnames=list(c('PC&S','PCacao','PSuiker','Pnoten'),1:57))
>  for (i in 1:dim(x)[1])
+  {
+  	for (j in 1:dim(x)[2])
+  	{
+  		y[i,j] <- as.numeric(x[i,j])
+  	}
+  }
> par4 = 'no'
> par3 = ''
> par2 = 'none'
> 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(party)
Loading required package: survival
Loading required package: splines
Loading required package: grid
Loading required package: modeltools
Loading required package: stats4
Loading required package: coin
Loading required package: mvtnorm
Loading required package: zoo

Attaching package: 'zoo'


	The following object(s) are masked from package:base :

	 as.Date.numeric 

Loading required package: sandwich
Loading required package: strucchange
Loading required package: vcd
Loading required package: MASS
Loading required package: colorspace
> library(Hmisc)

Attaching package: 'Hmisc'


	The following object(s) are masked from package:survival :

	 untangle.specials 


	The following object(s) are masked from package:base :

	 format.pval,
	 round.POSIXt,
	 trunc.POSIXt,
	 units 

> par1 <- as.numeric(par1)
> par3 <- as.numeric(par3)
> x <- data.frame(t(y))
> is.data.frame(x)
[1] TRUE
> x <- x[!is.na(x[,par1]),]
> k <- length(x[1,])
> n <- length(x[,1])
> colnames(x)[par1]
[1] "PC.S"
> x[,par1]
 [1] 105.31 105.63 106.02 105.85 106.57 106.48 106.60 106.75 106.69 106.69
[11] 106.93 107.21 107.88 108.84 108.96 109.52 108.45 108.67 108.96 108.76
[21] 107.85 108.78 107.51 108.83 111.54 111.74 112.04 111.74 111.81 111.86
[31] 114.23 114.80 115.17 115.11 114.43 114.66 115.11 117.74 118.18 118.56
[41] 117.63 117.71 117.46 117.37 117.34 117.09 116.65 116.71 116.82 117.33
[51] 117.95 123.53 124.91 125.99 126.29 125.68 125.52
> 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])

105.31 105.63 105.85 106.02 106.48 106.57  106.6 106.69 106.75 106.93 107.21 
     1      1      1      1      1      1      1      2      1      1      1 
107.51 107.85 107.88 108.45 108.67 108.76 108.78 108.83 108.84 108.96 109.52 
     1      1      1      1      1      1      1      1      1      2      1 
111.54 111.74 111.81 111.86 112.04 114.23 114.43 114.66  114.8 115.11 115.17 
     1      2      1      1      1      1      1      1      1      2      1 
116.65 116.71 116.82 117.09 117.33 117.34 117.37 117.46 117.63 117.71 117.74 
     1      1      1      1      1      1      1      1      1      1      1 
117.95 118.18 118.56 123.53 124.91 125.52 125.68 125.99 126.29 
     1      1      1      1      1      1      1      1      1 
> colnames(x)
[1] "PC.S"    "PCacao"  "PSuiker" "Pnoten" 
> colnames(x)[par1]
[1] "PC.S"
> x[,par1]
 [1] 105.31 105.63 106.02 105.85 106.57 106.48 106.60 106.75 106.69 106.69
[11] 106.93 107.21 107.88 108.84 108.96 109.52 108.45 108.67 108.96 108.76
[21] 107.85 108.78 107.51 108.83 111.54 111.74 112.04 111.74 111.81 111.86
[31] 114.23 114.80 115.17 115.11 114.43 114.66 115.11 117.74 118.18 118.56
[41] 117.63 117.71 117.46 117.37 117.34 117.09 116.65 116.71 116.82 117.33
[51] 117.95 123.53 124.91 125.99 126.29 125.68 125.52
> if (par2 == 'none') {
+ m <- ctree(as.formula(paste(colnames(x)[par1],' ~ .',sep='')),data = x)
+ }
> 
> #Note: the /var/www/html/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/html/freestat/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/www/html/freestat/rcomp/tmp/1yg7f1292936292.tab") 
+ }
+ }
> m

	 Conditional inference tree with 4 terminal nodes

Response:  PC.S 
Inputs:  PCacao, PSuiker, Pnoten 
Number of observations:  57 

1) PCacao <= 2216.38; criterion = 1, statistic = 40.16
  2) PCacao <= 1702.16; criterion = 1, statistic = 15.594
    3)*  weights = 13 
  2) PCacao > 1702.16
    4)*  weights = 13 
1) PCacao > 2216.38
  5) PSuiker <= 27.54; criterion = 0.999, statistic = 13.566
    6)*  weights = 22 
  5) PSuiker > 27.54
    7)*  weights = 9 
> postscript(file="/var/www/html/freestat/rcomp/tmp/2yg7f1292936292.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) 
> plot(m)
> dev.off()
null device 
          1 
> postscript(file="/var/www/html/freestat/rcomp/tmp/39p601292936292.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()
null device 
          1 
> 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)
+ }
   Actuals Forecasts   Residuals
1   105.31  106.5085 -1.19846154
2   105.63  106.5085 -0.87846154
3   106.02  106.5085 -0.48846154
4   105.85  106.5085 -0.65846154
5   106.57  106.5085  0.06153846
6   106.48  106.5085 -0.02846154
7   106.60  106.5085  0.09153846
8   106.75  106.5085  0.24153846
9   106.69  106.5085  0.18153846
10  106.69  106.5085  0.18153846
11  106.93  106.5085  0.42153846
12  107.21  106.5085  0.70153846
13  107.88  106.5085  1.37153846
14  108.84  109.3154 -0.47538462
15  108.96  109.3154 -0.35538462
16  109.52  109.3154  0.20461538
17  108.45  109.3154 -0.86538462
18  108.67  109.3154 -0.64538462
19  108.96  109.3154 -0.35538462
20  108.76  109.3154 -0.55538462
21  107.85  109.3154 -1.46538462
22  108.78  109.3154 -0.53538462
23  107.51  109.3154 -1.80538462
24  108.83  109.3154 -0.48538462
25  111.54  109.3154  2.22461538
26  111.74  113.1667 -1.42666667
27  112.04  113.1667 -1.12666667
28  111.74  113.1667 -1.42666667
29  111.81  113.1667 -1.35666667
30  111.86  113.1667 -1.30666667
31  114.23  113.1667  1.06333333
32  114.80  113.1667  1.63333333
33  115.17  113.1667  2.00333333
34  115.11  113.1667  1.94333333
35  114.43  109.3154  5.11461538
36  114.66  119.3741 -4.71409091
37  115.11  119.3741 -4.26409091
38  117.74  119.3741 -1.63409091
39  118.18  119.3741 -1.19409091
40  118.56  119.3741 -0.81409091
41  117.63  119.3741 -1.74409091
42  117.71  119.3741 -1.66409091
43  117.46  119.3741 -1.91409091
44  117.37  119.3741 -2.00409091
45  117.34  119.3741 -2.03409091
46  117.09  119.3741 -2.28409091
47  116.65  119.3741 -2.72409091
48  116.71  119.3741 -2.66409091
49  116.82  119.3741 -2.55409091
50  117.33  119.3741 -2.04409091
51  117.95  119.3741 -1.42409091
52  123.53  119.3741  4.15590909
53  124.91  119.3741  5.53590909
54  125.99  119.3741  6.61590909
55  126.29  119.3741  6.91590909
56  125.68  119.3741  6.30590909
57  125.52  119.3741  6.14590909
> if (par2 != 'none') {
+ print(cbind(as.factor(x[,par1]),predict(m)))
+ myt <- table(as.factor(x[,par1]),predict(m))
+ print(myt)
+ }
> postscript(file="/var/www/html/freestat/rcomp/tmp/4kgol1292936292.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()
null device 
          1 
> 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/www/html/freestat/rcomp/tmp/5nzmr1292936292.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/www/html/freestat/rcomp/tmp/68ilf1292936292.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/www/html/freestat/rcomp/tmp/7j9k01292936292.tab") 
+ }
> try(system("convert tmp/2yg7f1292936292.ps tmp/2yg7f1292936292.png",intern=TRUE))
character(0)
> try(system("convert tmp/39p601292936292.ps tmp/39p601292936292.png",intern=TRUE))
character(0)
> try(system("convert tmp/4kgol1292936292.ps tmp/4kgol1292936292.png",intern=TRUE))
character(0)
> 
> 
> proc.time()
   user  system elapsed 
  4.155   0.834  20.054