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