R version 2.15.2 (2012-10-26) -- "Trick or Treat" Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i686-pc-linux-gnu (32-bit) 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(1 + ,1 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,1 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,1 + ,1 + ,0 + ,1 + ,1 + ,0 + ,1 + ,1 + ,1 + ,1 + ,1 + ,1 + ,1 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,1 + ,1 + ,1 + ,1 + ,1 + ,1 + ,0 + ,0 + ,0 + ,1 + ,0 + ,1 + ,0 + ,1 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,1 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,1 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,0 + ,1 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,1 + ,1 + ,1 + ,1 + ,1 + ,1 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,1 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,1 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,1 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,1 + ,1 + ,0 + ,1 + ,1 + ,1 + ,0 + ,1 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,0 + ,1 + ,0 + ,0 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0 + ,0 + ,1 + ,1 + ,1 + ,0 + ,0 + ,0 + ,0 + ,0) + ,dim=c(6 + ,86) + ,dimnames=list(c('UseLimit' + ,'T40' + ,'Used' + ,'CorrectAnalysis' + ,'Useful' + ,'Outcome') + ,1:86)) > y <- array(NA,dim=c(6,86),dimnames=list(c('UseLimit','T40','Used','CorrectAnalysis','Useful','Outcome'),1:86)) > for (i in 1:dim(x)[1]) + { + for (j in 1:dim(x)[2]) + { + y[i,j] <- as.numeric(x[i,j]) + } + } > par3 = 'No Linear Trend' > par2 = 'Do not include Seasonal Dummies' > par1 = '6' > par3 <- 'No Linear Trend' > par2 <- 'Do not include Seasonal Dummies' > par1 <- '6' > #'GNU S' R Code compiled by R2WASP v. 1.0.44 () > #Author: Prof. Dr. P. Wessa > #To cite this work: AUTHOR(S), (YEAR), 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: Office for Research, Development, and Education > #Technical description: Write here your technical program description (don't use hard returns!) > library(lattice) > library(lmtest) Loading required package: zoo Attaching package: 'zoo' The following object(s) are masked from 'package:base': as.Date, as.Date.numeric > n25 <- 25 #minimum number of obs. for Goldfeld-Quandt test > par1 <- as.numeric(par1) > 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] > x <- x1 > if (par3 == 'First Differences'){ + x2 <- array(0, dim=c(n-1,k), dimnames=list(1:(n-1), paste('(1-B)',colnames(x),sep=''))) + for (i in 1:n-1) { + for (j in 1:k) { + x2[i,j] <- x[i+1,j] - x[i,j] + } + } + x <- x2 + } > if (par2 == 'Include Monthly Dummies'){ + x2 <- array(0, dim=c(n,11), dimnames=list(1:n, paste('M', seq(1:11), sep =''))) + for (i in 1:11){ + x2[seq(i,n,12),i] <- 1 + } + x <- cbind(x, x2) + } > if (par2 == 'Include Quarterly Dummies'){ + x2 <- array(0, dim=c(n,3), dimnames=list(1:n, paste('Q', seq(1:3), sep =''))) + for (i in 1:3){ + x2[seq(i,n,4),i] <- 1 + } + x <- cbind(x, x2) + } > k <- length(x[1,]) > if (par3 == 'Linear Trend'){ + x <- cbind(x, c(1:n)) + colnames(x)[k+1] <- 't' + } > x Outcome UseLimit T40 Used CorrectAnalysis Useful 1 1 1 1 0 0 0 2 0 0 0 0 0 0 3 0 0 0 0 0 0 4 0 0 0 0 0 0 5 0 0 0 0 0 0 6 1 1 0 0 0 1 7 0 0 0 0 0 0 8 0 0 1 0 0 0 9 1 0 0 0 0 0 10 0 1 0 0 0 0 11 0 1 1 0 0 0 12 0 0 0 0 0 0 13 0 0 0 1 0 1 14 0 1 1 0 0 0 15 1 0 0 1 0 1 16 1 0 1 1 0 1 17 0 1 1 1 1 1 18 0 1 1 0 0 0 19 1 0 0 0 0 0 20 1 0 1 1 1 1 21 0 1 0 0 0 1 22 1 1 0 1 0 1 23 1 0 0 0 0 1 24 1 1 0 0 0 1 25 1 0 1 1 0 0 26 0 0 0 1 0 1 27 1 1 0 0 0 0 28 0 0 0 1 0 0 29 1 0 0 0 0 0 30 0 0 0 0 0 1 31 0 0 0 0 0 0 32 0 1 0 0 0 0 33 0 1 0 0 0 1 34 1 0 1 0 0 0 35 0 0 0 0 0 0 36 0 0 0 0 0 0 37 0 1 1 1 0 1 38 1 0 0 1 0 0 39 1 0 0 0 0 1 40 0 0 1 0 0 1 41 1 0 0 1 1 1 42 1 0 0 1 0 0 43 1 1 0 0 0 1 44 0 1 1 0 0 0 45 0 0 0 0 0 1 46 1 0 0 0 0 1 47 0 0 0 0 0 0 48 1 0 0 0 0 0 49 1 0 0 0 0 1 50 0 0 0 0 0 0 51 0 0 1 1 0 0 52 0 1 1 1 1 1 53 1 0 0 0 0 0 54 0 0 0 1 1 0 55 0 0 0 0 0 0 56 1 0 1 1 0 0 57 1 0 0 1 0 1 58 1 0 0 0 0 0 59 1 0 0 0 0 0 60 1 1 1 1 1 1 61 1 1 1 0 0 0 62 0 0 0 1 0 1 63 0 0 0 0 0 0 64 1 1 1 0 0 0 65 0 0 0 0 0 0 66 0 0 0 0 0 0 67 0 0 1 1 1 1 68 0 1 0 0 0 0 69 1 0 0 0 0 0 70 0 0 0 1 0 0 71 0 0 0 0 0 0 72 1 0 0 0 0 0 73 1 0 0 1 0 0 74 0 1 0 1 0 0 75 1 0 0 0 0 0 76 1 0 1 0 0 1 77 1 0 0 0 0 0 78 1 0 0 1 0 1 79 1 0 1 1 1 0 80 0 0 1 0 0 1 81 0 0 0 0 0 0 82 1 1 0 1 0 0 83 0 0 0 0 0 0 84 0 0 0 1 1 0 85 1 0 0 0 0 1 86 0 1 0 0 0 0 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) UseLimit T40 Used 0.40730 -0.08160 0.03589 0.10570 CorrectAnalysis Useful -0.16898 0.15555 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -0.6685 -0.4073 -0.3257 0.5142 0.6743 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 0.40730 0.08244 4.941 4.2e-06 *** UseLimit -0.08160 0.12779 -0.639 0.525 T40 0.03589 0.13526 0.265 0.791 Used 0.10570 0.13749 0.769 0.444 CorrectAnalysis -0.16898 0.21291 -0.794 0.430 Useful 0.15555 0.12059 1.290 0.201 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 0.5071 on 80 degrees of freedom Multiple R-squared: 0.03833, Adjusted R-squared: -0.02178 F-statistic: 0.6377 on 5 and 80 DF, p-value: 0.6716 > if (n > n25) { + kp3 <- k + 3 + nmkm3 <- n - k - 3 + gqarr <- array(NA, dim=c(nmkm3-kp3+1,3)) + numgqtests <- 0 + numsignificant1 <- 0 + numsignificant5 <- 0 + numsignificant10 <- 0 + for (mypoint in kp3:nmkm3) { + j <- 0 + numgqtests <- numgqtests + 1 + for (myalt in c('greater', 'two.sided', 'less')) { + j <- j + 1 + gqarr[mypoint-kp3+1,j] <- gqtest(mylm, point=mypoint, alternative=myalt)$p.value + } + if (gqarr[mypoint-kp3+1,2] < 0.01) numsignificant1 <- numsignificant1 + 1 + if (gqarr[mypoint-kp3+1,2] < 0.05) numsignificant5 <- numsignificant5 + 1 + if (gqarr[mypoint-kp3+1,2] < 0.10) numsignificant10 <- numsignificant10 + 1 + } + gqarr + } [,1] [,2] [,3] [1,] 0.6343092 0.7313815 0.3656908 [2,] 0.6855594 0.6288812 0.3144406 [3,] 0.6527098 0.6945805 0.3472902 [4,] 0.5392698 0.9214603 0.4607302 [5,] 0.4335347 0.8670694 0.5664653 [6,] 0.3691832 0.7383663 0.6308168 [7,] 0.4586751 0.9173502 0.5413249 [8,] 0.3913483 0.7826967 0.6086517 [9,] 0.3074734 0.6149469 0.6925266 [10,] 0.2557090 0.5114179 0.7442910 [11,] 0.3927380 0.7854760 0.6072620 [12,] 0.4052760 0.8105521 0.5947240 [13,] 0.4663825 0.9327649 0.5336175 [14,] 0.4589131 0.9178262 0.5410869 [15,] 0.4115946 0.8231892 0.5884054 [16,] 0.3806304 0.7612608 0.6193696 [17,] 0.3962830 0.7925660 0.6037170 [18,] 0.4700828 0.9401656 0.5299172 [19,] 0.5551998 0.8896003 0.4448002 [20,] 0.5223146 0.9553707 0.4776854 [21,] 0.5633025 0.8733950 0.4366975 [22,] 0.5825866 0.8348268 0.4174134 [23,] 0.5441669 0.9116661 0.4558331 [24,] 0.4965855 0.9931710 0.5034145 [25,] 0.4901268 0.9802536 0.5098732 [26,] 0.4984118 0.9968237 0.5015882 [27,] 0.4627736 0.9255473 0.5372264 [28,] 0.4285165 0.8570330 0.5714835 [29,] 0.4675322 0.9350643 0.5324678 [30,] 0.4761600 0.9523200 0.5238400 [31,] 0.4559821 0.9119641 0.5440179 [32,] 0.4791355 0.9582711 0.5208645 [33,] 0.4679689 0.9359379 0.5320311 [34,] 0.4575852 0.9151704 0.5424148 [35,] 0.4563153 0.9126305 0.5436847 [36,] 0.4393099 0.8786197 0.5606901 [37,] 0.4467756 0.8935512 0.5532244 [38,] 0.4307101 0.8614202 0.5692899 [39,] 0.4088398 0.8176795 0.5911602 [40,] 0.4241849 0.8483697 0.5758151 [41,] 0.4081041 0.8162082 0.5918959 [42,] 0.3865653 0.7731305 0.6134347 [43,] 0.4461552 0.8923105 0.5538448 [44,] 0.4326414 0.8652828 0.5673586 [45,] 0.4507791 0.9015583 0.5492209 [46,] 0.4024383 0.8048765 0.5975617 [47,] 0.3792165 0.7584329 0.6207835 [48,] 0.3414845 0.6829690 0.6585155 [49,] 0.3168306 0.6336611 0.6831694 [50,] 0.3304709 0.6609418 0.6695291 [51,] 0.3502730 0.7005460 0.6497270 [52,] 0.3636518 0.7273035 0.6363482 [53,] 0.3451791 0.6903583 0.6548209 [54,] 0.3608949 0.7217898 0.6391051 [55,] 0.3332631 0.6665262 0.6667369 [56,] 0.3481344 0.6962689 0.6518656 [57,] 0.3255811 0.6511622 0.6744189 [58,] 0.3126864 0.6253729 0.6873136 [59,] 0.3023010 0.6046020 0.6976990 [60,] 0.2421571 0.4843143 0.7578429 [61,] 0.2379821 0.4759642 0.7620179 [62,] 0.3067064 0.6134127 0.6932936 [63,] 0.2998442 0.5996885 0.7001558 [64,] 0.2777383 0.5554765 0.7222617 [65,] 0.2035987 0.4071974 0.7964013 [66,] 0.2173947 0.4347893 0.7826053 [67,] 0.2181508 0.4363016 0.7818492 [68,] 0.1677194 0.3354388 0.8322806 [69,] 0.2615927 0.5231855 0.7384073 > postscript(file="/var/fisher/rcomp/tmp/1lh7o1355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index') > points(x[,1]-mysum$resid) > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/24jj41355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index') > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/3etqy1355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals') > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/4wz4z1355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals') > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/5fy4m1355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > qqnorm(mysum$resid, main='Residual Normal Q-Q Plot') > qqline(mysum$resid) > grid() > dev.off() null device 1 > (myerror <- as.ts(mysum$resid)) Time Series: Start = 1 End = 86 Frequency = 1 1 2 3 4 5 6 7 0.6384151 -0.4072979 -0.4072979 -0.4072979 -0.4072979 0.5187501 -0.4072979 8 9 10 11 12 13 14 -0.4431845 0.5927021 -0.3256984 -0.3615849 -0.4072979 -0.6685516 -0.3615849 15 16 17 18 19 20 21 0.3314484 0.2955618 -0.4538574 -0.3615849 0.5927021 0.4645430 -0.4812499 22 23 24 25 26 27 28 0.4130479 0.4371505 0.5187501 0.4511134 -0.6685516 0.6743016 -0.5130001 29 30 31 32 33 34 35 0.5927021 -0.5628495 -0.4072979 -0.3256984 -0.4812499 0.5568155 -0.4072979 36 37 38 39 40 41 42 -0.4072979 -0.6228386 0.4869999 0.4371505 -0.5987361 0.5004295 0.4869999 43 44 45 46 47 48 49 0.5187501 -0.3615849 -0.5628495 0.4371505 -0.4072979 0.5927021 0.4371505 50 51 52 53 54 55 56 -0.4072979 -0.5488866 -0.4538574 0.5927021 -0.3440189 -0.4072979 0.4511134 57 58 59 60 61 62 63 0.3314484 0.5927021 0.5927021 0.5461426 0.6384151 -0.6685516 -0.4072979 64 65 66 67 68 69 70 0.6384151 -0.4072979 -0.4072979 -0.5354570 -0.3256984 0.5927021 -0.5130001 71 72 73 74 75 76 77 -0.4072979 0.5927021 0.4869999 -0.4314005 0.5927021 0.4012639 0.5927021 78 79 80 81 82 83 84 0.3314484 0.6200946 -0.5987361 -0.4072979 0.5685995 -0.4072979 -0.3440189 85 86 0.4371505 -0.3256984 > postscript(file="/var/fisher/rcomp/tmp/6q1931355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > dum <- cbind(lag(myerror,k=1),myerror) > dum Time Series: Start = 0 End = 86 Frequency = 1 lag(myerror, k = 1) myerror 0 0.6384151 NA 1 -0.4072979 0.6384151 2 -0.4072979 -0.4072979 3 -0.4072979 -0.4072979 4 -0.4072979 -0.4072979 5 0.5187501 -0.4072979 6 -0.4072979 0.5187501 7 -0.4431845 -0.4072979 8 0.5927021 -0.4431845 9 -0.3256984 0.5927021 10 -0.3615849 -0.3256984 11 -0.4072979 -0.3615849 12 -0.6685516 -0.4072979 13 -0.3615849 -0.6685516 14 0.3314484 -0.3615849 15 0.2955618 0.3314484 16 -0.4538574 0.2955618 17 -0.3615849 -0.4538574 18 0.5927021 -0.3615849 19 0.4645430 0.5927021 20 -0.4812499 0.4645430 21 0.4130479 -0.4812499 22 0.4371505 0.4130479 23 0.5187501 0.4371505 24 0.4511134 0.5187501 25 -0.6685516 0.4511134 26 0.6743016 -0.6685516 27 -0.5130001 0.6743016 28 0.5927021 -0.5130001 29 -0.5628495 0.5927021 30 -0.4072979 -0.5628495 31 -0.3256984 -0.4072979 32 -0.4812499 -0.3256984 33 0.5568155 -0.4812499 34 -0.4072979 0.5568155 35 -0.4072979 -0.4072979 36 -0.6228386 -0.4072979 37 0.4869999 -0.6228386 38 0.4371505 0.4869999 39 -0.5987361 0.4371505 40 0.5004295 -0.5987361 41 0.4869999 0.5004295 42 0.5187501 0.4869999 43 -0.3615849 0.5187501 44 -0.5628495 -0.3615849 45 0.4371505 -0.5628495 46 -0.4072979 0.4371505 47 0.5927021 -0.4072979 48 0.4371505 0.5927021 49 -0.4072979 0.4371505 50 -0.5488866 -0.4072979 51 -0.4538574 -0.5488866 52 0.5927021 -0.4538574 53 -0.3440189 0.5927021 54 -0.4072979 -0.3440189 55 0.4511134 -0.4072979 56 0.3314484 0.4511134 57 0.5927021 0.3314484 58 0.5927021 0.5927021 59 0.5461426 0.5927021 60 0.6384151 0.5461426 61 -0.6685516 0.6384151 62 -0.4072979 -0.6685516 63 0.6384151 -0.4072979 64 -0.4072979 0.6384151 65 -0.4072979 -0.4072979 66 -0.5354570 -0.4072979 67 -0.3256984 -0.5354570 68 0.5927021 -0.3256984 69 -0.5130001 0.5927021 70 -0.4072979 -0.5130001 71 0.5927021 -0.4072979 72 0.4869999 0.5927021 73 -0.4314005 0.4869999 74 0.5927021 -0.4314005 75 0.4012639 0.5927021 76 0.5927021 0.4012639 77 0.3314484 0.5927021 78 0.6200946 0.3314484 79 -0.5987361 0.6200946 80 -0.4072979 -0.5987361 81 0.5685995 -0.4072979 82 -0.4072979 0.5685995 83 -0.3440189 -0.4072979 84 0.4371505 -0.3440189 85 -0.3256984 0.4371505 86 NA -0.3256984 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] -0.4072979 0.6384151 [2,] -0.4072979 -0.4072979 [3,] -0.4072979 -0.4072979 [4,] -0.4072979 -0.4072979 [5,] 0.5187501 -0.4072979 [6,] -0.4072979 0.5187501 [7,] -0.4431845 -0.4072979 [8,] 0.5927021 -0.4431845 [9,] -0.3256984 0.5927021 [10,] -0.3615849 -0.3256984 [11,] -0.4072979 -0.3615849 [12,] -0.6685516 -0.4072979 [13,] -0.3615849 -0.6685516 [14,] 0.3314484 -0.3615849 [15,] 0.2955618 0.3314484 [16,] -0.4538574 0.2955618 [17,] -0.3615849 -0.4538574 [18,] 0.5927021 -0.3615849 [19,] 0.4645430 0.5927021 [20,] -0.4812499 0.4645430 [21,] 0.4130479 -0.4812499 [22,] 0.4371505 0.4130479 [23,] 0.5187501 0.4371505 [24,] 0.4511134 0.5187501 [25,] -0.6685516 0.4511134 [26,] 0.6743016 -0.6685516 [27,] -0.5130001 0.6743016 [28,] 0.5927021 -0.5130001 [29,] -0.5628495 0.5927021 [30,] -0.4072979 -0.5628495 [31,] -0.3256984 -0.4072979 [32,] -0.4812499 -0.3256984 [33,] 0.5568155 -0.4812499 [34,] -0.4072979 0.5568155 [35,] -0.4072979 -0.4072979 [36,] -0.6228386 -0.4072979 [37,] 0.4869999 -0.6228386 [38,] 0.4371505 0.4869999 [39,] -0.5987361 0.4371505 [40,] 0.5004295 -0.5987361 [41,] 0.4869999 0.5004295 [42,] 0.5187501 0.4869999 [43,] -0.3615849 0.5187501 [44,] -0.5628495 -0.3615849 [45,] 0.4371505 -0.5628495 [46,] -0.4072979 0.4371505 [47,] 0.5927021 -0.4072979 [48,] 0.4371505 0.5927021 [49,] -0.4072979 0.4371505 [50,] -0.5488866 -0.4072979 [51,] -0.4538574 -0.5488866 [52,] 0.5927021 -0.4538574 [53,] -0.3440189 0.5927021 [54,] -0.4072979 -0.3440189 [55,] 0.4511134 -0.4072979 [56,] 0.3314484 0.4511134 [57,] 0.5927021 0.3314484 [58,] 0.5927021 0.5927021 [59,] 0.5461426 0.5927021 [60,] 0.6384151 0.5461426 [61,] -0.6685516 0.6384151 [62,] -0.4072979 -0.6685516 [63,] 0.6384151 -0.4072979 [64,] -0.4072979 0.6384151 [65,] -0.4072979 -0.4072979 [66,] -0.5354570 -0.4072979 [67,] -0.3256984 -0.5354570 [68,] 0.5927021 -0.3256984 [69,] -0.5130001 0.5927021 [70,] -0.4072979 -0.5130001 [71,] 0.5927021 -0.4072979 [72,] 0.4869999 0.5927021 [73,] -0.4314005 0.4869999 [74,] 0.5927021 -0.4314005 [75,] 0.4012639 0.5927021 [76,] 0.5927021 0.4012639 [77,] 0.3314484 0.5927021 [78,] 0.6200946 0.3314484 [79,] -0.5987361 0.6200946 [80,] -0.4072979 -0.5987361 [81,] 0.5685995 -0.4072979 [82,] -0.4072979 0.5685995 [83,] -0.3440189 -0.4072979 [84,] 0.4371505 -0.3440189 [85,] -0.3256984 0.4371505 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 -0.4072979 0.6384151 2 -0.4072979 -0.4072979 3 -0.4072979 -0.4072979 4 -0.4072979 -0.4072979 5 0.5187501 -0.4072979 6 -0.4072979 0.5187501 7 -0.4431845 -0.4072979 8 0.5927021 -0.4431845 9 -0.3256984 0.5927021 10 -0.3615849 -0.3256984 11 -0.4072979 -0.3615849 12 -0.6685516 -0.4072979 13 -0.3615849 -0.6685516 14 0.3314484 -0.3615849 15 0.2955618 0.3314484 16 -0.4538574 0.2955618 17 -0.3615849 -0.4538574 18 0.5927021 -0.3615849 19 0.4645430 0.5927021 20 -0.4812499 0.4645430 21 0.4130479 -0.4812499 22 0.4371505 0.4130479 23 0.5187501 0.4371505 24 0.4511134 0.5187501 25 -0.6685516 0.4511134 26 0.6743016 -0.6685516 27 -0.5130001 0.6743016 28 0.5927021 -0.5130001 29 -0.5628495 0.5927021 30 -0.4072979 -0.5628495 31 -0.3256984 -0.4072979 32 -0.4812499 -0.3256984 33 0.5568155 -0.4812499 34 -0.4072979 0.5568155 35 -0.4072979 -0.4072979 36 -0.6228386 -0.4072979 37 0.4869999 -0.6228386 38 0.4371505 0.4869999 39 -0.5987361 0.4371505 40 0.5004295 -0.5987361 41 0.4869999 0.5004295 42 0.5187501 0.4869999 43 -0.3615849 0.5187501 44 -0.5628495 -0.3615849 45 0.4371505 -0.5628495 46 -0.4072979 0.4371505 47 0.5927021 -0.4072979 48 0.4371505 0.5927021 49 -0.4072979 0.4371505 50 -0.5488866 -0.4072979 51 -0.4538574 -0.5488866 52 0.5927021 -0.4538574 53 -0.3440189 0.5927021 54 -0.4072979 -0.3440189 55 0.4511134 -0.4072979 56 0.3314484 0.4511134 57 0.5927021 0.3314484 58 0.5927021 0.5927021 59 0.5461426 0.5927021 60 0.6384151 0.5461426 61 -0.6685516 0.6384151 62 -0.4072979 -0.6685516 63 0.6384151 -0.4072979 64 -0.4072979 0.6384151 65 -0.4072979 -0.4072979 66 -0.5354570 -0.4072979 67 -0.3256984 -0.5354570 68 0.5927021 -0.3256984 69 -0.5130001 0.5927021 70 -0.4072979 -0.5130001 71 0.5927021 -0.4072979 72 0.4869999 0.5927021 73 -0.4314005 0.4869999 74 0.5927021 -0.4314005 75 0.4012639 0.5927021 76 0.5927021 0.4012639 77 0.3314484 0.5927021 78 0.6200946 0.3314484 79 -0.5987361 0.6200946 80 -0.4072979 -0.5987361 81 0.5685995 -0.4072979 82 -0.4072979 0.5685995 83 -0.3440189 -0.4072979 84 0.4371505 -0.3440189 85 -0.3256984 0.4371505 > plot(z,main=paste('Residual Lag plot, lowess, and regression line'), ylab='values of Residuals', xlab='lagged values of Residuals') > lines(lowess(z)) > abline(lm(z)) > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/7lwtz1355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function') > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/8bcpb1355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function') > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/9giwe1355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0)) > plot(mylm, las = 1, sub='Residual Diagnostics') > par(opar) > dev.off() null device 1 > if (n > n25) { + postscript(file="/var/fisher/rcomp/tmp/10sgt41355495564.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) + plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint') + grid() + dev.off() + } null device 1 > > #Note: the /var/fisher/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/fisher/rcomp/createtable") > > a<-table.start() > a<-table.row.start(a) > a<-table.element(a, 'Multiple Linear Regression - Estimated Regression Equation', 1, TRUE) > a<-table.row.end(a) > myeq <- colnames(x)[1] > myeq <- paste(myeq, '[t] = ', sep='') > for (i in 1:k){ + if (mysum$coefficients[i,1] > 0) myeq <- paste(myeq, '+', '') + myeq <- paste(myeq, mysum$coefficients[i,1], sep=' ') + if (rownames(mysum$coefficients)[i] != '(Intercept)') { + myeq <- paste(myeq, rownames(mysum$coefficients)[i], sep='') + if (rownames(mysum$coefficients)[i] != 't') myeq <- paste(myeq, '[t]', sep='') + } + } > myeq <- paste(myeq, ' + e[t]') > a<-table.row.start(a) > a<-table.element(a, myeq) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/fisher/rcomp/tmp/11zags1355495565.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,hyperlink('http://www.xycoon.com/ols1.htm','Multiple Linear Regression - Ordinary Least Squares',''), 6, TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Variable',header=TRUE) > a<-table.element(a,'Parameter',header=TRUE) > a<-table.element(a,'S.D.',header=TRUE) > a<-table.element(a,'T-STAT
H0: parameter = 0',header=TRUE) > a<-table.element(a,'2-tail p-value',header=TRUE) > a<-table.element(a,'1-tail p-value',header=TRUE) > a<-table.row.end(a) > for (i in 1:k){ + a<-table.row.start(a) + a<-table.element(a,rownames(mysum$coefficients)[i],header=TRUE) + a<-table.element(a,mysum$coefficients[i,1]) + a<-table.element(a, round(mysum$coefficients[i,2],6)) + a<-table.element(a, round(mysum$coefficients[i,3],4)) + a<-table.element(a, round(mysum$coefficients[i,4],6)) + a<-table.element(a, round(mysum$coefficients[i,4]/2,6)) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/fisher/rcomp/tmp/12wuys1355495565.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a, 'Multiple Linear Regression - Regression Statistics', 2, TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'Multiple R',1,TRUE) > a<-table.element(a, sqrt(mysum$r.squared)) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'R-squared',1,TRUE) > a<-table.element(a, mysum$r.squared) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'Adjusted R-squared',1,TRUE) > a<-table.element(a, mysum$adj.r.squared) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'F-TEST (value)',1,TRUE) > a<-table.element(a, mysum$fstatistic[1]) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE) > a<-table.element(a, mysum$fstatistic[2]) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE) > a<-table.element(a, mysum$fstatistic[3]) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'p-value',1,TRUE) > a<-table.element(a, 1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3])) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'Multiple Linear Regression - Residual Statistics', 2, TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'Residual Standard Deviation',1,TRUE) > a<-table.element(a, mysum$sigma) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'Sum Squared Residuals',1,TRUE) > a<-table.element(a, sum(myerror*myerror)) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/fisher/rcomp/tmp/13g0lc1355495565.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a, 'Multiple Linear Regression - Actuals, Interpolation, and Residuals', 4, TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a, 'Time or Index', 1, TRUE) > a<-table.element(a, 'Actuals', 1, TRUE) > a<-table.element(a, 'Interpolation
Forecast', 1, TRUE) > a<-table.element(a, 'Residuals
Prediction Error', 1, TRUE) > a<-table.row.end(a) > for (i in 1:n) { + a<-table.row.start(a) + a<-table.element(a,i, 1, TRUE) + a<-table.element(a,x[i]) + a<-table.element(a,x[i]-mysum$resid[i]) + a<-table.element(a,mysum$resid[i]) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/fisher/rcomp/tmp/14x5cx1355495565.tab") > if (n > n25) { + a<-table.start() + a<-table.row.start(a) + a<-table.element(a,'Goldfeld-Quandt test for Heteroskedasticity',4,TRUE) + a<-table.row.end(a) + a<-table.row.start(a) + a<-table.element(a,'p-values',header=TRUE) + a<-table.element(a,'Alternative Hypothesis',3,header=TRUE) + a<-table.row.end(a) + a<-table.row.start(a) + a<-table.element(a,'breakpoint index',header=TRUE) + a<-table.element(a,'greater',header=TRUE) + a<-table.element(a,'2-sided',header=TRUE) + a<-table.element(a,'less',header=TRUE) + a<-table.row.end(a) + for (mypoint in kp3:nmkm3) { + a<-table.row.start(a) + a<-table.element(a,mypoint,header=TRUE) + a<-table.element(a,gqarr[mypoint-kp3+1,1]) + a<-table.element(a,gqarr[mypoint-kp3+1,2]) + a<-table.element(a,gqarr[mypoint-kp3+1,3]) + a<-table.row.end(a) + } + a<-table.end(a) + table.save(a,file="/var/fisher/rcomp/tmp/15bxhv1355495565.tab") + a<-table.start() + a<-table.row.start(a) + a<-table.element(a,'Meta Analysis of Goldfeld-Quandt test for Heteroskedasticity',4,TRUE) + a<-table.row.end(a) + a<-table.row.start(a) + a<-table.element(a,'Description',header=TRUE) + a<-table.element(a,'# significant tests',header=TRUE) + a<-table.element(a,'% significant tests',header=TRUE) + a<-table.element(a,'OK/NOK',header=TRUE) + a<-table.row.end(a) + a<-table.row.start(a) + a<-table.element(a,'1% type I error level',header=TRUE) + a<-table.element(a,numsignificant1) + a<-table.element(a,numsignificant1/numgqtests) + if (numsignificant1/numgqtests < 0.01) dum <- 'OK' else dum <- 'NOK' + a<-table.element(a,dum) + a<-table.row.end(a) + a<-table.row.start(a) + a<-table.element(a,'5% type I error level',header=TRUE) + a<-table.element(a,numsignificant5) + a<-table.element(a,numsignificant5/numgqtests) + if (numsignificant5/numgqtests < 0.05) dum <- 'OK' else dum <- 'NOK' + a<-table.element(a,dum) + a<-table.row.end(a) + a<-table.row.start(a) + a<-table.element(a,'10% type I error level',header=TRUE) + a<-table.element(a,numsignificant10) + a<-table.element(a,numsignificant10/numgqtests) + if (numsignificant10/numgqtests < 0.1) dum <- 'OK' else dum <- 'NOK' + a<-table.element(a,dum) + a<-table.row.end(a) + a<-table.end(a) + table.save(a,file="/var/fisher/rcomp/tmp/16g92t1355495565.tab") + } > > try(system("convert tmp/1lh7o1355495564.ps tmp/1lh7o1355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/24jj41355495564.ps tmp/24jj41355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/3etqy1355495564.ps tmp/3etqy1355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/4wz4z1355495564.ps tmp/4wz4z1355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/5fy4m1355495564.ps tmp/5fy4m1355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/6q1931355495564.ps tmp/6q1931355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/7lwtz1355495564.ps tmp/7lwtz1355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/8bcpb1355495564.ps tmp/8bcpb1355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/9giwe1355495564.ps tmp/9giwe1355495564.png",intern=TRUE)) character(0) > try(system("convert tmp/10sgt41355495564.ps tmp/10sgt41355495564.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 6.311 1.584 7.903