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 + ,21221 + ,6086 + ,7 + ,2 + ,7106 + ,2273 + ,6 + ,3 + ,4794 + ,2601 + ,6 + ,4 + ,6776 + ,2708 + ,6 + ,5 + ,7122 + ,3321 + ,6 + ,6 + ,12025 + ,2106 + ,6 + ,7 + ,18945 + ,6353 + ,5 + ,8 + ,24380 + ,7912 + ,5 + ,9 + ,6942 + ,1905 + ,5 + ,10 + ,7370 + ,3841 + ,5 + ,11 + ,13143 + ,2737 + ,5 + ,12 + ,11831 + ,4449 + ,5 + ,13 + ,3102 + ,1253 + ,5 + ,14 + ,7355 + ,2372 + ,5 + ,15 + ,7749 + ,3303 + ,5 + ,16 + ,10306 + ,2876 + ,5 + ,17 + ,7860 + ,2939 + ,4 + ,18 + ,9138 + ,2149 + ,4 + ,19 + ,6545 + ,2545 + ,4 + ,20 + ,9829 + ,3734 + ,4 + ,21 + ,9460 + ,1911 + ,4 + ,22 + ,30631 + ,875 + ,4 + ,23 + ,4807 + ,1466 + ,4 + ,24 + ,4618 + ,1562 + ,4 + ,25 + ,7304 + ,1203 + ,4 + ,26 + ,12621 + ,2860 + ,4 + ,27 + ,9320 + ,3245 + ,4 + ,28 + ,9569 + ,2668 + ,4 + ,29 + ,8480 + ,2111 + ,4 + ,30 + ,4718 + ,942 + ,4 + ,31 + ,7342 + ,2690 + ,4 + ,32 + ,8230 + ,3675 + ,4 + ,33 + ,10007 + ,3798 + ,4 + ,34 + ,9176 + ,672 + ,4 + ,35 + ,2334 + ,523 + ,4 + ,36 + ,5459 + ,2352 + ,3 + ,37 + ,5440 + ,1476 + ,3 + ,38 + ,9027 + ,2646 + ,3 + ,39 + ,9852 + ,2500 + ,3 + ,40 + ,3122 + ,720 + ,3 + ,41 + ,8863 + ,2483 + ,3 + ,42 + ,4286 + ,1550 + ,3 + ,43 + ,9119 + ,2551 + ,3 + ,44 + ,7297 + ,2544 + ,3 + ,45 + ,8009 + ,2870 + ,3 + ,46 + ,5049 + ,1288 + ,3 + ,47 + ,5057 + ,1474 + ,3 + ,48 + ,4372 + ,1521 + ,3 + ,49 + ,3393 + ,1331 + ,3 + ,50 + ,7586 + ,1685 + ,3 + ,51 + ,7440 + ,2755 + ,3 + ,52 + ,8636 + ,1268 + ,3 + ,53 + ,10846 + ,1266 + ,3 + ,54 + ,6843 + ,1564 + ,3 + ,55 + ,5102 + ,1062 + ,3 + ,56 + ,4183 + ,1306 + ,3 + ,57 + ,10083 + ,4566 + ,3 + ,58 + ,5273 + ,1655 + ,3 + ,59 + ,6433 + ,1423 + ,3 + ,60 + ,6116 + ,474 + ,3 + ,61 + ,5887 + ,725 + ,3 + ,62 + ,7415 + ,1061 + ,2 + ,63 + ,7953 + ,687 + ,2 + ,64 + ,8447 + ,2411 + ,2 + ,65 + ,4733 + ,1705 + ,2 + ,66 + ,11478 + ,2422 + ,2 + ,67 + ,4234 + ,1900 + ,2 + ,68 + ,11089 + ,2848 + ,2 + ,69 + ,6390 + ,1398 + ,2 + ,70 + ,7737 + ,2625 + ,2) + ,dim=c(4 + ,70) + ,dimnames=list(c('Ranking' + ,'Characters' + ,'Revisions' + ,'Hours ') + ,1:70)) > y <- array(NA,dim=c(4,70),dimnames=list(c('Ranking','Characters','Revisions','Hours '),1:70)) > 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 = '1' > 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 Ranking Characters Revisions Hours\r 1 1 21221 6086 7 2 2 7106 2273 6 3 3 4794 2601 6 4 4 6776 2708 6 5 5 7122 3321 6 6 6 12025 2106 6 7 7 18945 6353 5 8 8 24380 7912 5 9 9 6942 1905 5 10 10 7370 3841 5 11 11 13143 2737 5 12 12 11831 4449 5 13 13 3102 1253 5 14 14 7355 2372 5 15 15 7749 3303 5 16 16 10306 2876 5 17 17 7860 2939 4 18 18 9138 2149 4 19 19 6545 2545 4 20 20 9829 3734 4 21 21 9460 1911 4 22 22 30631 875 4 23 23 4807 1466 4 24 24 4618 1562 4 25 25 7304 1203 4 26 26 12621 2860 4 27 27 9320 3245 4 28 28 9569 2668 4 29 29 8480 2111 4 30 30 4718 942 4 31 31 7342 2690 4 32 32 8230 3675 4 33 33 10007 3798 4 34 34 9176 672 4 35 35 2334 523 4 36 36 5459 2352 3 37 37 5440 1476 3 38 38 9027 2646 3 39 39 9852 2500 3 40 40 3122 720 3 41 41 8863 2483 3 42 42 4286 1550 3 43 43 9119 2551 3 44 44 7297 2544 3 45 45 8009 2870 3 46 46 5049 1288 3 47 47 5057 1474 3 48 48 4372 1521 3 49 49 3393 1331 3 50 50 7586 1685 3 51 51 7440 2755 3 52 52 8636 1268 3 53 53 10846 1266 3 54 54 6843 1564 3 55 55 5102 1062 3 56 56 4183 1306 3 57 57 10083 4566 3 58 58 5273 1655 3 59 59 6433 1423 3 60 60 6116 474 3 61 61 5887 725 3 62 62 7415 1061 2 63 63 7953 687 2 64 64 8447 2411 2 65 65 4733 1705 2 66 66 11478 2422 2 67 67 4234 1900 2 68 68 11089 2848 2 69 69 6390 1398 2 70 70 7737 2625 2 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) Characters Revisions `Hours\\r` 9.649e+01 -1.141e-04 -3.876e-04 -1.598e+01 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -13.5173 -4.7695 -0.3675 4.2772 21.1773 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 9.649e+01 2.921e+00 33.037 <2e-16 *** Characters -1.141e-04 2.242e-04 -0.509 0.613 Revisions -3.876e-04 8.230e-04 -0.471 0.639 `Hours\\r` -1.598e+01 8.145e-01 -19.623 <2e-16 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 7.077 on 66 degrees of freedom Multiple R-squared: 0.8843, Adjusted R-squared: 0.8791 F-statistic: 168.2 on 3 and 66 DF, p-value: < 2.2e-16 > 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,] 2.275577e-02 4.551154e-02 0.977244231 [2,] 4.970907e-03 9.941814e-03 0.995029093 [3,] 1.244508e-03 2.489016e-03 0.998755492 [4,] 1.034344e-03 2.068688e-03 0.998965656 [5,] 3.527939e-04 7.055878e-04 0.999647206 [6,] 4.800418e-04 9.600836e-04 0.999519958 [7,] 4.074384e-04 8.148768e-04 0.999592562 [8,] 4.400187e-04 8.800375e-04 0.999559981 [9,] 8.518153e-04 1.703631e-03 0.999148185 [10,] 1.370307e-03 2.740613e-03 0.998629693 [11,] 6.898755e-04 1.379751e-03 0.999310124 [12,] 3.446193e-04 6.892386e-04 0.999655381 [13,] 1.883773e-04 3.767545e-04 0.999811623 [14,] 1.265984e-04 2.531968e-04 0.999873402 [15,] 7.300153e-05 1.460031e-04 0.999926998 [16,] 3.284060e-05 6.568119e-05 0.999967159 [17,] 3.730359e-05 7.460718e-05 0.999962696 [18,] 5.154466e-05 1.030893e-04 0.999948455 [19,] 6.933661e-05 1.386732e-04 0.999930663 [20,] 1.665333e-04 3.330667e-04 0.999833467 [21,] 4.697770e-04 9.395540e-04 0.999530223 [22,] 1.006258e-03 2.012516e-03 0.998993742 [23,] 1.889053e-03 3.778106e-03 0.998110947 [24,] 2.850131e-03 5.700262e-03 0.997149869 [25,] 5.856680e-03 1.171336e-02 0.994143320 [26,] 1.257658e-02 2.515316e-02 0.987423418 [27,] 2.394821e-02 4.789641e-02 0.976051795 [28,] 3.355121e-02 6.710242e-02 0.966448788 [29,] 4.218799e-02 8.437598e-02 0.957812012 [30,] 4.941092e-02 9.882184e-02 0.950589081 [31,] 6.904510e-02 1.380902e-01 0.930954896 [32,] 8.868656e-02 1.773731e-01 0.911313438 [33,] 1.273517e-01 2.547033e-01 0.872648340 [34,] 1.828943e-01 3.657886e-01 0.817105677 [35,] 2.544048e-01 5.088096e-01 0.745595202 [36,] 3.426536e-01 6.853071e-01 0.657346428 [37,] 4.564098e-01 9.128196e-01 0.543590193 [38,] 5.753945e-01 8.492110e-01 0.424605504 [39,] 7.099182e-01 5.801635e-01 0.290081761 [40,] 8.034370e-01 3.931261e-01 0.196563046 [41,] 8.787510e-01 2.424979e-01 0.121248967 [42,] 9.333189e-01 1.333623e-01 0.066681126 [43,] 9.728300e-01 5.433999e-02 0.027169996 [44,] 9.864386e-01 2.712281e-02 0.013561407 [45,] 9.951480e-01 9.704045e-03 0.004852023 [46,] 9.965314e-01 6.937279e-03 0.003468639 [47,] 9.968986e-01 6.202886e-03 0.003101443 [48,] 9.976235e-01 4.753007e-03 0.002376504 [49,] 9.976616e-01 4.676809e-03 0.002338404 [50,] 9.974023e-01 5.195433e-03 0.002597717 [51,] 9.983774e-01 3.245108e-03 0.001622554 [52,] 9.983445e-01 3.311068e-03 0.001655534 [53,] 9.978446e-01 4.310858e-03 0.002155429 [54,] 9.942542e-01 1.149153e-02 0.005745766 [55,] 9.847022e-01 3.059562e-02 0.015297810 [56,] 9.655083e-01 6.898344e-02 0.034491722 [57,] 9.112089e-01 1.775821e-01 0.088791060 > postscript(file="/var/fisher/rcomp/tmp/1j19m1352156420.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/2s0fi1352156420.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/31ser1352156420.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/432xn1352156420.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/50gjg1352156420.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 = 70 Frequency = 1 1 2 3 4 5 6 21.1772877 3.1055307 3.9688639 5.2364707 6.5135259 7.6020424 7 8 9 10 11 12 -4.9459317 -2.7216070 -6.0393147 -4.2401531 -3.0093537 -1.4955330 13 14 15 16 17 18 -2.7301344 -0.8111997 0.5945783 1.7208294 -13.5173406 -12.6777042 19 20 21 22 23 24 -11.8200772 -9.9845714 -9.7332063 -6.7192113 -8.4365593 -7.4209170 25 26 27 28 29 30 -6.2535931 -4.0047501 -3.2321659 -2.4273821 -1.7675065 -1.6497987 31 32 33 34 35 36 0.3270539 1.8101237 3.0605419 2.7541957 2.9158075 -11.0022945 37 38 39 40 41 42 -10.3439705 -8.4812575 -7.4437135 -7.9014440 -5.5631425 -5.4469567 43 44 45 46 47 48 -3.5075795 -2.7181743 -1.5105917 -1.4614443 -0.3884442 0.5516161 49 50 51 52 53 54 1.3662790 2.9818796 4.3799178 4.9400646 6.1914403 6.8502112 55 56 57 58 59 60 7.4570127 8.4467252 11.3833550 10.7063499 11.7487852 12.3448163 61 62 63 64 65 66 13.4159677 -1.2629815 -0.3465480 1.3779798 1.6806075 3.7280662 67 68 69 70 3.6992494 5.8487864 5.7506807 7.3799113 > postscript(file="/var/fisher/rcomp/tmp/6w3nc1352156420.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 = 70 Frequency = 1 lag(myerror, k = 1) myerror 0 21.1772877 NA 1 3.1055307 21.1772877 2 3.9688639 3.1055307 3 5.2364707 3.9688639 4 6.5135259 5.2364707 5 7.6020424 6.5135259 6 -4.9459317 7.6020424 7 -2.7216070 -4.9459317 8 -6.0393147 -2.7216070 9 -4.2401531 -6.0393147 10 -3.0093537 -4.2401531 11 -1.4955330 -3.0093537 12 -2.7301344 -1.4955330 13 -0.8111997 -2.7301344 14 0.5945783 -0.8111997 15 1.7208294 0.5945783 16 -13.5173406 1.7208294 17 -12.6777042 -13.5173406 18 -11.8200772 -12.6777042 19 -9.9845714 -11.8200772 20 -9.7332063 -9.9845714 21 -6.7192113 -9.7332063 22 -8.4365593 -6.7192113 23 -7.4209170 -8.4365593 24 -6.2535931 -7.4209170 25 -4.0047501 -6.2535931 26 -3.2321659 -4.0047501 27 -2.4273821 -3.2321659 28 -1.7675065 -2.4273821 29 -1.6497987 -1.7675065 30 0.3270539 -1.6497987 31 1.8101237 0.3270539 32 3.0605419 1.8101237 33 2.7541957 3.0605419 34 2.9158075 2.7541957 35 -11.0022945 2.9158075 36 -10.3439705 -11.0022945 37 -8.4812575 -10.3439705 38 -7.4437135 -8.4812575 39 -7.9014440 -7.4437135 40 -5.5631425 -7.9014440 41 -5.4469567 -5.5631425 42 -3.5075795 -5.4469567 43 -2.7181743 -3.5075795 44 -1.5105917 -2.7181743 45 -1.4614443 -1.5105917 46 -0.3884442 -1.4614443 47 0.5516161 -0.3884442 48 1.3662790 0.5516161 49 2.9818796 1.3662790 50 4.3799178 2.9818796 51 4.9400646 4.3799178 52 6.1914403 4.9400646 53 6.8502112 6.1914403 54 7.4570127 6.8502112 55 8.4467252 7.4570127 56 11.3833550 8.4467252 57 10.7063499 11.3833550 58 11.7487852 10.7063499 59 12.3448163 11.7487852 60 13.4159677 12.3448163 61 -1.2629815 13.4159677 62 -0.3465480 -1.2629815 63 1.3779798 -0.3465480 64 1.6806075 1.3779798 65 3.7280662 1.6806075 66 3.6992494 3.7280662 67 5.8487864 3.6992494 68 5.7506807 5.8487864 69 7.3799113 5.7506807 70 NA 7.3799113 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] 3.1055307 21.1772877 [2,] 3.9688639 3.1055307 [3,] 5.2364707 3.9688639 [4,] 6.5135259 5.2364707 [5,] 7.6020424 6.5135259 [6,] -4.9459317 7.6020424 [7,] -2.7216070 -4.9459317 [8,] -6.0393147 -2.7216070 [9,] -4.2401531 -6.0393147 [10,] -3.0093537 -4.2401531 [11,] -1.4955330 -3.0093537 [12,] -2.7301344 -1.4955330 [13,] -0.8111997 -2.7301344 [14,] 0.5945783 -0.8111997 [15,] 1.7208294 0.5945783 [16,] -13.5173406 1.7208294 [17,] -12.6777042 -13.5173406 [18,] -11.8200772 -12.6777042 [19,] -9.9845714 -11.8200772 [20,] -9.7332063 -9.9845714 [21,] -6.7192113 -9.7332063 [22,] -8.4365593 -6.7192113 [23,] -7.4209170 -8.4365593 [24,] -6.2535931 -7.4209170 [25,] -4.0047501 -6.2535931 [26,] -3.2321659 -4.0047501 [27,] -2.4273821 -3.2321659 [28,] -1.7675065 -2.4273821 [29,] -1.6497987 -1.7675065 [30,] 0.3270539 -1.6497987 [31,] 1.8101237 0.3270539 [32,] 3.0605419 1.8101237 [33,] 2.7541957 3.0605419 [34,] 2.9158075 2.7541957 [35,] -11.0022945 2.9158075 [36,] -10.3439705 -11.0022945 [37,] -8.4812575 -10.3439705 [38,] -7.4437135 -8.4812575 [39,] -7.9014440 -7.4437135 [40,] -5.5631425 -7.9014440 [41,] -5.4469567 -5.5631425 [42,] -3.5075795 -5.4469567 [43,] -2.7181743 -3.5075795 [44,] -1.5105917 -2.7181743 [45,] -1.4614443 -1.5105917 [46,] -0.3884442 -1.4614443 [47,] 0.5516161 -0.3884442 [48,] 1.3662790 0.5516161 [49,] 2.9818796 1.3662790 [50,] 4.3799178 2.9818796 [51,] 4.9400646 4.3799178 [52,] 6.1914403 4.9400646 [53,] 6.8502112 6.1914403 [54,] 7.4570127 6.8502112 [55,] 8.4467252 7.4570127 [56,] 11.3833550 8.4467252 [57,] 10.7063499 11.3833550 [58,] 11.7487852 10.7063499 [59,] 12.3448163 11.7487852 [60,] 13.4159677 12.3448163 [61,] -1.2629815 13.4159677 [62,] -0.3465480 -1.2629815 [63,] 1.3779798 -0.3465480 [64,] 1.6806075 1.3779798 [65,] 3.7280662 1.6806075 [66,] 3.6992494 3.7280662 [67,] 5.8487864 3.6992494 [68,] 5.7506807 5.8487864 [69,] 7.3799113 5.7506807 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 3.1055307 21.1772877 2 3.9688639 3.1055307 3 5.2364707 3.9688639 4 6.5135259 5.2364707 5 7.6020424 6.5135259 6 -4.9459317 7.6020424 7 -2.7216070 -4.9459317 8 -6.0393147 -2.7216070 9 -4.2401531 -6.0393147 10 -3.0093537 -4.2401531 11 -1.4955330 -3.0093537 12 -2.7301344 -1.4955330 13 -0.8111997 -2.7301344 14 0.5945783 -0.8111997 15 1.7208294 0.5945783 16 -13.5173406 1.7208294 17 -12.6777042 -13.5173406 18 -11.8200772 -12.6777042 19 -9.9845714 -11.8200772 20 -9.7332063 -9.9845714 21 -6.7192113 -9.7332063 22 -8.4365593 -6.7192113 23 -7.4209170 -8.4365593 24 -6.2535931 -7.4209170 25 -4.0047501 -6.2535931 26 -3.2321659 -4.0047501 27 -2.4273821 -3.2321659 28 -1.7675065 -2.4273821 29 -1.6497987 -1.7675065 30 0.3270539 -1.6497987 31 1.8101237 0.3270539 32 3.0605419 1.8101237 33 2.7541957 3.0605419 34 2.9158075 2.7541957 35 -11.0022945 2.9158075 36 -10.3439705 -11.0022945 37 -8.4812575 -10.3439705 38 -7.4437135 -8.4812575 39 -7.9014440 -7.4437135 40 -5.5631425 -7.9014440 41 -5.4469567 -5.5631425 42 -3.5075795 -5.4469567 43 -2.7181743 -3.5075795 44 -1.5105917 -2.7181743 45 -1.4614443 -1.5105917 46 -0.3884442 -1.4614443 47 0.5516161 -0.3884442 48 1.3662790 0.5516161 49 2.9818796 1.3662790 50 4.3799178 2.9818796 51 4.9400646 4.3799178 52 6.1914403 4.9400646 53 6.8502112 6.1914403 54 7.4570127 6.8502112 55 8.4467252 7.4570127 56 11.3833550 8.4467252 57 10.7063499 11.3833550 58 11.7487852 10.7063499 59 12.3448163 11.7487852 60 13.4159677 12.3448163 61 -1.2629815 13.4159677 62 -0.3465480 -1.2629815 63 1.3779798 -0.3465480 64 1.6806075 1.3779798 65 3.7280662 1.6806075 66 3.6992494 3.7280662 67 5.8487864 3.6992494 68 5.7506807 5.8487864 69 7.3799113 5.7506807 > 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/70w8g1352156420.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/8o2ls1352156420.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/9l8fr1352156420.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/104r3f1352156420.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/11wlww1352156420.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/12x2l31352156420.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/133qlv1352156421.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/14dcdx1352156421.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/150zpm1352156421.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/16l1hp1352156421.tab") + } > > try(system("convert tmp/1j19m1352156420.ps tmp/1j19m1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/2s0fi1352156420.ps tmp/2s0fi1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/31ser1352156420.ps tmp/31ser1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/432xn1352156420.ps tmp/432xn1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/50gjg1352156420.ps tmp/50gjg1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/6w3nc1352156420.ps tmp/6w3nc1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/70w8g1352156420.ps tmp/70w8g1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/8o2ls1352156420.ps tmp/8o2ls1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/9l8fr1352156420.ps tmp/9l8fr1352156420.png",intern=TRUE)) character(0) > try(system("convert tmp/104r3f1352156420.ps tmp/104r3f1352156420.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 6.405 1.201 7.606