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(100 + ,0 + ,108.1560276 + ,0 + ,114.0150276 + ,0 + ,102.1880309 + ,0 + ,110.3672031 + ,0 + ,96.8602511 + ,0 + ,94.1944583 + ,0 + ,99.51621961 + ,0 + ,94.06333487 + ,0 + ,97.5541476 + ,0 + ,78.15062422 + ,0 + ,81.2434643 + ,0 + ,92.36262465 + ,0 + ,96.06324371 + ,0 + ,114.0523777 + ,0 + ,110.6616666 + ,0 + ,104.9171949 + ,0 + ,90.00187193 + ,0 + ,95.7008067 + ,0 + ,86.02741157 + ,0 + ,84.85287668 + ,0 + ,100.04328 + ,0 + ,80.91713823 + ,0 + ,74.06539709 + ,0 + ,77.30281369 + ,0 + ,97.23043249 + ,0 + ,90.75515676 + ,0 + ,100.5614455 + ,0 + ,92.01293267 + ,0 + ,99.24012138 + ,0 + ,105.8672755 + ,0 + ,90.9920463 + ,0 + ,93.30624423 + ,0 + ,91.17419413 + ,0 + ,77.33295039 + ,0 + ,91.1277721 + ,0 + ,85.01249943 + ,0 + ,83.90390242 + ,0 + ,104.8626302 + ,0 + ,110.9039108 + ,0 + ,95.43714373 + ,0 + ,111.6238727 + ,0 + ,108.8925403 + ,0 + ,96.17511682 + ,0 + ,101.9740205 + ,0 + ,99.11953031 + ,0 + ,86.78158147 + ,0 + ,118.4195003 + ,0 + ,118.7441447 + ,0 + ,106.5296192 + ,0 + ,134.7772694 + ,0 + ,104.6778714 + ,0 + ,105.2954304 + ,0 + ,139.4139849 + ,0 + ,103.6060491 + ,0 + ,99.78182974 + ,0 + ,103.4610301 + ,0 + ,120.0594945 + ,0 + ,96.71377168 + ,0 + ,107.1308929 + ,0 + ,105.3608372 + ,0 + ,111.6942359 + ,0 + ,132.0519998 + ,0 + ,126.8037879 + ,0 + ,154.4824253 + ,0 + ,141.5570984 + ,0 + ,109.9506882 + ,0 + ,127.904198 + ,0 + ,133.0888617 + ,0 + ,120.0796299 + ,0 + ,117.5557142 + ,0 + ,143.0362309 + ,0 + ,159.982927 + ,1 + ,128.5991124 + ,1 + ,149.7373327 + ,1 + ,126.8169313 + ,1 + ,140.9639674 + ,1 + ,137.6691981 + ,1 + ,117.9402337 + ,1 + ,122.3095247 + ,1 + ,127.7804207 + ,1 + ,136.1677176 + ,1 + ,116.2405856 + ,1 + ,123.1576893 + ,1 + ,116.3400234 + ,1 + ,108.6119282 + ,1 + ,125.8982264 + ,1 + ,112.8003105 + ,1 + ,107.5182447 + ,1 + ,135.0955413 + ,1 + ,115.5096488 + ,1 + ,115.8640759 + ,1 + ,104.5883906 + ,1 + ,163.7213386 + ,1 + ,113.4482275 + ,1 + ,98.0428844 + ,1 + ,116.7868521 + ,1 + ,126.5330444 + ,1 + ,113.0336597 + ,1 + ,124.3392163 + ,1 + ,109.8298759 + ,1 + ,124.4434777 + ,1 + ,111.5039454 + ,1 + ,102.0350019 + ,1 + ,116.8726598 + ,1 + ,112.2073122 + ,1 + ,101.1513902 + ,1 + ,124.4255108 + ,1) + ,dim=c(2 + ,108) + ,dimnames=list(c('Y' + ,'X') + ,1:108)) > y <- array(NA,dim=c(2,108),dimnames=list(c('Y','X'),1:108)) > 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 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 1 100.00000 0 2 108.15603 0 3 114.01503 0 4 102.18803 0 5 110.36720 0 6 96.86025 0 7 94.19446 0 8 99.51622 0 9 94.06333 0 10 97.55415 0 11 78.15062 0 12 81.24346 0 13 92.36262 0 14 96.06324 0 15 114.05238 0 16 110.66167 0 17 104.91719 0 18 90.00187 0 19 95.70081 0 20 86.02741 0 21 84.85288 0 22 100.04328 0 23 80.91714 0 24 74.06540 0 25 77.30281 0 26 97.23043 0 27 90.75516 0 28 100.56145 0 29 92.01293 0 30 99.24012 0 31 105.86728 0 32 90.99205 0 33 93.30624 0 34 91.17419 0 35 77.33295 0 36 91.12777 0 37 85.01250 0 38 83.90390 0 39 104.86263 0 40 110.90391 0 41 95.43714 0 42 111.62387 0 43 108.89254 0 44 96.17512 0 45 101.97402 0 46 99.11953 0 47 86.78158 0 48 118.41950 0 49 118.74414 0 50 106.52962 0 51 134.77727 0 52 104.67787 0 53 105.29543 0 54 139.41398 0 55 103.60605 0 56 99.78183 0 57 103.46103 0 58 120.05949 0 59 96.71377 0 60 107.13089 0 61 105.36084 0 62 111.69424 0 63 132.05200 0 64 126.80379 0 65 154.48243 0 66 141.55710 0 67 109.95069 0 68 127.90420 0 69 133.08886 0 70 120.07963 0 71 117.55571 0 72 143.03623 0 73 159.98293 1 74 128.59911 1 75 149.73733 1 76 126.81693 1 77 140.96397 1 78 137.66920 1 79 117.94023 1 80 122.30952 1 81 127.78042 1 82 136.16772 1 83 116.24059 1 84 123.15769 1 85 116.34002 1 86 108.61193 1 87 125.89823 1 88 112.80031 1 89 107.51824 1 90 135.09554 1 91 115.50965 1 92 115.86408 1 93 104.58839 1 94 163.72134 1 95 113.44823 1 96 98.04288 1 97 116.78685 1 98 126.53304 1 99 113.03366 1 100 124.33922 1 101 109.82988 1 102 124.44348 1 103 111.50395 1 104 102.03500 1 105 116.87266 1 106 112.20731 1 107 101.15139 1 108 124.42551 1 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) X 103.89 18.00 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -29.820 -9.962 -3.583 6.727 50.597 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 103.885 1.913 54.315 < 2e-16 *** X 18.003 3.313 5.434 3.54e-07 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 16.23 on 106 degrees of freedom Multiple R-squared: 0.2179, Adjusted R-squared: 0.2105 F-statistic: 29.53 on 1 and 106 DF, p-value: 3.538e-07 > 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.078762146 0.157524292 0.92123785 [2,] 0.060867330 0.121734660 0.93913267 [3,] 0.052497037 0.104994075 0.94750296 [4,] 0.024200770 0.048401540 0.97579923 [5,] 0.017098886 0.034197771 0.98290111 [6,] 0.008066857 0.016133714 0.99193314 [7,] 0.052957189 0.105914379 0.94704281 [8,] 0.080376905 0.160753811 0.91962309 [9,] 0.053135576 0.106271153 0.94686442 [10,] 0.031333300 0.062666600 0.96866670 [11,] 0.038828183 0.077656366 0.96117182 [12,] 0.033481121 0.066962242 0.96651888 [13,] 0.021419004 0.042838009 0.97858100 [14,] 0.016760307 0.033520613 0.98323969 [15,] 0.010128745 0.020257490 0.98987126 [16,] 0.010168943 0.020337887 0.98983106 [17,] 0.010748919 0.021497838 0.98925108 [18,] 0.006384375 0.012768751 0.99361562 [19,] 0.009504669 0.019009337 0.99049533 [20,] 0.025325783 0.050651567 0.97467422 [21,] 0.040306599 0.080613199 0.95969340 [22,] 0.028186995 0.056373990 0.97181300 [23,] 0.021286525 0.042573049 0.97871348 [24,] 0.014964578 0.029929156 0.98503542 [25,] 0.010846016 0.021692032 0.98915398 [26,] 0.007348986 0.014697971 0.99265101 [27,] 0.005855997 0.011711994 0.99414400 [28,] 0.004396668 0.008793336 0.99560333 [29,] 0.003078564 0.006157127 0.99692144 [30,] 0.002319431 0.004638862 0.99768057 [31,] 0.005239836 0.010479673 0.99476016 [32,] 0.004228311 0.008456622 0.99577169 [33,] 0.004898210 0.009796420 0.99510179 [34,] 0.006413451 0.012826902 0.99358655 [35,] 0.005529711 0.011059421 0.99447029 [36,] 0.006236417 0.012472833 0.99376358 [37,] 0.005056771 0.010113542 0.99494323 [38,] 0.005740889 0.011481779 0.99425911 [39,] 0.005460221 0.010920442 0.99453978 [40,] 0.004549382 0.009098764 0.99545062 [41,] 0.003667369 0.007334737 0.99633263 [42,] 0.003019606 0.006039211 0.99698039 [43,] 0.004621405 0.009242811 0.99537859 [44,] 0.007715458 0.015430916 0.99228454 [45,] 0.011631834 0.023263669 0.98836817 [46,] 0.010315742 0.020631484 0.98968426 [47,] 0.047257186 0.094514372 0.95274281 [48,] 0.041664725 0.083329450 0.95833528 [49,] 0.037073601 0.074147203 0.96292640 [50,] 0.133790707 0.267581415 0.86620929 [51,] 0.122781688 0.245563377 0.87721831 [52,] 0.122062048 0.244124096 0.87793795 [53,] 0.117982546 0.235965091 0.88201745 [54,] 0.123970415 0.247940830 0.87602958 [55,] 0.147593333 0.295186667 0.85240667 [56,] 0.148935581 0.297871162 0.85106442 [57,] 0.161496700 0.322993399 0.83850330 [58,] 0.169790134 0.339580268 0.83020987 [59,] 0.228722367 0.457444734 0.77127763 [60,] 0.254764908 0.509529817 0.74523509 [61,] 0.582584177 0.834831646 0.41741582 [62,] 0.695816515 0.608366970 0.30418349 [63,] 0.689278458 0.621443085 0.31072154 [64,] 0.685551509 0.628896981 0.31444849 [65,] 0.701387759 0.597224481 0.29861224 [66,] 0.678362179 0.643275643 0.32163782 [67,] 0.686176818 0.627646363 0.31382318 [68,] 0.726914675 0.546170650 0.27308533 [69,] 0.868582870 0.262834260 0.13141713 [70,] 0.858357231 0.283285537 0.14164277 [71,] 0.915328865 0.169342269 0.08467113 [72,] 0.900722356 0.198555288 0.09927764 [73,] 0.918872673 0.162254654 0.08112733 [74,] 0.927924441 0.144151118 0.07207556 [75,] 0.912493649 0.175012702 0.08750635 [76,] 0.890100065 0.219799870 0.10989994 [77,] 0.870059935 0.259880131 0.12994007 [78,] 0.882007781 0.235984439 0.11799222 [79,] 0.853369914 0.293260173 0.14663009 [80,] 0.819114940 0.361770121 0.18088506 [81,] 0.776925065 0.446149870 0.22307493 [82,] 0.752012419 0.495975162 0.24798758 [83,] 0.709017808 0.581964384 0.29098219 [84,] 0.655995467 0.688009065 0.34400453 [85,] 0.623548262 0.752903476 0.37645174 [86,] 0.636939273 0.726121453 0.36306073 [87,] 0.564906252 0.870187497 0.43509375 [88,] 0.487886389 0.975772778 0.51211361 [89,] 0.462806317 0.925612635 0.53719368 [90,] 0.982103363 0.035793274 0.01789664 [91,] 0.967971722 0.064056556 0.03202828 [92,] 0.981747948 0.036504104 0.01825205 [93,] 0.964890125 0.070219751 0.03510988 [94,] 0.962759830 0.074480339 0.03724017 [95,] 0.929030183 0.141939634 0.07096982 [96,] 0.914925681 0.170148638 0.08507432 [97,] 0.847977571 0.304044859 0.15202243 [98,] 0.831214362 0.337571275 0.16878564 [99,] 0.688484830 0.623030340 0.31151517 > postscript(file="/var/www/html/rcomp/tmp/1ml911258617967.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/2hxnv1258617967.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/3vbx31258617967.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/4n3lw1258617967.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/5burc1258617967.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 = 108 Frequency = 1 1 2 3 4 5 6 -3.8852973 4.2707303 10.1297303 -1.6972664 6.4819058 -7.0250462 7 8 9 10 11 12 -9.6908390 -4.3690777 -9.8219625 -6.3311497 -25.7346731 -22.6418330 13 14 15 16 17 18 -11.5226727 -7.8220536 10.1670804 6.7763693 1.0318976 -13.8834254 19 20 21 22 23 24 -8.1844906 -17.8578858 -19.0324207 -3.8420173 -22.9681591 -29.8199003 25 26 27 28 29 30 -26.5824837 -6.6548649 -13.1301406 -3.3238518 -11.8723647 -4.6451760 31 32 33 34 35 36 1.9819782 -12.8932510 -10.5790531 -12.7111032 -26.5523470 -12.7575252 37 38 39 40 41 42 -18.8727979 -19.9813949 0.9773329 7.0186135 -8.4481536 7.7385754 43 44 45 46 47 48 5.0072430 -7.7101805 -1.9112768 -4.7657670 -17.1037159 14.5342030 49 50 51 52 53 54 14.8588474 2.6443219 30.8919721 0.7925741 1.4101331 35.5286876 55 56 57 58 59 60 -0.2792482 -4.1034676 -0.4242672 16.1741972 -7.1715257 3.2455956 61 62 63 64 65 66 1.4755399 7.8089386 28.1667025 22.9184906 50.5971280 37.6718011 67 68 69 70 71 72 6.0653909 24.0189007 29.2035644 16.1943326 13.6704169 39.1509336 73 74 75 76 77 78 38.0949707 6.7111561 27.8493764 4.9289750 19.0760111 15.7812418 79 80 81 82 83 84 -3.9477226 0.4215684 5.8924644 14.2797613 -5.6473707 1.2697330 85 86 87 88 89 90 -5.5479329 -13.2760281 4.0102701 -9.0876458 -14.3697116 13.2075850 91 92 93 94 95 96 -6.3783075 -6.0238804 -17.2995657 41.8333823 -8.4397288 -23.8450719 97 98 99 100 101 102 -5.1011042 4.6450881 -8.8542966 2.4512600 -12.0580804 2.5555214 103 104 105 106 107 108 -10.3840109 -19.8529544 -5.0152965 -9.6806441 -20.7365661 2.5375545 > postscript(file="/var/www/html/rcomp/tmp/6tmk41258617967.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 = 108 Frequency = 1 lag(myerror, k = 1) myerror 0 -3.8852973 NA 1 4.2707303 -3.8852973 2 10.1297303 4.2707303 3 -1.6972664 10.1297303 4 6.4819058 -1.6972664 5 -7.0250462 6.4819058 6 -9.6908390 -7.0250462 7 -4.3690777 -9.6908390 8 -9.8219625 -4.3690777 9 -6.3311497 -9.8219625 10 -25.7346731 -6.3311497 11 -22.6418330 -25.7346731 12 -11.5226727 -22.6418330 13 -7.8220536 -11.5226727 14 10.1670804 -7.8220536 15 6.7763693 10.1670804 16 1.0318976 6.7763693 17 -13.8834254 1.0318976 18 -8.1844906 -13.8834254 19 -17.8578858 -8.1844906 20 -19.0324207 -17.8578858 21 -3.8420173 -19.0324207 22 -22.9681591 -3.8420173 23 -29.8199003 -22.9681591 24 -26.5824837 -29.8199003 25 -6.6548649 -26.5824837 26 -13.1301406 -6.6548649 27 -3.3238518 -13.1301406 28 -11.8723647 -3.3238518 29 -4.6451760 -11.8723647 30 1.9819782 -4.6451760 31 -12.8932510 1.9819782 32 -10.5790531 -12.8932510 33 -12.7111032 -10.5790531 34 -26.5523470 -12.7111032 35 -12.7575252 -26.5523470 36 -18.8727979 -12.7575252 37 -19.9813949 -18.8727979 38 0.9773329 -19.9813949 39 7.0186135 0.9773329 40 -8.4481536 7.0186135 41 7.7385754 -8.4481536 42 5.0072430 7.7385754 43 -7.7101805 5.0072430 44 -1.9112768 -7.7101805 45 -4.7657670 -1.9112768 46 -17.1037159 -4.7657670 47 14.5342030 -17.1037159 48 14.8588474 14.5342030 49 2.6443219 14.8588474 50 30.8919721 2.6443219 51 0.7925741 30.8919721 52 1.4101331 0.7925741 53 35.5286876 1.4101331 54 -0.2792482 35.5286876 55 -4.1034676 -0.2792482 56 -0.4242672 -4.1034676 57 16.1741972 -0.4242672 58 -7.1715257 16.1741972 59 3.2455956 -7.1715257 60 1.4755399 3.2455956 61 7.8089386 1.4755399 62 28.1667025 7.8089386 63 22.9184906 28.1667025 64 50.5971280 22.9184906 65 37.6718011 50.5971280 66 6.0653909 37.6718011 67 24.0189007 6.0653909 68 29.2035644 24.0189007 69 16.1943326 29.2035644 70 13.6704169 16.1943326 71 39.1509336 13.6704169 72 38.0949707 39.1509336 73 6.7111561 38.0949707 74 27.8493764 6.7111561 75 4.9289750 27.8493764 76 19.0760111 4.9289750 77 15.7812418 19.0760111 78 -3.9477226 15.7812418 79 0.4215684 -3.9477226 80 5.8924644 0.4215684 81 14.2797613 5.8924644 82 -5.6473707 14.2797613 83 1.2697330 -5.6473707 84 -5.5479329 1.2697330 85 -13.2760281 -5.5479329 86 4.0102701 -13.2760281 87 -9.0876458 4.0102701 88 -14.3697116 -9.0876458 89 13.2075850 -14.3697116 90 -6.3783075 13.2075850 91 -6.0238804 -6.3783075 92 -17.2995657 -6.0238804 93 41.8333823 -17.2995657 94 -8.4397288 41.8333823 95 -23.8450719 -8.4397288 96 -5.1011042 -23.8450719 97 4.6450881 -5.1011042 98 -8.8542966 4.6450881 99 2.4512600 -8.8542966 100 -12.0580804 2.4512600 101 2.5555214 -12.0580804 102 -10.3840109 2.5555214 103 -19.8529544 -10.3840109 104 -5.0152965 -19.8529544 105 -9.6806441 -5.0152965 106 -20.7365661 -9.6806441 107 2.5375545 -20.7365661 108 NA 2.5375545 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] 4.2707303 -3.8852973 [2,] 10.1297303 4.2707303 [3,] -1.6972664 10.1297303 [4,] 6.4819058 -1.6972664 [5,] -7.0250462 6.4819058 [6,] -9.6908390 -7.0250462 [7,] -4.3690777 -9.6908390 [8,] -9.8219625 -4.3690777 [9,] -6.3311497 -9.8219625 [10,] -25.7346731 -6.3311497 [11,] -22.6418330 -25.7346731 [12,] -11.5226727 -22.6418330 [13,] -7.8220536 -11.5226727 [14,] 10.1670804 -7.8220536 [15,] 6.7763693 10.1670804 [16,] 1.0318976 6.7763693 [17,] -13.8834254 1.0318976 [18,] -8.1844906 -13.8834254 [19,] -17.8578858 -8.1844906 [20,] -19.0324207 -17.8578858 [21,] -3.8420173 -19.0324207 [22,] -22.9681591 -3.8420173 [23,] -29.8199003 -22.9681591 [24,] -26.5824837 -29.8199003 [25,] -6.6548649 -26.5824837 [26,] -13.1301406 -6.6548649 [27,] -3.3238518 -13.1301406 [28,] -11.8723647 -3.3238518 [29,] -4.6451760 -11.8723647 [30,] 1.9819782 -4.6451760 [31,] -12.8932510 1.9819782 [32,] -10.5790531 -12.8932510 [33,] -12.7111032 -10.5790531 [34,] -26.5523470 -12.7111032 [35,] -12.7575252 -26.5523470 [36,] -18.8727979 -12.7575252 [37,] -19.9813949 -18.8727979 [38,] 0.9773329 -19.9813949 [39,] 7.0186135 0.9773329 [40,] -8.4481536 7.0186135 [41,] 7.7385754 -8.4481536 [42,] 5.0072430 7.7385754 [43,] -7.7101805 5.0072430 [44,] -1.9112768 -7.7101805 [45,] -4.7657670 -1.9112768 [46,] -17.1037159 -4.7657670 [47,] 14.5342030 -17.1037159 [48,] 14.8588474 14.5342030 [49,] 2.6443219 14.8588474 [50,] 30.8919721 2.6443219 [51,] 0.7925741 30.8919721 [52,] 1.4101331 0.7925741 [53,] 35.5286876 1.4101331 [54,] -0.2792482 35.5286876 [55,] -4.1034676 -0.2792482 [56,] -0.4242672 -4.1034676 [57,] 16.1741972 -0.4242672 [58,] -7.1715257 16.1741972 [59,] 3.2455956 -7.1715257 [60,] 1.4755399 3.2455956 [61,] 7.8089386 1.4755399 [62,] 28.1667025 7.8089386 [63,] 22.9184906 28.1667025 [64,] 50.5971280 22.9184906 [65,] 37.6718011 50.5971280 [66,] 6.0653909 37.6718011 [67,] 24.0189007 6.0653909 [68,] 29.2035644 24.0189007 [69,] 16.1943326 29.2035644 [70,] 13.6704169 16.1943326 [71,] 39.1509336 13.6704169 [72,] 38.0949707 39.1509336 [73,] 6.7111561 38.0949707 [74,] 27.8493764 6.7111561 [75,] 4.9289750 27.8493764 [76,] 19.0760111 4.9289750 [77,] 15.7812418 19.0760111 [78,] -3.9477226 15.7812418 [79,] 0.4215684 -3.9477226 [80,] 5.8924644 0.4215684 [81,] 14.2797613 5.8924644 [82,] -5.6473707 14.2797613 [83,] 1.2697330 -5.6473707 [84,] -5.5479329 1.2697330 [85,] -13.2760281 -5.5479329 [86,] 4.0102701 -13.2760281 [87,] -9.0876458 4.0102701 [88,] -14.3697116 -9.0876458 [89,] 13.2075850 -14.3697116 [90,] -6.3783075 13.2075850 [91,] -6.0238804 -6.3783075 [92,] -17.2995657 -6.0238804 [93,] 41.8333823 -17.2995657 [94,] -8.4397288 41.8333823 [95,] -23.8450719 -8.4397288 [96,] -5.1011042 -23.8450719 [97,] 4.6450881 -5.1011042 [98,] -8.8542966 4.6450881 [99,] 2.4512600 -8.8542966 [100,] -12.0580804 2.4512600 [101,] 2.5555214 -12.0580804 [102,] -10.3840109 2.5555214 [103,] -19.8529544 -10.3840109 [104,] -5.0152965 -19.8529544 [105,] -9.6806441 -5.0152965 [106,] -20.7365661 -9.6806441 [107,] 2.5375545 -20.7365661 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 4.2707303 -3.8852973 2 10.1297303 4.2707303 3 -1.6972664 10.1297303 4 6.4819058 -1.6972664 5 -7.0250462 6.4819058 6 -9.6908390 -7.0250462 7 -4.3690777 -9.6908390 8 -9.8219625 -4.3690777 9 -6.3311497 -9.8219625 10 -25.7346731 -6.3311497 11 -22.6418330 -25.7346731 12 -11.5226727 -22.6418330 13 -7.8220536 -11.5226727 14 10.1670804 -7.8220536 15 6.7763693 10.1670804 16 1.0318976 6.7763693 17 -13.8834254 1.0318976 18 -8.1844906 -13.8834254 19 -17.8578858 -8.1844906 20 -19.0324207 -17.8578858 21 -3.8420173 -19.0324207 22 -22.9681591 -3.8420173 23 -29.8199003 -22.9681591 24 -26.5824837 -29.8199003 25 -6.6548649 -26.5824837 26 -13.1301406 -6.6548649 27 -3.3238518 -13.1301406 28 -11.8723647 -3.3238518 29 -4.6451760 -11.8723647 30 1.9819782 -4.6451760 31 -12.8932510 1.9819782 32 -10.5790531 -12.8932510 33 -12.7111032 -10.5790531 34 -26.5523470 -12.7111032 35 -12.7575252 -26.5523470 36 -18.8727979 -12.7575252 37 -19.9813949 -18.8727979 38 0.9773329 -19.9813949 39 7.0186135 0.9773329 40 -8.4481536 7.0186135 41 7.7385754 -8.4481536 42 5.0072430 7.7385754 43 -7.7101805 5.0072430 44 -1.9112768 -7.7101805 45 -4.7657670 -1.9112768 46 -17.1037159 -4.7657670 47 14.5342030 -17.1037159 48 14.8588474 14.5342030 49 2.6443219 14.8588474 50 30.8919721 2.6443219 51 0.7925741 30.8919721 52 1.4101331 0.7925741 53 35.5286876 1.4101331 54 -0.2792482 35.5286876 55 -4.1034676 -0.2792482 56 -0.4242672 -4.1034676 57 16.1741972 -0.4242672 58 -7.1715257 16.1741972 59 3.2455956 -7.1715257 60 1.4755399 3.2455956 61 7.8089386 1.4755399 62 28.1667025 7.8089386 63 22.9184906 28.1667025 64 50.5971280 22.9184906 65 37.6718011 50.5971280 66 6.0653909 37.6718011 67 24.0189007 6.0653909 68 29.2035644 24.0189007 69 16.1943326 29.2035644 70 13.6704169 16.1943326 71 39.1509336 13.6704169 72 38.0949707 39.1509336 73 6.7111561 38.0949707 74 27.8493764 6.7111561 75 4.9289750 27.8493764 76 19.0760111 4.9289750 77 15.7812418 19.0760111 78 -3.9477226 15.7812418 79 0.4215684 -3.9477226 80 5.8924644 0.4215684 81 14.2797613 5.8924644 82 -5.6473707 14.2797613 83 1.2697330 -5.6473707 84 -5.5479329 1.2697330 85 -13.2760281 -5.5479329 86 4.0102701 -13.2760281 87 -9.0876458 4.0102701 88 -14.3697116 -9.0876458 89 13.2075850 -14.3697116 90 -6.3783075 13.2075850 91 -6.0238804 -6.3783075 92 -17.2995657 -6.0238804 93 41.8333823 -17.2995657 94 -8.4397288 41.8333823 95 -23.8450719 -8.4397288 96 -5.1011042 -23.8450719 97 4.6450881 -5.1011042 98 -8.8542966 4.6450881 99 2.4512600 -8.8542966 100 -12.0580804 2.4512600 101 2.5555214 -12.0580804 102 -10.3840109 2.5555214 103 -19.8529544 -10.3840109 104 -5.0152965 -19.8529544 105 -9.6806441 -5.0152965 106 -20.7365661 -9.6806441 107 2.5375545 -20.7365661 > 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/79wdz1258617967.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/8rd6e1258617967.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/9oa4w1258617967.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/107wta1258617967.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/114j351258617967.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/12uvu11258617967.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/13yoat1258617967.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/142nnb1258617967.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/15dm5v1258617967.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/16xrgt1258617967.tab") + } > system("convert tmp/1ml911258617967.ps tmp/1ml911258617967.png") > system("convert tmp/2hxnv1258617967.ps tmp/2hxnv1258617967.png") > system("convert tmp/3vbx31258617967.ps tmp/3vbx31258617967.png") > system("convert tmp/4n3lw1258617967.ps tmp/4n3lw1258617967.png") > system("convert tmp/5burc1258617967.ps tmp/5burc1258617967.png") > system("convert tmp/6tmk41258617967.ps tmp/6tmk41258617967.png") > system("convert tmp/79wdz1258617967.ps tmp/79wdz1258617967.png") > system("convert tmp/8rd6e1258617967.ps tmp/8rd6e1258617967.png") > system("convert tmp/9oa4w1258617967.ps tmp/9oa4w1258617967.png") > system("convert tmp/107wta1258617967.ps tmp/107wta1258617967.png") > > > proc.time() user system elapsed 2.982 1.598 4.102