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 <- c(40.7819 + ,39.5915 + ,38.8859 + ,39.9068 + ,41.47 + ,41.5613 + ,41.6005 + ,41.4113 + ,41.84 + ,42.2892 + ,43.1521 + ,43.5998 + ,43.116 + ,42.4185 + ,42.3687 + ,42.2975 + ,42.8528 + ,43.535 + ,44.7265 + ,45.7293 + ,45.7585 + ,46.1685 + ,46.5075 + ,46.527 + ,46.601 + ,46.4607 + ,46.7135 + ,46.4113 + ,45.55 + ,44.6081 + ,44.4395 + ,44.9847 + ,45.7558 + ,45.3942 + ,45.697 + ,45.5664 + ,46.0205 + ,45.9195 + ,45.8005 + ,45.535 + ,45.4977 + ,45.5782 + ,45.7697 + ,45.2445 + ,45.0615 + ,45.2865 + ,44.791 + ,44.7625 + ,44.7644 + ,44.9973 + ,44.7265 + ,45.1465 + ,44.7465 + ,45.1795 + ,45.6515 + ,45.492 + ,45.2775 + ,45.2115 + ,45.411 + ,45.4005 + ,44.7692 + ,44.8913 + ,45.032 + ,44.879 + ,44.833 + ,44.8257 + ,44.7815 + ,44.479 + ,44.6317 + ,44.5043 + ,44.3217 + ,44.1005 + ,44.047 + ,43.6835 + ,43.7864 + ,44.1807 + ,43.9595 + ,43.937 + ,43.991 + ,43.865 + ,43.671 + ,43.93 + ,43.863 + ,43.7095 + ,43.9435 + ,43.736 + ,43.6295 + ,43.598 + ,43.8726 + ,43.8935 + ,43.5957 + ,43.7155 + ,43.528 + ,43.3415 + ,43.3374 + ,43.332 + ,43.3869 + ,43.5016 + ,43.4875 + ,43.6023 + ,43.3886 + ,43.3105 + ,43.4455 + ,43.5185 + ,43.5755 + ,43.6217 + ,43.644 + ,43.5789 + ,43.5215 + ,43.5033 + ,43.632 + ,43.263 + ,43.3717 + ,43.2745 + ,43.2647 + ,43.324 + ,43.4455 + ,43.4098 + ,43.41 + ,43.93 + ,43.8104 + ,43.54 + ,43.858 + ,43.8375 + ,43.881 + ,43.887 + ,43.8009 + ,43.7877 + ,43.811 + ,44.0625 + ,44.125 + ,44.52 + ,45.4005 + ,45.89 + ,45.189 + ,44.9035 + ,44.9351 + ,44.801 + ,43.98 + ,44.11 + ,44.2661 + ,44.361 + ,44.099 + ,43.8435 + ,43.8914 + ,44.217 + ,44.506 + ,44.54 + ,44.4465 + ,44.842 + ,44.8946 + ,44.951 + ,45.445 + ,45.0035 + ,45.769 + ,46.09 + ,45.412 + ,45.12 + ,45.48 + ,45.105 + ,45.056 + ,45.22 + ,45.39 + ,45.041 + ,44.9399 + ,44.9315 + ,45.1935 + ,45.3466 + ,45.4645 + ,45.5685 + ,45.3921 + ,45.34 + ,45.1308 + ,45.1005 + ,45.37 + ,45.2 + ,44.9614 + ,44.8015 + ,44.9152 + ,45.095 + ,44.9271 + ,44.6026 + ,44.5 + ,44.54 + ,44.5532 + ,44.407 + ,44.259 + ,44.1365 + ,44.112 + ,43.8814 + ,43.98 + ,43.7294 + ,43.9119 + ,43.955 + ,43.9 + ,43.7065 + ,43.6939 + ,43.6587 + ,43.5885 + ,43.8885 + ,43.8216 + ,43.751 + ,43.699 + ,43.7425 + ,43.639 + ,43.589 + ,43.606 + ,43.5325 + ,43.385 + ,43.3745 + ,43.236 + ,43.1957 + ,43.01 + ,43.1401 + ,43.0487 + ,43.1972 + ,43.2461 + ,43.0866 + ,43.0865 + ,43.0194 + ,43.08 + ,43.007 + ,42.9278 + ,42.9545 + ,42.7995 + ,42.9048 + ,42.9468 + ,43.08 + ,43.1274 + ,43.1625 + ,43.45 + ,43.831 + ,43.7769 + ,43.98 + ,43.92 + ,44.11 + ,44.03 + ,44.1582 + ,44.14 + ,45.07 + ,44.8737 + ,44.8505 + ,44.373 + ,44.075 + ,43.9725 + ,44.094 + ,44.191 + ,43.9685 + ,43.79 + ,43.6041 + ,43.1707 + ,42.71 + ,42.755 + ,43.3316 + ,43.5 + ,43.154 + ,43.16 + ,43.1 + ,42.85 + ,42.6175 + ,42.5 + ,42.6285 + ,42.6974 + ,43.04 + ,42.673 + ,42.5015 + ,42.538 + ,42.3735 + ,42.014 + ,41.8618 + ,42.1824 + ,42.605 + ,42.7345 + ,42.615 + ,42.465 + ,42.34 + ,42.251 + ,42.0475 + ,41.86 + ,41.685 + ,41.735 + ,41.706 + ,41.764 + ,41.58 + ,41.373 + ,41.088 + ,41.137 + ,41.1587 + ,41.185 + ,40.819 + ,40.633 + ,40.858 + ,40.794 + ,40.69 + ,40.595 + ,40.7305 + ,40.5471 + ,40.5145 + ,40.7 + ,40.7 + ,40.522 + ,40.6165 + ,40.3985 + ,40.2815 + ,40.245 + ,40.3055 + ,40.2696 + ,40.251 + ,40.127 + ,39.95 + ,39.675 + ,39.954 + ,39.8828 + ,39.62 + ,39.5415 + ,39.525 + ,39.8145 + ,39.6675 + ,39.695 + ,39.5985 + ,39.2735 + ,39.1435 + ,39.1742 + ,39.2025 + ,39.3946 + ,39.5025 + ,39.4845 + ,39.33 + ,39.295 + ,39.2675 + ,39.2535 + ,38.9845 + ,38.9285 + ,38.8592 + ,38.77 + ,38.79 + ,38.8205 + ,38.7577 + ,38.839 + ,38.78 + ,38.54 + ,38.511 + ,38.615 + ,38.898 + ,38.8691 + ,38.384 + ,38.0277 + ,37.72 + ,37.7325 + ,37.626 + ,37.603 + ,37.78 + ,38.559 + ,39.0459 + ,38.45 + ,38.505 + ,38.2885 + ,37.795 + ,37.92 + ,38.034 + ,38.029 + ,38.063 + ,37.9828 + ,37.745 + ,37.969 + ,38.007 + ,38.0615 + ,38.0912 + ,38.091 + ,38.431 + ,38.48 + ,38.35 + ,38.214 + ,38.384 + ,38.1375 + ,38.0075 + ,38.0524 + ,38.235 + ,38.31 + ,38.2615 + ,38.13 + ,38.282 + ,38.581 + ,39.0801 + ,39.0387 + ,39.1015 + ,39.1503 + ,39.14 + ,39.0275 + ,38.7665 + ,38.691 + ,38.849 + ,39.1644 + ,39.4907 + ,39.5095 + ,39.2795 + ,39.0437 + ,39.1355 + ,39.143 + ,39.185 + ,39.355 + ,39.297 + ,39.4514 + ,39.4173 + ,39.4305 + ,39.384 + ,39.3261 + ,39.301 + ,39.35 + ,39.64 + ,39.4723 + ,39.3685 + ,39.1906 + ,39.1183 + ,39.1325 + ,39.1144 + ,39.1614 + ,39.0908 + ,38.9199 + ,38.913 + ,38.9655 + ,39.029 + ,39.089 + ,39.07 + ,39.0046 + ,39.1038 + ,39.3572 + ,39.388 + ,39.382 + ,39.4398 + ,39.2537 + ,39.2301 + ,39.2763 + ,39.282 + ,39.3325 + ,39.557 + ,40.1 + ,40.5875 + ,40.485 + ,40.55 + ,40.7955 + ,41.456 + ,41.3557 + ,41.374 + ,41.2235 + ,41.15 + ,41.3725 + ,41.6923 + ,41.8 + ,41.8045 + ,41.64 + ,41.36 + ,41.5745 + ,41.593 + ,41.575 + ,41.68 + ,42.0055 + ,42.3188 + ,42.565 + ,42.3575 + ,42.29 + ,42.695 + ,43.0028 + ,42.4507 + ,42.4705 + ,42.2875 + ,42.3172 + ,42.55 + ,42.7523 + ,42.8993 + ,43.1555 + ,43.1885 + ,43.43 + ,43.31 + ,42.815 + ,42.7017 + ,42.28 + ,41.922 + ,42.17 + ,42.1962 + ,42.3215 + ,42.3173 + ,42.391 + ,42.463 + ,42.4125 + ,42.304 + ,41.813 + ,41.651 + ,41.539 + ,41.1575 + ,40.9545) > par9 = '1' > par8 = '2' > par7 = '1' > par6 = '3' > par5 = '1' > par4 = '0' > par3 = '1' > par2 = '-1.8' > par1 = 'FALSE' > library(lattice) > if (par1 == 'TRUE') par1 <- TRUE > if (par1 == 'FALSE') par1 <- FALSE > par2 <- as.numeric(par2) #Box-Cox lambda transformation parameter > par3 <- as.numeric(par3) #degree of non-seasonal differencing > par4 <- as.numeric(par4) #degree of seasonal differencing > par5 <- 5 #seasonal period > par6 <- as.numeric(par6) #degree (p) of the non-seasonal AR(p) polynomial > par7 <- as.numeric(par7) #degree (q) of the non-seasonal MA(q) polynomial > par8 <- as.numeric(par8) #degree (P) of the seasonal AR(P) polynomial > par9 <- as.numeric(par9) #degree (Q) of the seasonal MA(Q) polynomial > armaGR <- function(arima.out, names, n){ + try1 <- arima.out$coef + try2 <- sqrt(diag(arima.out$var.coef)) + try.data.frame <- data.frame(matrix(NA,ncol=4,nrow=length(names))) + dimnames(try.data.frame) <- list(names,c('coef','std','tstat','pv')) + try.data.frame[,1] <- try1 + for(i in 1:length(try2)) try.data.frame[which(rownames(try.data.frame)==names(try2)[i]),2] <- try2[i] + try.data.frame[,3] <- try.data.frame[,1] / try.data.frame[,2] + try.data.frame[,4] <- round((1-pt(abs(try.data.frame[,3]),df=n-(length(try2)+1)))*2,5) + vector <- rep(NA,length(names)) + vector[is.na(try.data.frame[,4])] <- 0 + maxi <- which.max(try.data.frame[,4]) + continue <- max(try.data.frame[,4],na.rm=TRUE) > .05 + vector[maxi] <- 0 + list(summary=try.data.frame,next.vector=vector,continue=continue) + } > arimaSelect <- function(series, order=c(13,0,0), seasonal=list(order=c(2,0,0),period=12), include.mean=F){ + nrc <- order[1]+order[3]+seasonal$order[1]+seasonal$order[3] + coeff <- matrix(NA, nrow=nrc*2, ncol=nrc) + pval <- matrix(NA, nrow=nrc*2, ncol=nrc) + mylist <- rep(list(NULL), nrc) + names <- NULL + if(order[1] > 0) names <- paste('ar',1:order[1],sep='') + if(order[3] > 0) names <- c( names , paste('ma',1:order[3],sep='') ) + if(seasonal$order[1] > 0) names <- c(names, paste('sar',1:seasonal$order[1],sep='')) + if(seasonal$order[3] > 0) names <- c(names, paste('sma',1:seasonal$order[3],sep='')) + arima.out <- arima(series, order=order, seasonal=seasonal, include.mean=include.mean, method='ML') + mylist[[1]] <- arima.out + last.arma <- armaGR(arima.out, names, length(series)) + mystop <- FALSE + i <- 1 + coeff[i,] <- last.arma[[1]][,1] + pval [i,] <- last.arma[[1]][,4] + i <- 2 + aic <- arima.out$aic + while(!mystop){ + mylist[[i]] <- arima.out + arima.out <- arima(series, order=order, seasonal=seasonal, include.mean=include.mean, method='ML', fixed=last.arma$next.vector) + aic <- c(aic, arima.out$aic) + last.arma <- armaGR(arima.out, names, length(series)) + mystop <- !last.arma$continue + coeff[i,] <- last.arma[[1]][,1] + pval [i,] <- last.arma[[1]][,4] + i <- i+1 + } + list(coeff, pval, mylist, aic=aic) + } > arimaSelectplot <- function(arimaSelect.out,noms,choix){ + noms <- names(arimaSelect.out[[3]][[1]]$coef) + coeff <- arimaSelect.out[[1]] + k <- min(which(is.na(coeff[,1])))-1 + coeff <- coeff[1:k,] + pval <- arimaSelect.out[[2]][1:k,] + aic <- arimaSelect.out$aic[1:k] + coeff[coeff==0] <- NA + n <- ncol(coeff) + if(missing(choix)) choix <- k + layout(matrix(c(1,1,1,2, + 3,3,3,2, + 3,3,3,4, + 5,6,7,7),nr=4), + widths=c(10,35,45,15), + heights=c(30,30,15,15)) + couleurs <- rainbow(75)[1:50]#(50) + ticks <- pretty(coeff) + par(mar=c(1,1,3,1)) + plot(aic,k:1-.5,type='o',pch=21,bg='blue',cex=2,axes=F,lty=2,xpd=NA) + points(aic[choix],k-choix+.5,pch=21,cex=4,bg=2,xpd=NA) + title('aic',line=2) + par(mar=c(3,0,0,0)) + plot(0,axes=F,xlab='',ylab='',xlim=range(ticks),ylim=c(.1,1)) + rect(xleft = min(ticks) + (0:49)/50*(max(ticks)-min(ticks)), + xright = min(ticks) + (1:50)/50*(max(ticks)-min(ticks)), + ytop = rep(1,50), + ybottom= rep(0,50),col=couleurs,border=NA) + axis(1,ticks) + rect(xleft=min(ticks),xright=max(ticks),ytop=1,ybottom=0) + text(mean(coeff,na.rm=T),.5,'coefficients',cex=2,font=2) + par(mar=c(1,1,3,1)) + image(1:n,1:k,t(coeff[k:1,]),axes=F,col=couleurs,zlim=range(ticks)) + for(i in 1:n) for(j in 1:k) if(!is.na(coeff[j,i])) { + if(pval[j,i]<.01) symb = 'green' + else if( (pval[j,i]<.05) & (pval[j,i]>=.01)) symb = 'orange' + else if( (pval[j,i]<.1) & (pval[j,i]>=.05)) symb = 'red' + else symb = 'black' + polygon(c(i+.5 ,i+.2 ,i+.5 ,i+.5), + c(k-j+0.5,k-j+0.5,k-j+0.8,k-j+0.5), + col=symb) + if(j==choix) { + rect(xleft=i-.5, + xright=i+.5, + ybottom=k-j+1.5, + ytop=k-j+.5, + lwd=4) + text(i, + k-j+1, + round(coeff[j,i],2), + cex=1.2, + font=2) + } + else{ + rect(xleft=i-.5,xright=i+.5,ybottom=k-j+1.5,ytop=k-j+.5) + text(i,k-j+1,round(coeff[j,i],2),cex=1.2,font=1) + } + } + axis(3,1:n,noms) + par(mar=c(0.5,0,0,0.5)) + plot(0,axes=F,xlab='',ylab='',type='n',xlim=c(0,8),ylim=c(-.2,.8)) + cols <- c('green','orange','red','black') + niv <- c('0','0.01','0.05','0.1') + for(i in 0:3){ + polygon(c(1+2*i ,1+2*i ,1+2*i-.5 ,1+2*i), + c(.4 ,.7 , .4 , .4), + col=cols[i+1]) + text(2*i,0.5,niv[i+1],cex=1.5) + } + text(8,.5,1,cex=1.5) + text(4,0,'p-value',cex=2) + box() + residus <- arimaSelect.out[[3]][[choix]]$res + par(mar=c(1,2,4,1)) + acf(residus,main='') + title('acf',line=.5) + par(mar=c(1,2,4,1)) + pacf(residus,main='') + title('pacf',line=.5) + par(mar=c(2,2,4,1)) + qqnorm(residus,main='') + title('qq-norm',line=.5) + qqline(residus) + residus + } > if (par2 == 0) x <- log(x) > if (par2 != 0) x <- x^par2 > (selection <- arimaSelect(x, order=c(par6,par3,par7), seasonal=list(order=c(par8,par4,par9), period=par5))) [[1]] [,1] [,2] [,3] [,4] [,5] [,6] [1,] 0.1415161 -0.14615919 -0.1038974 0.1752222 0.4182325 0.09580098 [2,] 0.0000000 -0.10129709 -0.1329508 0.3170740 0.4550993 0.09472927 [3,] 0.0000000 -0.09318170 -0.1298711 0.3192973 0.0000000 0.06444705 [4,] 0.0000000 -0.08905061 -0.1296174 0.3186951 0.0000000 0.06165557 [5,] 0.0000000 -0.08714675 -0.1230937 0.3178127 0.0000000 0.00000000 [6,] 0.0000000 0.00000000 -0.1249547 0.3481668 0.0000000 0.00000000 [7,] NA NA NA NA NA NA [8,] NA NA NA NA NA NA [9,] NA NA NA NA NA NA [10,] NA NA NA NA NA NA [11,] NA NA NA NA NA NA [12,] NA NA NA NA NA NA [13,] NA NA NA NA NA NA [14,] NA NA NA NA NA NA [,7] [1,] -0.46620426 [2,] -0.50133048 [3,] -0.04823813 [4,] 0.00000000 [5,] 0.00000000 [6,] 0.00000000 [7,] NA [8,] NA [9,] NA [10,] NA [11,] NA [12,] NA [13,] NA [14,] NA [[2]] [,1] [,2] [,3] [,4] [,5] [,6] [,7] [1,] 0.87037 0.60443 0.58266 0.84016 0.36064 0.05985 0.30999 [2,] NA 0.04254 0.00535 0.00000 0.34655 0.06007 0.30315 [3,] NA 0.05808 0.00617 0.00000 NA 0.19719 0.32531 [4,] NA 0.06956 0.00635 0.00000 NA 0.21437 NA [5,] NA 0.07579 0.00913 0.00000 NA NA NA [6,] NA NA 0.00791 0.00000 NA NA NA [7,] NA NA NA NA NA NA NA [8,] NA NA NA NA NA NA NA [9,] NA NA NA NA NA NA NA [10,] NA NA NA NA NA NA NA [11,] NA NA NA NA NA NA NA [12,] NA NA NA NA NA NA NA [13,] NA NA NA NA NA NA NA [14,] NA NA NA NA NA NA NA [[3]] [[3]][[1]] Call: arima(x = series, order = order, seasonal = seasonal, include.mean = include.mean, method = "ML") Coefficients: ar1 ar2 ar3 ma1 sar1 sar2 sma1 0.1415 -0.1462 -0.1039 0.1752 0.4182 0.0958 -0.4662 s.e. 0.8668 0.2820 0.1889 0.8683 0.4571 0.0508 0.4587 sigma^2 estimated as 1.686e-10: log likelihood = 4817.87, aic = -9619.75 [[3]][[2]] Call: arima(x = series, order = order, seasonal = seasonal, include.mean = include.mean, method = "ML") Coefficients: ar1 ar2 ar3 ma1 sar1 sar2 sma1 0.1415 -0.1462 -0.1039 0.1752 0.4182 0.0958 -0.4662 s.e. 0.8668 0.2820 0.1889 0.8683 0.4571 0.0508 0.4587 sigma^2 estimated as 1.686e-10: log likelihood = 4817.87, aic = -9619.75 [[3]][[3]] Call: arima(x = series, order = order, seasonal = seasonal, include.mean = include.mean, fixed = last.arma$next.vector, method = "ML") Coefficients: ar1 ar2 ar3 ma1 sar1 sar2 sma1 0 -0.1013 -0.1330 0.3171 0.4551 0.0947 -0.5013 s.e. 0 0.0498 0.0475 0.0458 0.4830 0.0503 0.4864 sigma^2 estimated as 1.686e-10: log likelihood = 4817.88, aic = -9621.76 [[3]][[4]] Call: arima(x = series, order = order, seasonal = seasonal, include.mean = include.mean, fixed = last.arma$next.vector, method = "ML") Coefficients: ar1 ar2 ar3 ma1 sar1 sar2 sma1 0 -0.0932 -0.1299 0.3193 0 0.0644 -0.0482 s.e. 0 0.0491 0.0472 0.0457 0 0.0499 0.0490 sigma^2 estimated as 1.69e-10: log likelihood = 4817.36, aic = -9622.73 [[3]][[5]] Call: arima(x = series, order = order, seasonal = seasonal, include.mean = include.mean, fixed = last.arma$next.vector, method = "ML") Coefficients: ar1 ar2 ar3 ma1 sar1 sar2 sma1 0 -0.0891 -0.1296 0.3187 0 0.0617 0 s.e. 0 0.0490 0.0473 0.0458 0 0.0496 0 sigma^2 estimated as 1.693e-10: log likelihood = 4816.88, aic = -9623.76 [[3]][[6]] Call: arima(x = series, order = order, seasonal = seasonal, include.mean = include.mean, fixed = last.arma$next.vector, method = "ML") Coefficients: ar1 ar2 ar3 ma1 sar1 sar2 sma1 0 -0.0871 -0.1231 0.3178 0 0 0 s.e. 0 0.0490 0.0470 0.0457 0 0 0 sigma^2 estimated as 1.699e-10: log likelihood = 4816.11, aic = -9624.22 [[3]][[7]] NULL $aic [1] -9619.749 -9621.757 -9622.730 -9623.762 -9624.222 -9623.123 Warning messages: 1: In arima(series, order = order, seasonal = seasonal, include.mean = include.mean, : some AR parameters were fixed: setting transform.pars = FALSE 2: In arima(series, order = order, seasonal = seasonal, include.mean = include.mean, : some AR parameters were fixed: setting transform.pars = FALSE 3: In arima(series, order = order, seasonal = seasonal, include.mean = include.mean, : some AR parameters were fixed: setting transform.pars = FALSE 4: In arima(series, order = order, seasonal = seasonal, include.mean = include.mean, : some AR parameters were fixed: setting transform.pars = FALSE 5: In arima(series, order = order, seasonal = seasonal, include.mean = include.mean, : some AR parameters were fixed: setting transform.pars = FALSE > postscript(file="/var/www/html/freestat/rcomp/tmp/1lw8l1293574286.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > resid <- arimaSelectplot(selection) > dev.off() null device 1 > resid Time Series: Start = 1 End = 491 Frequency = 1 [1] 1.262288e-06 6.493993e-05 2.442521e-05 -6.245155e-05 -5.556105e-05 [6] 1.274528e-05 -2.147786e-05 5.641940e-06 -2.512361e-05 -1.434474e-05 [11] -3.839216e-05 -1.356333e-05 2.051390e-05 2.047595e-05 -4.622975e-06 [16] 1.080241e-05 -2.645645e-05 -2.333706e-05 -4.777279e-05 -3.283948e-05 [21] 6.339331e-07 -2.674222e-05 -9.960899e-06 8.447171e-07 -6.275419e-06 [26] 5.706285e-06 -1.185703e-05 1.550561e-05 2.919085e-05 3.018971e-05 [31] 2.170712e-06 -1.649089e-05 -2.111539e-05 2.032870e-05 -2.451063e-05 [36] 1.046433e-05 -2.087706e-05 9.597877e-06 7.828787e-07 8.626078e-06 [41] -3.022119e-07 -1.670599e-06 -5.778647e-06 2.326549e-05 -8.098763e-07 [46] -8.236807e-06 2.684312e-05 -7.185489e-06 2.865292e-06 -8.154312e-06 [51] 1.428725e-05 -2.325002e-05 2.414894e-05 -2.616327e-05 -1.185747e-05 [56] 1.076983e-05 1.480527e-06 4.513338e-07 -6.846335e-06 3.938929e-06 [61] 2.492730e-05 -1.412228e-05 8.913413e-07 9.020456e-06 -2.066019e-06 [66] 7.998918e-07 2.608032e-06 1.253313e-05 -1.042026e-05 1.022144e-05 [71] 5.796626e-06 7.653273e-06 1.345748e-06 1.792925e-05 -8.992834e-06 [76] -1.318637e-05 1.573158e-05 -6.112534e-06 -1.819202e-06 7.588945e-06 [81] 6.369126e-06 -1.364253e-05 8.848090e-06 4.258646e-06 -1.319999e-05 [86] 1.462184e-05 1.820357e-07 9.074352e-07 -1.127869e-05 3.365400e-06 [91] 1.164703e-05 -1.084857e-05 1.317106e-05 5.721455e-06 -1.550691e-06 [96] 2.570671e-06 -2.304460e-06 -4.575143e-06 1.916714e-06 -6.716248e-06 [101] 1.147391e-05 -3.622546e-07 -6.001618e-06 5.186975e-08 -2.759292e-06 [106] -2.333141e-06 -9.354355e-07 2.791811e-06 1.421117e-06 5.294543e-07 [111] -5.525653e-06 1.937649e-05 -1.168837e-05 9.053956e-06 -7.410765e-07 [116] -2.790258e-06 -4.202644e-06 2.817050e-06 -1.744734e-06 -2.392246e-05 [121] 1.323900e-05 6.141152e-06 -1.902422e-05 8.730706e-06 -4.493497e-06 [126] -5.602279e-07 4.037311e-06 -9.477924e-07 -4.539220e-07 -1.069473e-05 [131] 5.821181e-07 -1.873861e-05 -3.302688e-05 -1.126200e-05 2.685758e-05 [136] -2.825132e-06 -3.950868e-07 1.042339e-05 3.413260e-05 -1.635335e-05 [141] 2.098053e-06 -9.245763e-07 1.059953e-05 6.940134e-06 -3.883705e-06 [146] -1.094456e-05 -8.014291e-06 -4.759856e-07 1.325439e-06 -1.922009e-05 [151] 4.038754e-06 -4.668518e-06 -2.145611e-05 2.474947e-05 -4.157372e-05 [156] -5.442702e-07 2.706404e-05 -1.462248e-06 -1.367303e-05 2.434989e-05 [161] -5.483775e-06 -5.616554e-06 -3.173326e-06 1.522785e-05 -2.022506e-06 [166] 1.400286e-06 -9.317520e-06 -2.848687e-06 -4.875443e-06 -4.624412e-06 [171] 7.504276e-06 -1.200705e-06 9.207965e-06 -5.748938e-07 -1.002542e-05 [176] 1.143444e-05 5.584614e-06 4.270977e-06 -4.464265e-06 -4.353575e-06 [181] 8.895770e-06 9.840636e-06 1.018018e-06 2.038354e-08 1.522496e-06 [186] 6.301450e-06 4.259477e-06 4.580930e-06 9.940396e-07 1.135101e-05 [191] -7.302057e-06 1.475318e-05 -1.211323e-05 2.343971e-06 2.421506e-06 [196] 6.862559e-06 -1.625616e-06 3.210652e-06 3.355233e-06 -1.459399e-05 [201] 8.160838e-06 -1.688041e-07 1.010978e-06 -1.660543e-06 5.886962e-06 [206] 5.569279e-07 -7.938035e-07 4.441669e-06 5.681865e-06 -1.113816e-06 [211] 7.894695e-06 2.893395e-07 9.397795e-06 -8.236497e-06 7.988384e-06 [216] -9.066554e-06 1.808631e-07 7.442250e-06 -3.433312e-06 4.678743e-06 [221] -3.455128e-06 4.875265e-06 2.405598e-06 -2.104439e-06 8.936961e-06 [226] -7.590337e-06 8.806395e-07 -6.193827e-06 -1.098422e-06 -2.127787e-06 [231] -1.384647e-05 -1.358321e-05 5.398831e-06 -1.413161e-05 5.255783e-06 [236] -1.071968e-05 6.095571e-06 -8.088409e-06 2.642824e-06 -4.122271e-05 [241] 2.078225e-05 -9.029343e-06 1.932314e-05 8.201767e-06 3.928181e-06 [246] -3.009023e-06 -1.342667e-06 1.048922e-05 3.718804e-06 7.690800e-06 [251] 1.979229e-05 1.766819e-05 -4.994708e-06 -2.159871e-05 1.533295e-06 [256] 1.310666e-05 -8.536978e-06 6.019805e-06 1.209351e-05 7.734755e-06 [261] 4.751390e-06 -5.390185e-06 2.348536e-07 -1.651190e-05 2.196638e-05 [266] -3.805740e-07 -2.180767e-06 1.180992e-05 1.533999e-05 3.457020e-06 [271] -1.492898e-05 -1.345466e-05 -2.546466e-06 2.810158e-06 3.378351e-06 [276] 4.899245e-06 4.289389e-06 1.043284e-05 7.492979e-06 8.187044e-06 [281] -3.104067e-06 4.484779e-06 -3.560099e-06 1.059165e-05 7.556707e-06 [286] 1.346628e-05 -4.802220e-06 3.043180e-06 -7.308852e-07 1.988798e-05 [291] 3.812625e-06 -1.220399e-05 1.081265e-05 2.559092e-06 3.293369e-06 [296] -7.714123e-06 1.395163e-05 -2.591564e-06 -9.673974e-06 4.505770e-06 [301] 7.918386e-06 -9.146734e-06 1.615899e-05 2.354522e-06 1.780484e-06 [306] -1.947363e-06 3.699439e-06 -1.461758e-07 6.994860e-06 8.508056e-06 [311] 1.444481e-05 -1.942527e-05 1.309173e-05 1.213752e-05 -7.814810e-07 [316] 3.140509e-06 -1.607534e-05 1.458660e-05 -7.685257e-06 6.883143e-06 [321] 1.864369e-05 2.464625e-06 -2.499739e-07 1.469804e-06 -1.151147e-05 [326] -3.328975e-06 9.034926e-07 7.151665e-06 -8.284401e-07 2.925794e-06 [331] 1.292422e-06 1.684379e-05 -1.521289e-06 6.466850e-06 6.034352e-06 [336] -2.379899e-06 -1.564173e-07 4.666214e-06 -7.019280e-06 6.117835e-06 [341] 1.362686e-05 -2.747915e-06 -4.085325e-06 -1.480955e-05 6.189359e-06 [346] 2.708089e-05 1.314624e-05 1.987999e-05 -1.234077e-06 1.256872e-05 [351] 1.329427e-07 -1.179966e-05 -4.744807e-05 -1.697853e-05 3.764633e-05 [356] -2.471997e-05 2.163642e-05 3.093496e-05 -1.760849e-05 2.511476e-06 [361] 2.905471e-06 -4.955796e-06 6.086596e-06 1.421705e-05 -1.969936e-05 [366] 5.768761e-06 -4.856581e-06 -2.579086e-06 1.935772e-07 -2.333965e-05 [371] 3.953973e-06 5.334050e-06 4.281809e-06 -1.230348e-05 2.217306e-05 [376] 1.862663e-06 -3.593103e-06 -8.340197e-06 -1.525426e-06 2.266034e-06 [381] 6.136165e-06 -1.244958e-05 -1.455058e-05 -2.707929e-05 8.243293e-06 [386] -1.176649e-05 -3.013920e-06 1.578034e-06 5.807282e-06 1.443645e-05 [391] 9.674455e-07 -8.134052e-06 -1.485245e-05 -1.569095e-05 8.582776e-07 [396] 9.641644e-06 9.075420e-06 -7.553867e-06 4.952298e-06 -2.885058e-06 [401] -1.036105e-05 6.584415e-06 -1.282411e-05 5.178756e-06 -2.841207e-06 [406] 2.768297e-06 2.870767e-06 7.862546e-07 -2.610851e-06 -1.627162e-05 [411] 1.527543e-05 -4.071920e-07 9.850556e-06 3.190963e-06 -1.603168e-07 [416] 2.932378e-06 -3.391965e-06 5.484798e-06 8.926481e-06 -2.375394e-06 [421] -1.090577e-06 -2.297554e-06 -3.278283e-06 1.475409e-06 2.826939e-06 [426] -7.497844e-06 -1.284904e-05 2.153241e-06 -2.455246e-06 -4.865764e-06 [431] 1.280791e-05 -2.869501e-06 -1.389947e-06 1.627554e-06 -3.705812e-06 [436] -1.293733e-05 -2.853058e-05 -2.050962e-05 7.818858e-06 -1.259194e-05 [441] -1.271915e-05 -3.151227e-05 1.371559e-05 -1.016795e-05 7.367699e-06 [446] 2.215559e-06 -1.211680e-05 -1.173647e-05 -2.445578e-06 -2.408996e-06 [451] 6.790018e-06 1.198592e-05 -1.449468e-05 5.983631e-06 -1.173484e-07 [456] -6.975052e-06 -1.469775e-05 -1.159561e-05 -1.073253e-05 1.027324e-05 [461] -2.901931e-06 -1.980060e-05 -7.072781e-06 2.791223e-05 -1.363321e-05 [466] 1.400371e-05 -2.706414e-06 -1.006900e-05 -5.744891e-06 -6.511396e-06 [471] -1.251051e-05 5.622765e-07 -1.349866e-05 8.268866e-06 1.987059e-05 [476] -1.695321e-06 2.423976e-05 1.393857e-05 -1.460936e-05 7.483219e-06 [481] -7.556976e-06 9.354574e-07 -4.704901e-06 -2.855111e-06 3.128339e-06 [486] 3.674167e-06 2.370476e-05 1.711932e-06 8.216433e-06 2.166569e-05 [491] 5.770437e-06 > postscript(file="/var/www/html/freestat/rcomp/tmp/2lw8l1293574286.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > acf(resid,length(resid)/2, main='Residual Autocorrelation Function') > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/3w5p61293574286.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > pacf(resid,length(resid)/2, main='Residual Partial Autocorrelation Function') > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/4w5p61293574286.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > cpgram(resid, main='Residual Cumulative Periodogram') > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/5w5p61293574286.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > hist(resid, main='Residual Histogram', xlab='values of Residuals') > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/6w5p61293574286.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > densityplot(~resid,col='black',main='Residual Density Plot', xlab='values of Residuals') > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/7w5p61293574286.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > qqnorm(resid, main='Residual Normal Q-Q Plot') > qqline(resid) > dev.off() null device 1 > ncols <- length(selection[[1]][1,]) > nrows <- length(selection[[2]][,1])-1 > > #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") > > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'ARIMA Parameter Estimation and Backward Selection', ncols+1,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Iteration', header=TRUE) > for (i in 1:ncols) { + a<-table.element(a,names(selection[[3]][[1]]$coef)[i],header=TRUE) + } > a<-table.row.end(a) > for (j in 1:nrows) { + a<-table.row.start(a) + mydum <- 'Estimates (' + mydum <- paste(mydum,j) + mydum <- paste(mydum,')') + a<-table.element(a,mydum, header=TRUE) + for (i in 1:ncols) { + a<-table.element(a,round(selection[[1]][j,i],4)) + } + a<-table.row.end(a) + a<-table.row.start(a) + a<-table.element(a,'(p-val)', header=TRUE) + for (i in 1:ncols) { + mydum <- '(' + mydum <- paste(mydum,round(selection[[2]][j,i],4),sep='') + mydum <- paste(mydum,')') + a<-table.element(a,mydum) + } + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/freestat/rcomp/tmp/83on01293574286.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Estimated ARIMA Residuals', 1,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Value', 1,TRUE) > a<-table.row.end(a) > for (i in (par4*par5+par3):length(resid)) { + a<-table.row.start(a) + a<-table.element(a,resid[i]) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/freestat/rcomp/tmp/9vxml1293574286.tab") > > try(system("convert tmp/1lw8l1293574286.ps tmp/1lw8l1293574286.png",intern=TRUE)) character(0) > try(system("convert tmp/2lw8l1293574286.ps tmp/2lw8l1293574286.png",intern=TRUE)) character(0) > try(system("convert tmp/3w5p61293574286.ps tmp/3w5p61293574286.png",intern=TRUE)) character(0) > try(system("convert tmp/4w5p61293574286.ps tmp/4w5p61293574286.png",intern=TRUE)) character(0) > try(system("convert tmp/5w5p61293574286.ps tmp/5w5p61293574286.png",intern=TRUE)) character(0) > try(system("convert tmp/6w5p61293574286.ps tmp/6w5p61293574286.png",intern=TRUE)) character(0) > try(system("convert tmp/7w5p61293574286.ps tmp/7w5p61293574286.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 5.469 1.650 5.678