R version 2.13.0 (2011-04-13) Copyright (C) 2011 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i486-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(210907 + ,56 + ,79 + ,94 + ,0 + ,2 + ,179321 + ,89 + ,108 + ,103 + ,0 + ,4 + ,149061 + ,44 + ,43 + ,93 + ,0 + ,0 + ,237213 + ,84 + ,78 + ,123 + ,0 + ,0 + ,173326 + ,88 + ,86 + ,148 + ,0 + ,-4 + ,133131 + ,55 + ,44 + ,90 + ,0 + ,4 + ,258873 + ,60 + ,104 + ,124 + ,0 + ,4 + ,324799 + ,154 + ,158 + ,168 + ,0 + ,0 + ,230964 + ,53 + ,102 + ,115 + ,0 + ,-1 + ,236785 + ,119 + ,77 + ,71 + ,0 + ,0 + ,344297 + ,75 + ,80 + ,108 + ,0 + ,1 + ,174724 + ,92 + ,123 + ,120 + ,0 + ,0 + ,174415 + ,100 + ,73 + ,114 + ,0 + ,3 + ,223632 + ,73 + ,105 + ,120 + ,0 + ,-1 + ,294424 + ,77 + ,107 + ,124 + ,0 + ,4 + ,325107 + ,99 + ,84 + ,126 + ,0 + ,3 + ,106408 + ,30 + ,33 + ,37 + ,0 + ,1 + ,96560 + ,76 + ,42 + ,38 + ,1 + ,0 + ,265769 + ,146 + ,96 + ,120 + ,0 + ,-2 + ,269651 + ,67 + ,106 + ,93 + ,0 + ,-3 + ,149112 + ,56 + ,56 + ,95 + ,0 + ,-4 + ,152871 + ,58 + ,59 + ,90 + ,0 + ,2 + ,362301 + ,119 + ,76 + ,110 + ,0 + ,2 + ,183167 + ,66 + ,91 + ,138 + ,0 + ,-4 + ,277965 + ,89 + ,115 + ,133 + ,0 + ,3 + ,218946 + ,41 + ,76 + ,96 + ,0 + ,2 + ,244052 + ,68 + ,101 + ,164 + ,0 + ,2 + ,341570 + ,168 + ,94 + ,78 + ,1 + ,0 + ,233328 + ,132 + ,92 + ,102 + ,0 + ,5 + ,206161 + ,71 + ,75 + ,99 + ,0 + ,-2 + ,311473 + ,112 + ,128 + ,129 + ,0 + ,0 + ,207176 + ,70 + ,56 + ,114 + ,0 + ,-2 + ,196553 + ,57 + ,41 + ,99 + ,0 + ,-3 + ,143246 + ,103 + ,67 + ,104 + ,0 + ,2 + ,182192 + ,52 + ,77 + ,138 + ,0 + ,2 + ,194979 + ,62 + ,66 + ,151 + ,0 + ,2 + ,167488 + ,45 + ,69 + ,72 + ,0 + ,0 + ,143756 + ,46 + ,105 + ,120 + ,0 + ,4 + ,275541 + ,63 + ,116 + ,115 + ,0 + ,4 + ,152299 + ,53 + ,62 + ,98 + ,0 + ,2 + ,193339 + ,78 + ,100 + ,71 + ,0 + ,2 + ,130585 + ,46 + ,67 + ,107 + ,0 + ,-4 + ,112611 + ,41 + ,46 + ,73 + ,1 + ,3 + ,148446 + ,91 + ,135 + ,129 + ,0 + ,3 + ,182079 + ,63 + ,124 + ,118 + ,0 + ,2 + ,243060 + ,63 + ,58 + ,104 + ,0 + ,-1 + ,162765 + ,32 + ,68 + ,107 + ,0 + ,-3 + ,85574 + ,34 + ,37 + ,36 + ,1 + ,0 + ,225060 + ,93 + ,93 + ,139 + ,0 + ,1 + ,133328 + ,55 + ,56 + ,56 + ,1 + ,-3 + ,100750 + ,72 + ,83 + ,93 + ,0 + ,3 + ,101523 + ,42 + ,59 + ,87 + ,1 + ,0 + ,243511 + ,71 + ,133 + ,110 + ,0 + ,0 + ,152474 + ,65 + ,106 + ,83 + ,0 + ,0 + ,132487 + ,41 + ,71 + ,98 + ,0 + ,3 + ,317394 + ,86 + ,116 + ,82 + ,0 + ,-3 + ,244749 + ,95 + ,98 + ,115 + ,0 + ,0 + ,184510 + ,49 + ,64 + ,140 + ,0 + ,-4 + ,128423 + ,64 + ,32 + ,120 + ,0 + ,2 + ,97839 + ,38 + ,25 + ,66 + ,0 + ,-1 + ,172494 + ,52 + ,46 + ,139 + ,0 + ,3 + ,229242 + ,247 + ,63 + ,119 + ,0 + ,2 + ,351619 + ,139 + ,95 + ,141 + ,0 + ,5 + ,324598 + ,110 + ,113 + ,133 + ,0 + ,2 + ,195838 + ,67 + ,111 + ,98 + ,0 + ,-2 + ,254488 + ,83 + ,120 + ,117 + ,0 + ,0 + ,199476 + ,70 + ,87 + ,105 + ,0 + ,3 + ,92499 + ,32 + ,25 + ,55 + ,1 + ,-2 + ,224330 + ,83 + ,131 + ,132 + ,0 + ,0 + ,181633 + ,70 + ,47 + ,73 + ,0 + ,6 + ,271856 + ,103 + ,109 + ,86 + ,0 + ,-3 + ,95227 + ,34 + ,37 + ,48 + ,0 + ,3 + ,98146 + ,40 + ,15 + ,48 + ,1 + ,0 + ,118612 + ,46 + ,54 + ,43 + ,1 + ,-2 + ,65475 + ,18 + ,16 + ,46 + ,1 + ,1 + ,108446 + ,60 + ,22 + ,65 + ,1 + ,0 + ,121848 + ,39 + ,37 + ,52 + ,1 + ,2 + ,76302 + ,31 + ,29 + ,68 + ,1 + ,2 + ,98104 + ,54 + ,55 + ,47 + ,1 + ,-3 + ,30989 + ,14 + ,5 + ,41 + ,1 + ,-2 + ,31774 + ,23 + ,0 + ,47 + ,1 + ,1 + ,150580 + ,77 + ,27 + ,71 + ,1 + ,-4 + ,54157 + ,19 + ,37 + ,30 + ,1 + ,0 + ,59382 + ,49 + ,29 + ,24 + ,1 + ,1 + ,84105 + ,20 + ,17 + ,63 + ,1 + ,0) + ,dim=c(6 + ,85) + ,dimnames=list(c('time' + ,'logins' + ,'BC' + ,'LFM' + ,'Course' + ,'Totaal') + ,1:85)) > y <- array(NA,dim=c(6,85),dimnames=list(c('time','logins','BC','LFM','Course','Totaal'),1:85)) > 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' > #'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 > 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 time logins BC LFM Course Totaal 1 210907 56 79 94 0 2 2 179321 89 108 103 0 4 3 149061 44 43 93 0 0 4 237213 84 78 123 0 0 5 173326 88 86 148 0 -4 6 133131 55 44 90 0 4 7 258873 60 104 124 0 4 8 324799 154 158 168 0 0 9 230964 53 102 115 0 -1 10 236785 119 77 71 0 0 11 344297 75 80 108 0 1 12 174724 92 123 120 0 0 13 174415 100 73 114 0 3 14 223632 73 105 120 0 -1 15 294424 77 107 124 0 4 16 325107 99 84 126 0 3 17 106408 30 33 37 0 1 18 96560 76 42 38 1 0 19 265769 146 96 120 0 -2 20 269651 67 106 93 0 -3 21 149112 56 56 95 0 -4 22 152871 58 59 90 0 2 23 362301 119 76 110 0 2 24 183167 66 91 138 0 -4 25 277965 89 115 133 0 3 26 218946 41 76 96 0 2 27 244052 68 101 164 0 2 28 341570 168 94 78 1 0 29 233328 132 92 102 0 5 30 206161 71 75 99 0 -2 31 311473 112 128 129 0 0 32 207176 70 56 114 0 -2 33 196553 57 41 99 0 -3 34 143246 103 67 104 0 2 35 182192 52 77 138 0 2 36 194979 62 66 151 0 2 37 167488 45 69 72 0 0 38 143756 46 105 120 0 4 39 275541 63 116 115 0 4 40 152299 53 62 98 0 2 41 193339 78 100 71 0 2 42 130585 46 67 107 0 -4 43 112611 41 46 73 1 3 44 148446 91 135 129 0 3 45 182079 63 124 118 0 2 46 243060 63 58 104 0 -1 47 162765 32 68 107 0 -3 48 85574 34 37 36 1 0 49 225060 93 93 139 0 1 50 133328 55 56 56 1 -3 51 100750 72 83 93 0 3 52 101523 42 59 87 1 0 53 243511 71 133 110 0 0 54 152474 65 106 83 0 0 55 132487 41 71 98 0 3 56 317394 86 116 82 0 -3 57 244749 95 98 115 0 0 58 184510 49 64 140 0 -4 59 128423 64 32 120 0 2 60 97839 38 25 66 0 -1 61 172494 52 46 139 0 3 62 229242 247 63 119 0 2 63 351619 139 95 141 0 5 64 324598 110 113 133 0 2 65 195838 67 111 98 0 -2 66 254488 83 120 117 0 0 67 199476 70 87 105 0 3 68 92499 32 25 55 1 -2 69 224330 83 131 132 0 0 70 181633 70 47 73 0 6 71 271856 103 109 86 0 -3 72 95227 34 37 48 0 3 73 98146 40 15 48 1 0 74 118612 46 54 43 1 -2 75 65475 18 16 46 1 1 76 108446 60 22 65 1 0 77 121848 39 37 52 1 2 78 76302 31 29 68 1 2 79 98104 54 55 47 1 -3 80 30989 14 5 41 1 -2 81 31774 23 0 47 1 1 82 150580 77 27 71 1 -4 83 54157 19 37 30 1 0 84 59382 49 29 24 1 1 85 84105 20 17 63 1 0 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) logins BC LFM Course Totaal 38756.5 818.9 796.8 376.5 -22382.5 -743.9 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -118740 -24618 -781 22061 140456 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 38756.5 25594.4 1.514 0.133952 logins 818.9 167.2 4.898 5.04e-06 *** BC 796.8 218.9 3.639 0.000486 *** LFM 376.5 247.6 1.521 0.132282 Course -22382.5 17971.1 -1.245 0.216641 Totaal -743.9 2127.6 -0.350 0.727544 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 46990 on 79 degrees of freedom Multiple R-squared: 0.6744, Adjusted R-squared: 0.6538 F-statistic: 32.72 on 5 and 79 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,] 0.5809897 8.380206e-01 4.190103e-01 [2,] 0.5474391 9.051219e-01 4.525609e-01 [3,] 0.9664984 6.700318e-02 3.350159e-02 [4,] 0.9832079 3.358413e-02 1.679207e-02 [5,] 0.9774350 4.512994e-02 2.256497e-02 [6,] 0.9602558 7.948837e-02 3.974419e-02 [7,] 0.9621088 7.578232e-02 3.789116e-02 [8,] 0.9839842 3.203150e-02 1.601575e-02 [9,] 0.9737979 5.240427e-02 2.620213e-02 [10,] 0.9596246 8.075078e-02 4.037539e-02 [11,] 0.9407992 1.184015e-01 5.920077e-02 [12,] 0.9501539 9.969218e-02 4.984609e-02 [13,] 0.9294479 1.411043e-01 7.055214e-02 [14,] 0.9079565 1.840870e-01 9.204348e-02 [15,] 0.9866349 2.673012e-02 1.336506e-02 [16,] 0.9824675 3.506508e-02 1.753254e-02 [17,] 0.9750933 4.981349e-02 2.490675e-02 [18,] 0.9743029 5.139417e-02 2.569708e-02 [19,] 0.9621302 7.573954e-02 3.786977e-02 [20,] 0.9830795 3.384093e-02 1.692047e-02 [21,] 0.9823293 3.534146e-02 1.767073e-02 [22,] 0.9738489 5.230220e-02 2.615110e-02 [23,] 0.9680197 6.396056e-02 3.198028e-02 [24,] 0.9568871 8.622582e-02 4.311291e-02 [25,] 0.9506484 9.870327e-02 4.935164e-02 [26,] 0.9699974 6.000512e-02 3.000256e-02 [27,] 0.9585112 8.297757e-02 4.148879e-02 [28,] 0.9417870 1.164260e-01 5.821300e-02 [29,] 0.9213727 1.572545e-01 7.862726e-02 [30,] 0.9348158 1.303684e-01 6.518422e-02 [31,] 0.9433607 1.132787e-01 5.663933e-02 [32,] 0.9245993 1.508014e-01 7.540071e-02 [33,] 0.9025384 1.949233e-01 9.746163e-02 [34,] 0.8997093 2.005814e-01 1.002907e-01 [35,] 0.8711727 2.576545e-01 1.288273e-01 [36,] 0.9727927 5.441466e-02 2.720733e-02 [37,] 0.9743598 5.128036e-02 2.564018e-02 [38,] 0.9853059 2.938825e-02 1.469412e-02 [39,] 0.9776306 4.473883e-02 2.236942e-02 [40,] 0.9669322 6.613554e-02 3.306777e-02 [41,] 0.9543744 9.125119e-02 4.562560e-02 [42,] 0.9350605 1.298791e-01 6.493954e-02 [43,] 0.9828755 3.424892e-02 1.712446e-02 [44,] 0.9818220 3.635608e-02 1.817804e-02 [45,] 0.9744297 5.114061e-02 2.557030e-02 [46,] 0.9836760 3.264801e-02 1.632400e-02 [47,] 0.9834170 3.316604e-02 1.658302e-02 [48,] 0.9969486 6.102783e-03 3.051391e-03 [49,] 0.9946058 1.078841e-02 5.394203e-03 [50,] 0.9906646 1.867085e-02 9.335425e-03 [51,] 0.9884110 2.317801e-02 1.158901e-02 [52,] 0.9814670 3.706596e-02 1.853298e-02 [53,] 0.9740997 5.180057e-02 2.590028e-02 [54,] 0.9999546 9.088765e-05 4.544383e-05 [55,] 0.9999125 1.750453e-04 8.752264e-05 [56,] 0.9999301 1.397846e-04 6.989228e-05 [57,] 0.9998305 3.389376e-04 1.694688e-04 [58,] 0.9996307 7.386935e-04 3.693467e-04 [59,] 0.9990340 1.932029e-03 9.660147e-04 [60,] 0.9979025 4.194921e-03 2.097460e-03 [61,] 0.9995610 8.779812e-04 4.389906e-04 [62,] 0.9987794 2.441150e-03 1.220575e-03 [63,] 0.9967243 6.551347e-03 3.275673e-03 [64,] 0.9909869 1.802626e-02 9.013131e-03 [65,] 0.9869072 2.618567e-02 1.309283e-02 [66,] 0.9703748 5.925033e-02 2.962517e-02 [67,] 0.9355264 1.289471e-01 6.447357e-02 [68,] 0.8484896 3.030207e-01 1.515104e-01 > postscript(file="/var/wessaorg/rcomp/tmp/1xrw11324634828.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/wessaorg/rcomp/tmp/2ldlt1324634828.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/wessaorg/rcomp/tmp/3oln01324634828.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/wessaorg/rcomp/tmp/45jmo1324634828.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/wessaorg/rcomp/tmp/50vs71324634828.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 = 85 Frequency = 1 1 2 3 4 5 6 29437.8096 -54180.7221 4991.8789 21203.0836 -64722.8459 -16638.0479 7 8 9 10 11 12 44400.2157 -29222.9032 23486.7512 12489.0932 140455.8728 -82562.8789 13 14 15 16 17 18 -45093.4926 -4496.9645 63638.9896 93134.4154 3601.7271 -29826.1876 19 20 21 22 23 24 -15715.0424 54317.6076 -18871.0697 -12794.3112 125604.6811 -37083.8146 25 26 27 28 29 30 26845.8207 51398.0922 8868.8296 83347.7921 -21519.0081 10736.7052 31 32 33 34 35 36 30434.6964 22061.4003 38940.4265 -70917.1241 -10975.5312 -2508.2807 37 38 39 40 41 42 9790.9475 -58542.3325 52438.8879 -14674.2766 -14218.7717 -42491.7513 43 44 45 46 47 48 753.2365 -118740.4986 -50014.7239 66593.6344 1100.4317 -1679.9650 49 50 51 52 53 54 -15552.7419 3975.2061 -95888.8446 -29014.9994 -780.6819 -55224.5105 55 56 57 58 59 60 -31086.2022 82674.9624 6807.5029 -1058.9439 -31938.9764 -17551.8495 61 62 63 64 65 66 4393.9905 -105308.6286 73964.1581 57130.8762 -24618.0983 8091.4831 67 68 69 70 71 72 -3230.5424 7802.2955 -36479.1551 25078.4884 27286.3870 -16696.2800 73 74 75 76 77 78 18989.1403 3862.0408 5034.8737 931.8578 25962.5271 -12682.3532 79 80 81 82 83 84 -26244.2442 -17759.7631 -20388.8566 19925.3276 -18553.7243 -28519.4155 85 14085.2115 > postscript(file="/var/wessaorg/rcomp/tmp/6uapt1324634828.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 = 85 Frequency = 1 lag(myerror, k = 1) myerror 0 29437.8096 NA 1 -54180.7221 29437.8096 2 4991.8789 -54180.7221 3 21203.0836 4991.8789 4 -64722.8459 21203.0836 5 -16638.0479 -64722.8459 6 44400.2157 -16638.0479 7 -29222.9032 44400.2157 8 23486.7512 -29222.9032 9 12489.0932 23486.7512 10 140455.8728 12489.0932 11 -82562.8789 140455.8728 12 -45093.4926 -82562.8789 13 -4496.9645 -45093.4926 14 63638.9896 -4496.9645 15 93134.4154 63638.9896 16 3601.7271 93134.4154 17 -29826.1876 3601.7271 18 -15715.0424 -29826.1876 19 54317.6076 -15715.0424 20 -18871.0697 54317.6076 21 -12794.3112 -18871.0697 22 125604.6811 -12794.3112 23 -37083.8146 125604.6811 24 26845.8207 -37083.8146 25 51398.0922 26845.8207 26 8868.8296 51398.0922 27 83347.7921 8868.8296 28 -21519.0081 83347.7921 29 10736.7052 -21519.0081 30 30434.6964 10736.7052 31 22061.4003 30434.6964 32 38940.4265 22061.4003 33 -70917.1241 38940.4265 34 -10975.5312 -70917.1241 35 -2508.2807 -10975.5312 36 9790.9475 -2508.2807 37 -58542.3325 9790.9475 38 52438.8879 -58542.3325 39 -14674.2766 52438.8879 40 -14218.7717 -14674.2766 41 -42491.7513 -14218.7717 42 753.2365 -42491.7513 43 -118740.4986 753.2365 44 -50014.7239 -118740.4986 45 66593.6344 -50014.7239 46 1100.4317 66593.6344 47 -1679.9650 1100.4317 48 -15552.7419 -1679.9650 49 3975.2061 -15552.7419 50 -95888.8446 3975.2061 51 -29014.9994 -95888.8446 52 -780.6819 -29014.9994 53 -55224.5105 -780.6819 54 -31086.2022 -55224.5105 55 82674.9624 -31086.2022 56 6807.5029 82674.9624 57 -1058.9439 6807.5029 58 -31938.9764 -1058.9439 59 -17551.8495 -31938.9764 60 4393.9905 -17551.8495 61 -105308.6286 4393.9905 62 73964.1581 -105308.6286 63 57130.8762 73964.1581 64 -24618.0983 57130.8762 65 8091.4831 -24618.0983 66 -3230.5424 8091.4831 67 7802.2955 -3230.5424 68 -36479.1551 7802.2955 69 25078.4884 -36479.1551 70 27286.3870 25078.4884 71 -16696.2800 27286.3870 72 18989.1403 -16696.2800 73 3862.0408 18989.1403 74 5034.8737 3862.0408 75 931.8578 5034.8737 76 25962.5271 931.8578 77 -12682.3532 25962.5271 78 -26244.2442 -12682.3532 79 -17759.7631 -26244.2442 80 -20388.8566 -17759.7631 81 19925.3276 -20388.8566 82 -18553.7243 19925.3276 83 -28519.4155 -18553.7243 84 14085.2115 -28519.4155 85 NA 14085.2115 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] -54180.7221 29437.8096 [2,] 4991.8789 -54180.7221 [3,] 21203.0836 4991.8789 [4,] -64722.8459 21203.0836 [5,] -16638.0479 -64722.8459 [6,] 44400.2157 -16638.0479 [7,] -29222.9032 44400.2157 [8,] 23486.7512 -29222.9032 [9,] 12489.0932 23486.7512 [10,] 140455.8728 12489.0932 [11,] -82562.8789 140455.8728 [12,] -45093.4926 -82562.8789 [13,] -4496.9645 -45093.4926 [14,] 63638.9896 -4496.9645 [15,] 93134.4154 63638.9896 [16,] 3601.7271 93134.4154 [17,] -29826.1876 3601.7271 [18,] -15715.0424 -29826.1876 [19,] 54317.6076 -15715.0424 [20,] -18871.0697 54317.6076 [21,] -12794.3112 -18871.0697 [22,] 125604.6811 -12794.3112 [23,] -37083.8146 125604.6811 [24,] 26845.8207 -37083.8146 [25,] 51398.0922 26845.8207 [26,] 8868.8296 51398.0922 [27,] 83347.7921 8868.8296 [28,] -21519.0081 83347.7921 [29,] 10736.7052 -21519.0081 [30,] 30434.6964 10736.7052 [31,] 22061.4003 30434.6964 [32,] 38940.4265 22061.4003 [33,] -70917.1241 38940.4265 [34,] -10975.5312 -70917.1241 [35,] -2508.2807 -10975.5312 [36,] 9790.9475 -2508.2807 [37,] -58542.3325 9790.9475 [38,] 52438.8879 -58542.3325 [39,] -14674.2766 52438.8879 [40,] -14218.7717 -14674.2766 [41,] -42491.7513 -14218.7717 [42,] 753.2365 -42491.7513 [43,] -118740.4986 753.2365 [44,] -50014.7239 -118740.4986 [45,] 66593.6344 -50014.7239 [46,] 1100.4317 66593.6344 [47,] -1679.9650 1100.4317 [48,] -15552.7419 -1679.9650 [49,] 3975.2061 -15552.7419 [50,] -95888.8446 3975.2061 [51,] -29014.9994 -95888.8446 [52,] -780.6819 -29014.9994 [53,] -55224.5105 -780.6819 [54,] -31086.2022 -55224.5105 [55,] 82674.9624 -31086.2022 [56,] 6807.5029 82674.9624 [57,] -1058.9439 6807.5029 [58,] -31938.9764 -1058.9439 [59,] -17551.8495 -31938.9764 [60,] 4393.9905 -17551.8495 [61,] -105308.6286 4393.9905 [62,] 73964.1581 -105308.6286 [63,] 57130.8762 73964.1581 [64,] -24618.0983 57130.8762 [65,] 8091.4831 -24618.0983 [66,] -3230.5424 8091.4831 [67,] 7802.2955 -3230.5424 [68,] -36479.1551 7802.2955 [69,] 25078.4884 -36479.1551 [70,] 27286.3870 25078.4884 [71,] -16696.2800 27286.3870 [72,] 18989.1403 -16696.2800 [73,] 3862.0408 18989.1403 [74,] 5034.8737 3862.0408 [75,] 931.8578 5034.8737 [76,] 25962.5271 931.8578 [77,] -12682.3532 25962.5271 [78,] -26244.2442 -12682.3532 [79,] -17759.7631 -26244.2442 [80,] -20388.8566 -17759.7631 [81,] 19925.3276 -20388.8566 [82,] -18553.7243 19925.3276 [83,] -28519.4155 -18553.7243 [84,] 14085.2115 -28519.4155 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 -54180.7221 29437.8096 2 4991.8789 -54180.7221 3 21203.0836 4991.8789 4 -64722.8459 21203.0836 5 -16638.0479 -64722.8459 6 44400.2157 -16638.0479 7 -29222.9032 44400.2157 8 23486.7512 -29222.9032 9 12489.0932 23486.7512 10 140455.8728 12489.0932 11 -82562.8789 140455.8728 12 -45093.4926 -82562.8789 13 -4496.9645 -45093.4926 14 63638.9896 -4496.9645 15 93134.4154 63638.9896 16 3601.7271 93134.4154 17 -29826.1876 3601.7271 18 -15715.0424 -29826.1876 19 54317.6076 -15715.0424 20 -18871.0697 54317.6076 21 -12794.3112 -18871.0697 22 125604.6811 -12794.3112 23 -37083.8146 125604.6811 24 26845.8207 -37083.8146 25 51398.0922 26845.8207 26 8868.8296 51398.0922 27 83347.7921 8868.8296 28 -21519.0081 83347.7921 29 10736.7052 -21519.0081 30 30434.6964 10736.7052 31 22061.4003 30434.6964 32 38940.4265 22061.4003 33 -70917.1241 38940.4265 34 -10975.5312 -70917.1241 35 -2508.2807 -10975.5312 36 9790.9475 -2508.2807 37 -58542.3325 9790.9475 38 52438.8879 -58542.3325 39 -14674.2766 52438.8879 40 -14218.7717 -14674.2766 41 -42491.7513 -14218.7717 42 753.2365 -42491.7513 43 -118740.4986 753.2365 44 -50014.7239 -118740.4986 45 66593.6344 -50014.7239 46 1100.4317 66593.6344 47 -1679.9650 1100.4317 48 -15552.7419 -1679.9650 49 3975.2061 -15552.7419 50 -95888.8446 3975.2061 51 -29014.9994 -95888.8446 52 -780.6819 -29014.9994 53 -55224.5105 -780.6819 54 -31086.2022 -55224.5105 55 82674.9624 -31086.2022 56 6807.5029 82674.9624 57 -1058.9439 6807.5029 58 -31938.9764 -1058.9439 59 -17551.8495 -31938.9764 60 4393.9905 -17551.8495 61 -105308.6286 4393.9905 62 73964.1581 -105308.6286 63 57130.8762 73964.1581 64 -24618.0983 57130.8762 65 8091.4831 -24618.0983 66 -3230.5424 8091.4831 67 7802.2955 -3230.5424 68 -36479.1551 7802.2955 69 25078.4884 -36479.1551 70 27286.3870 25078.4884 71 -16696.2800 27286.3870 72 18989.1403 -16696.2800 73 3862.0408 18989.1403 74 5034.8737 3862.0408 75 931.8578 5034.8737 76 25962.5271 931.8578 77 -12682.3532 25962.5271 78 -26244.2442 -12682.3532 79 -17759.7631 -26244.2442 80 -20388.8566 -17759.7631 81 19925.3276 -20388.8566 82 -18553.7243 19925.3276 83 -28519.4155 -18553.7243 84 14085.2115 -28519.4155 > 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/wessaorg/rcomp/tmp/7dp0w1324634828.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/wessaorg/rcomp/tmp/8ckoa1324634828.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/wessaorg/rcomp/tmp/9clbx1324634828.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/wessaorg/rcomp/tmp/10ywqn1324634828.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/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/wessaorg/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/wessaorg/rcomp/tmp/11henn1324634828.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/wessaorg/rcomp/tmp/12x9bx1324634828.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/wessaorg/rcomp/tmp/13545z1324634828.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/wessaorg/rcomp/tmp/14fuss1324634828.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/wessaorg/rcomp/tmp/15jtpw1324634828.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/wessaorg/rcomp/tmp/16ee6q1324634828.tab") + } > > try(system("convert tmp/1xrw11324634828.ps tmp/1xrw11324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/2ldlt1324634828.ps tmp/2ldlt1324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/3oln01324634828.ps tmp/3oln01324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/45jmo1324634828.ps tmp/45jmo1324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/50vs71324634828.ps tmp/50vs71324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/6uapt1324634828.ps tmp/6uapt1324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/7dp0w1324634828.ps tmp/7dp0w1324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/8ckoa1324634828.ps tmp/8ckoa1324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/9clbx1324634828.ps tmp/9clbx1324634828.png",intern=TRUE)) character(0) > try(system("convert tmp/10ywqn1324634828.ps tmp/10ywqn1324634828.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 3.620 0.529 4.206