R version 2.9.0 (2009-04-17) Copyright (C) 2009 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. 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(9487,1169,8700,2154,9627,2249,8947,2687,9283,4359,8829,5382,9947,4459,9628,6398,9318,4596,9605,3024,8640,1887,9214,2070,9567,1351,8547,2218,9185,2461,9470,3028,9123,4784,9278,4975,10170,4607,9434,6249,9655,4809,9429,3157,8739,1910,9552,2228,9687,1594,9019,2467,9672,2222,9206,3607,9069,4685,9788,4962,10312,5770,10105,5480,9863,5000,9656,3228,9295,1993,9946,2288,9701,1580,9049,2111,10190,2192,9706,3601,9765,4665,9893,4876,9994,5813,10433,5589,10073,5331,10112,3075,9266,2002,9820,2306,10097,1507,9115,1992,10411,2487,9678,3490,10408,4647,10153,5594,10368,5611,10581,5788,10597,6204,10680,3013,9738,1931,9556,2549),dim=c(2,60),dimnames=list(c('Y','X'),1:60)) > y <- array(NA,dim=c(2,60),dimnames=list(c('Y','X'),1:60)) > 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 = 'Include Monthly 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 Attaching package: 'zoo' The following object(s) are masked from package:base : 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 Y X M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 1 9487 1169 1 0 0 0 0 0 0 0 0 0 0 2 8700 2154 0 1 0 0 0 0 0 0 0 0 0 3 9627 2249 0 0 1 0 0 0 0 0 0 0 0 4 8947 2687 0 0 0 1 0 0 0 0 0 0 0 5 9283 4359 0 0 0 0 1 0 0 0 0 0 0 6 8829 5382 0 0 0 0 0 1 0 0 0 0 0 7 9947 4459 0 0 0 0 0 0 1 0 0 0 0 8 9628 6398 0 0 0 0 0 0 0 1 0 0 0 9 9318 4596 0 0 0 0 0 0 0 0 1 0 0 10 9605 3024 0 0 0 0 0 0 0 0 0 1 0 11 8640 1887 0 0 0 0 0 0 0 0 0 0 1 12 9214 2070 0 0 0 0 0 0 0 0 0 0 0 13 9567 1351 1 0 0 0 0 0 0 0 0 0 0 14 8547 2218 0 1 0 0 0 0 0 0 0 0 0 15 9185 2461 0 0 1 0 0 0 0 0 0 0 0 16 9470 3028 0 0 0 1 0 0 0 0 0 0 0 17 9123 4784 0 0 0 0 1 0 0 0 0 0 0 18 9278 4975 0 0 0 0 0 1 0 0 0 0 0 19 10170 4607 0 0 0 0 0 0 1 0 0 0 0 20 9434 6249 0 0 0 0 0 0 0 1 0 0 0 21 9655 4809 0 0 0 0 0 0 0 0 1 0 0 22 9429 3157 0 0 0 0 0 0 0 0 0 1 0 23 8739 1910 0 0 0 0 0 0 0 0 0 0 1 24 9552 2228 0 0 0 0 0 0 0 0 0 0 0 25 9687 1594 1 0 0 0 0 0 0 0 0 0 0 26 9019 2467 0 1 0 0 0 0 0 0 0 0 0 27 9672 2222 0 0 1 0 0 0 0 0 0 0 0 28 9206 3607 0 0 0 1 0 0 0 0 0 0 0 29 9069 4685 0 0 0 0 1 0 0 0 0 0 0 30 9788 4962 0 0 0 0 0 1 0 0 0 0 0 31 10312 5770 0 0 0 0 0 0 1 0 0 0 0 32 10105 5480 0 0 0 0 0 0 0 1 0 0 0 33 9863 5000 0 0 0 0 0 0 0 0 1 0 0 34 9656 3228 0 0 0 0 0 0 0 0 0 1 0 35 9295 1993 0 0 0 0 0 0 0 0 0 0 1 36 9946 2288 0 0 0 0 0 0 0 0 0 0 0 37 9701 1580 1 0 0 0 0 0 0 0 0 0 0 38 9049 2111 0 1 0 0 0 0 0 0 0 0 0 39 10190 2192 0 0 1 0 0 0 0 0 0 0 0 40 9706 3601 0 0 0 1 0 0 0 0 0 0 0 41 9765 4665 0 0 0 0 1 0 0 0 0 0 0 42 9893 4876 0 0 0 0 0 1 0 0 0 0 0 43 9994 5813 0 0 0 0 0 0 1 0 0 0 0 44 10433 5589 0 0 0 0 0 0 0 1 0 0 0 45 10073 5331 0 0 0 0 0 0 0 0 1 0 0 46 10112 3075 0 0 0 0 0 0 0 0 0 1 0 47 9266 2002 0 0 0 0 0 0 0 0 0 0 1 48 9820 2306 0 0 0 0 0 0 0 0 0 0 0 49 10097 1507 1 0 0 0 0 0 0 0 0 0 0 50 9115 1992 0 1 0 0 0 0 0 0 0 0 0 51 10411 2487 0 0 1 0 0 0 0 0 0 0 0 52 9678 3490 0 0 0 1 0 0 0 0 0 0 0 53 10408 4647 0 0 0 0 1 0 0 0 0 0 0 54 10153 5594 0 0 0 0 0 1 0 0 0 0 0 55 10368 5611 0 0 0 0 0 0 1 0 0 0 0 56 10581 5788 0 0 0 0 0 0 0 1 0 0 0 57 10597 6204 0 0 0 0 0 0 0 0 1 0 0 58 10680 3013 0 0 0 0 0 0 0 0 0 1 0 59 9738 1931 0 0 0 0 0 0 0 0 0 0 1 60 9556 2549 0 0 0 0 0 0 0 0 0 0 0 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) X M1 M2 M3 M4 9128.8690 0.2136 271.3222 -710.2840 192.1380 -428.5914 M5 M6 M7 M8 M9 M10 -587.7521 -642.3107 -92.4307 -353.0063 -335.7611 105.5378 M11 -408.6113 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -807.09 -272.19 -17.35 233.37 874.34 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 9128.8690 442.1685 20.646 < 2e-16 *** X 0.2136 0.1752 1.219 0.22900 M1 271.3222 302.5254 0.897 0.37437 M2 -710.2840 264.0903 -2.690 0.00987 ** M3 192.1380 263.5779 0.729 0.46964 M4 -428.5914 315.9193 -1.357 0.18137 M5 -587.7521 487.4070 -1.206 0.23390 M6 -642.3107 567.7357 -1.131 0.26365 M7 -92.4307 582.4081 -0.159 0.87458 M8 -353.0063 685.7349 -0.515 0.60912 M9 -335.7611 572.4287 -0.587 0.56031 M10 105.5378 299.4103 0.352 0.72605 M11 -408.6113 270.3026 -1.512 0.13731 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 416.6 on 47 degrees of freedom Multiple R-squared: 0.4859, Adjusted R-squared: 0.3546 F-statistic: 3.702 on 12 and 47 DF, p-value: 0.0005965 > 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.29786695 0.59573391 0.7021330 [2,] 0.20640599 0.41281198 0.7935940 [3,] 0.19024265 0.38048531 0.8097573 [4,] 0.12684516 0.25369032 0.8731548 [5,] 0.13726200 0.27452400 0.8627380 [6,] 0.11736186 0.23472373 0.8826381 [7,] 0.10356250 0.20712500 0.8964375 [8,] 0.08803031 0.17606062 0.9119697 [9,] 0.07528184 0.15056368 0.9247182 [10,] 0.04919124 0.09838249 0.9508088 [11,] 0.04748247 0.09496494 0.9525175 [12,] 0.05382550 0.10765099 0.9461745 [13,] 0.04387507 0.08775015 0.9561249 [14,] 0.12275334 0.24550668 0.8772467 [15,] 0.28478797 0.56957595 0.7152120 [16,] 0.22475443 0.44950886 0.7752456 [17,] 0.30301736 0.60603472 0.6969826 [18,] 0.29748232 0.59496465 0.7025177 [19,] 0.49449724 0.98899448 0.5055028 [20,] 0.54013608 0.91972785 0.4598639 [21,] 0.57737385 0.84525230 0.4226262 [22,] 0.55894501 0.88210998 0.4410550 [23,] 0.47715318 0.95430636 0.5228468 [24,] 0.48554013 0.97108026 0.5144599 [25,] 0.41451204 0.82902409 0.5854880 [26,] 0.57320855 0.85358291 0.4267915 [27,] 0.48667082 0.97334164 0.5133292 [28,] 0.45879656 0.91759313 0.5412034 [29,] 0.34689129 0.69378258 0.6531087 > postscript(file="/var/www/html/rcomp/tmp/1k66j1258906967.ps",horizontal=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/www/html/rcomp/tmp/2bvbc1258906967.ps",horizontal=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/www/html/rcomp/tmp/3d02p1258906967.ps",horizontal=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/www/html/rcomp/tmp/4do8t1258906967.ps",horizontal=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/www/html/rcomp/tmp/5pi6q1258906967.ps",horizontal=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 = 60 Frequency = 1 1 2 3 4 5 6 -162.875065 -178.652589 -174.365394 -327.187275 -189.144958 -807.086321 7 8 9 10 11 12 -41.825099 -514.395713 -456.756190 -275.295501 -483.297359 -356.995204 13 14 15 16 17 18 -121.747994 -345.322190 -661.645947 122.979382 -439.919653 -271.156202 19 20 21 22 23 24 149.563948 -676.571173 -165.250331 -479.702641 -389.209872 -52.742031 25 26 27 28 29 30 -53.649760 73.494518 -123.598531 -264.687790 -472.774488 241.620436 31 32 33 34 35 36 43.161666 158.677627 1.954453 -267.867355 149.062364 328.442718 37 38 39 40 41 42 -36.659535 179.531674 400.809095 236.593735 227.497262 364.988962 43 44 45 46 47 48 -284.022597 463.396588 141.256985 220.811535 118.140077 198.598142 49 50 51 52 53 54 374.932354 270.948588 558.800777 232.301949 874.341837 471.633125 55 56 57 58 59 60 133.122082 568.892672 478.795083 802.053961 605.304790 -117.303624 > postscript(file="/var/www/html/rcomp/tmp/6tnij1258906967.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > dum <- cbind(lag(myerror,k=1),myerror) > dum Time Series: Start = 0 End = 60 Frequency = 1 lag(myerror, k = 1) myerror 0 -162.875065 NA 1 -178.652589 -162.875065 2 -174.365394 -178.652589 3 -327.187275 -174.365394 4 -189.144958 -327.187275 5 -807.086321 -189.144958 6 -41.825099 -807.086321 7 -514.395713 -41.825099 8 -456.756190 -514.395713 9 -275.295501 -456.756190 10 -483.297359 -275.295501 11 -356.995204 -483.297359 12 -121.747994 -356.995204 13 -345.322190 -121.747994 14 -661.645947 -345.322190 15 122.979382 -661.645947 16 -439.919653 122.979382 17 -271.156202 -439.919653 18 149.563948 -271.156202 19 -676.571173 149.563948 20 -165.250331 -676.571173 21 -479.702641 -165.250331 22 -389.209872 -479.702641 23 -52.742031 -389.209872 24 -53.649760 -52.742031 25 73.494518 -53.649760 26 -123.598531 73.494518 27 -264.687790 -123.598531 28 -472.774488 -264.687790 29 241.620436 -472.774488 30 43.161666 241.620436 31 158.677627 43.161666 32 1.954453 158.677627 33 -267.867355 1.954453 34 149.062364 -267.867355 35 328.442718 149.062364 36 -36.659535 328.442718 37 179.531674 -36.659535 38 400.809095 179.531674 39 236.593735 400.809095 40 227.497262 236.593735 41 364.988962 227.497262 42 -284.022597 364.988962 43 463.396588 -284.022597 44 141.256985 463.396588 45 220.811535 141.256985 46 118.140077 220.811535 47 198.598142 118.140077 48 374.932354 198.598142 49 270.948588 374.932354 50 558.800777 270.948588 51 232.301949 558.800777 52 874.341837 232.301949 53 471.633125 874.341837 54 133.122082 471.633125 55 568.892672 133.122082 56 478.795083 568.892672 57 802.053961 478.795083 58 605.304790 802.053961 59 -117.303624 605.304790 60 NA -117.303624 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] -178.652589 -162.875065 [2,] -174.365394 -178.652589 [3,] -327.187275 -174.365394 [4,] -189.144958 -327.187275 [5,] -807.086321 -189.144958 [6,] -41.825099 -807.086321 [7,] -514.395713 -41.825099 [8,] -456.756190 -514.395713 [9,] -275.295501 -456.756190 [10,] -483.297359 -275.295501 [11,] -356.995204 -483.297359 [12,] -121.747994 -356.995204 [13,] -345.322190 -121.747994 [14,] -661.645947 -345.322190 [15,] 122.979382 -661.645947 [16,] -439.919653 122.979382 [17,] -271.156202 -439.919653 [18,] 149.563948 -271.156202 [19,] -676.571173 149.563948 [20,] -165.250331 -676.571173 [21,] -479.702641 -165.250331 [22,] -389.209872 -479.702641 [23,] -52.742031 -389.209872 [24,] -53.649760 -52.742031 [25,] 73.494518 -53.649760 [26,] -123.598531 73.494518 [27,] -264.687790 -123.598531 [28,] -472.774488 -264.687790 [29,] 241.620436 -472.774488 [30,] 43.161666 241.620436 [31,] 158.677627 43.161666 [32,] 1.954453 158.677627 [33,] -267.867355 1.954453 [34,] 149.062364 -267.867355 [35,] 328.442718 149.062364 [36,] -36.659535 328.442718 [37,] 179.531674 -36.659535 [38,] 400.809095 179.531674 [39,] 236.593735 400.809095 [40,] 227.497262 236.593735 [41,] 364.988962 227.497262 [42,] -284.022597 364.988962 [43,] 463.396588 -284.022597 [44,] 141.256985 463.396588 [45,] 220.811535 141.256985 [46,] 118.140077 220.811535 [47,] 198.598142 118.140077 [48,] 374.932354 198.598142 [49,] 270.948588 374.932354 [50,] 558.800777 270.948588 [51,] 232.301949 558.800777 [52,] 874.341837 232.301949 [53,] 471.633125 874.341837 [54,] 133.122082 471.633125 [55,] 568.892672 133.122082 [56,] 478.795083 568.892672 [57,] 802.053961 478.795083 [58,] 605.304790 802.053961 [59,] -117.303624 605.304790 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 -178.652589 -162.875065 2 -174.365394 -178.652589 3 -327.187275 -174.365394 4 -189.144958 -327.187275 5 -807.086321 -189.144958 6 -41.825099 -807.086321 7 -514.395713 -41.825099 8 -456.756190 -514.395713 9 -275.295501 -456.756190 10 -483.297359 -275.295501 11 -356.995204 -483.297359 12 -121.747994 -356.995204 13 -345.322190 -121.747994 14 -661.645947 -345.322190 15 122.979382 -661.645947 16 -439.919653 122.979382 17 -271.156202 -439.919653 18 149.563948 -271.156202 19 -676.571173 149.563948 20 -165.250331 -676.571173 21 -479.702641 -165.250331 22 -389.209872 -479.702641 23 -52.742031 -389.209872 24 -53.649760 -52.742031 25 73.494518 -53.649760 26 -123.598531 73.494518 27 -264.687790 -123.598531 28 -472.774488 -264.687790 29 241.620436 -472.774488 30 43.161666 241.620436 31 158.677627 43.161666 32 1.954453 158.677627 33 -267.867355 1.954453 34 149.062364 -267.867355 35 328.442718 149.062364 36 -36.659535 328.442718 37 179.531674 -36.659535 38 400.809095 179.531674 39 236.593735 400.809095 40 227.497262 236.593735 41 364.988962 227.497262 42 -284.022597 364.988962 43 463.396588 -284.022597 44 141.256985 463.396588 45 220.811535 141.256985 46 118.140077 220.811535 47 198.598142 118.140077 48 374.932354 198.598142 49 270.948588 374.932354 50 558.800777 270.948588 51 232.301949 558.800777 52 874.341837 232.301949 53 471.633125 874.341837 54 133.122082 471.633125 55 568.892672 133.122082 56 478.795083 568.892672 57 802.053961 478.795083 58 605.304790 802.053961 59 -117.303624 605.304790 > 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/www/html/rcomp/tmp/7j4cw1258906967.ps",horizontal=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/www/html/rcomp/tmp/8pd4x1258906967.ps",horizontal=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/www/html/rcomp/tmp/9e5ex1258906967.ps",horizontal=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/www/html/rcomp/tmp/10dv8q1258906967.ps",horizontal=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/www/html/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/html/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/www/html/rcomp/tmp/118gnk1258906967.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/www/html/rcomp/tmp/12geel1258906967.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/www/html/rcomp/tmp/13t1oo1258906967.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/www/html/rcomp/tmp/141cf61258906967.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/www/html/rcomp/tmp/15g7sp1258906967.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/www/html/rcomp/tmp/16ku7m1258906967.tab") + } > > system("convert tmp/1k66j1258906967.ps tmp/1k66j1258906967.png") > system("convert tmp/2bvbc1258906967.ps tmp/2bvbc1258906967.png") > system("convert tmp/3d02p1258906967.ps tmp/3d02p1258906967.png") > system("convert tmp/4do8t1258906967.ps tmp/4do8t1258906967.png") > system("convert tmp/5pi6q1258906967.ps tmp/5pi6q1258906967.png") > system("convert tmp/6tnij1258906967.ps tmp/6tnij1258906967.png") > system("convert tmp/7j4cw1258906967.ps tmp/7j4cw1258906967.png") > system("convert tmp/8pd4x1258906967.ps tmp/8pd4x1258906967.png") > system("convert tmp/9e5ex1258906967.ps tmp/9e5ex1258906967.png") > system("convert tmp/10dv8q1258906967.ps tmp/10dv8q1258906967.png") > > > proc.time() user system elapsed 2.332 1.497 3.107