R version 2.12.0 (2010-10-15)
Copyright (C) 2010 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(10554.27
+ ,2.08
+ ,83.9
+ ,61.2
+ ,11451
+ ,63.96
+ ,2.17
+ ,69
+ ,10532.54
+ ,2.09
+ ,85.6
+ ,62
+ ,11964
+ ,63.77
+ ,2.23
+ ,67
+ ,10324.31
+ ,2.07
+ ,87.5
+ ,65.1
+ ,12574
+ ,59.15
+ ,2.17
+ ,69
+ ,10695.25
+ ,2.04
+ ,88.5
+ ,63.2
+ ,13031
+ ,56.12
+ ,2.39
+ ,79
+ ,10827.81
+ ,2.35
+ ,91
+ ,66.3
+ ,13812
+ ,57.42
+ ,2.6
+ ,104
+ ,10872.48
+ ,2.33
+ ,90.6
+ ,61.9
+ ,14544
+ ,63.52
+ ,2.67
+ ,117
+ ,10971.19
+ ,2.37
+ ,91.2
+ ,62.1
+ ,14931
+ ,61.71
+ ,2.63
+ ,73
+ ,11145.65
+ ,2.59
+ ,93.2
+ ,66.3
+ ,14886
+ ,63.01
+ ,2.85
+ ,97
+ ,11234.68
+ ,2.62
+ ,90.1
+ ,72
+ ,16005
+ ,68.18
+ ,3.1
+ ,124
+ ,11333.88
+ ,2.6
+ ,95
+ ,65.3
+ ,17064
+ ,72.03
+ ,3.2
+ ,129
+ ,10997.97
+ ,2.83
+ ,95.4
+ ,67.6
+ ,15168
+ ,69.75
+ ,3.21
+ ,122
+ ,11036.89
+ ,2.78
+ ,93.7
+ ,70.5
+ ,16050
+ ,74.41
+ ,3.36
+ ,113
+ ,11257.35
+ ,3.01
+ ,93.9
+ ,74.2
+ ,15839
+ ,74.33
+ ,3.41
+ ,131
+ ,11533.59
+ ,3.06
+ ,92.5
+ ,77.8
+ ,15137
+ ,64.24
+ ,3.32
+ ,155
+ ,11963.12
+ ,3.33
+ ,89.2
+ ,78.5
+ ,14954
+ ,60.03
+ ,3.24
+ ,161
+ ,12185.15
+ ,3.32
+ ,93.3
+ ,77.8
+ ,15648
+ ,59.44
+ ,3.24
+ ,141
+ ,12377.62
+ ,3.6
+ ,93
+ ,81.4
+ ,15305
+ ,62.5
+ ,3.26
+ ,116
+ ,12512.89
+ ,3.57
+ ,96.1
+ ,84.5
+ ,15579
+ ,55.04
+ ,3.48
+ ,197
+ ,12631.48
+ ,3.57
+ ,96.7
+ ,88
+ ,16348
+ ,58.34
+ ,3.61
+ ,163
+ ,12268.53
+ ,3.83
+ ,97.6
+ ,93.9
+ ,15928
+ ,61.92
+ ,3.68
+ ,154
+ ,12754.8
+ ,3.84
+ ,102.6
+ ,98.9
+ ,16171
+ ,67.65
+ ,3.67
+ ,143
+ ,13407.75
+ ,3.8
+ ,107.6
+ ,96.7
+ ,15937
+ ,67.68
+ ,3.71
+ ,165
+ ,13480.21
+ ,4.07
+ ,103.5
+ ,98.9
+ ,15713
+ ,70.3
+ ,4.04
+ ,169
+ ,13673.28
+ ,4.05
+ ,100.8
+ ,102.2
+ ,15594
+ ,75.26
+ ,4.2
+ ,194
+ ,13239.71
+ ,4.272
+ ,94.5
+ ,105.4
+ ,15683
+ ,71.44
+ ,4.16
+ ,185
+ ,13557.69
+ ,3.858
+ ,100.1
+ ,105.1
+ ,16438
+ ,76.36
+ ,4.01
+ ,160
+ ,13901.28
+ ,4.067
+ ,97.4
+ ,116.6
+ ,17032
+ ,81.71
+ ,4.03
+ ,135
+ ,13200.58
+ ,3.964
+ ,103
+ ,112
+ ,17696
+ ,92.6
+ ,4
+ ,125
+ ,13406.97
+ ,3.782
+ ,100.2
+ ,108.8
+ ,17745
+ ,90.6
+ ,3.9
+ ,131
+ ,12538.12
+ ,4.114
+ ,100.2
+ ,106.9
+ ,19394
+ ,92.23
+ ,3.93
+ ,73
+ ,12419.57
+ ,4.009
+ ,99
+ ,109.5
+ ,20148
+ ,94.09
+ ,3.68
+ ,69
+ ,12193.88
+ ,4.025
+ ,102.4
+ ,106.7
+ ,20108
+ ,102.79
+ ,3.59
+ ,54
+ ,12656.63
+ ,4.082
+ ,99
+ ,118.9
+ ,18584
+ ,109.65
+ ,3.55
+ ,73
+ ,12812.48
+ ,4.044
+ ,103.7
+ ,117.5
+ ,18441
+ ,124.05
+ ,3.88
+ ,54
+ ,12056.67
+ ,3.916
+ ,103.4
+ ,113.7
+ ,18391
+ ,132.69
+ ,4.61
+ ,58
+ ,11322.38
+ ,4.289
+ ,95.3
+ ,119.6
+ ,19178
+ ,135.81
+ ,4.87
+ ,45
+ ,11530.75
+ ,4.296
+ ,93.6
+ ,120.6
+ ,18079
+ ,116.07
+ ,4.9
+ ,41
+ ,11114.08
+ ,4.193
+ ,102.4
+ ,117.5
+ ,18483
+ ,101.42
+ ,4.54
+ ,74
+ ,9181.73
+ ,3.48
+ ,110.5
+ ,120.3
+ ,19644
+ ,75.73
+ ,4.47
+ ,31
+ ,8614.55
+ ,2.934
+ ,109.1
+ ,119.8
+ ,19195
+ ,55.48
+ ,4.04
+ ,23
+ ,8595.56
+ ,2.221
+ ,100.9
+ ,108
+ ,19650
+ ,43.8
+ ,4.08
+ ,123
+ ,8396.2
+ ,1.211
+ ,108.1
+ ,98.8
+ ,20830
+ ,45.29
+ ,3.47
+ ,46
+ ,7690.5
+ ,1.28
+ ,105
+ ,94.6
+ ,23595
+ ,44.01
+ ,3.27
+ ,162
+ ,7235.47
+ ,0.96
+ ,111.5
+ ,84.6
+ ,22937
+ ,47.48
+ ,2.78
+ ,62
+ ,7992.12
+ ,0.5
+ ,109.5
+ ,84.4
+ ,21814
+ ,51.07
+ ,2.83
+ ,166
+ ,8398.37
+ ,0.687
+ ,110.5
+ ,79.1
+ ,21928
+ ,57.84
+ ,2.99
+ ,441
+ ,8593
+ ,0.344
+ ,114
+ ,73.3
+ ,21777
+ ,69.04
+ ,2.96
+ ,313
+ ,8679.75
+ ,0.346
+ ,108.2
+ ,74.3
+ ,21383
+ ,65.61
+ ,2.99
+ ,224
+ ,9374.63
+ ,0.334
+ ,110.3
+ ,67.8
+ ,21467
+ ,72.87
+ ,2.98
+ ,175
+ ,9634.97
+ ,0.34
+ ,111.8
+ ,64.8
+ ,22052
+ ,68.41
+ ,2.98
+ ,305
+ ,9857.34
+ ,0.328
+ ,107.5
+ ,66.5
+ ,22680
+ ,73.25
+ ,2.92
+ ,400
+ ,10238.83
+ ,0.344
+ ,114.1
+ ,57.7
+ ,24320
+ ,77.43
+ ,2.62
+ ,299
+ ,10433.44
+ ,0.341
+ ,113.8
+ ,53.8
+ ,24977
+ ,75.28
+ ,2.68
+ ,204
+ ,10471.24
+ ,0.32
+ ,114.5
+ ,51.8
+ ,25204
+ ,77.33
+ ,2.59
+ ,164
+ ,10214.51
+ ,0.314
+ ,114.8
+ ,50.9
+ ,25739
+ ,74.31
+ ,2.53
+ ,133
+ ,10677.52
+ ,0.325
+ ,117.8
+ ,49
+ ,26434
+ ,79.7
+ ,2.4
+ ,83
+ ,11052.15
+ ,0.339
+ ,116.7
+ ,48.1
+ ,27525
+ ,85.47
+ ,2.16
+ ,71
+ ,10500.19
+ ,0.329
+ ,122.8
+ ,42.6
+ ,30695
+ ,77.98
+ ,2.03
+ ,46
+ ,10159.27
+ ,0.48
+ ,122.3
+ ,40.9
+ ,32436
+ ,75.69
+ ,2.05
+ ,68
+ ,10222.24
+ ,0.399
+ ,115
+ ,43.3
+ ,30160
+ ,75.2
+ ,1.91
+ ,51
+ ,10350.4
+ ,0.37
+ ,118.5
+ ,43.7
+ ,30236
+ ,77.21
+ ,1.9
+ ,38)
+ ,dim=c(8
+ ,61)
+ ,dimnames=list(c('DowJones'
+ ,'Eonia'
+ ,'deposits'
+ ,'2JAAR'
+ ,'Goudkoers'
+ ,'Brent'
+ ,'gewrentevoet'
+ ,'kasbons')
+ ,1:61))
> y <- array(NA,dim=c(8,61),dimnames=list(c('DowJones','Eonia','deposits','2JAAR','Goudkoers','Brent','gewrentevoet','kasbons'),1:61))
> 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
DowJones Eonia deposits 2JAAR Goudkoers Brent gewrentevoet kasbons
1 10554.27 2.080 83.9 61.2 11451 63.96 2.17 69
2 10532.54 2.090 85.6 62.0 11964 63.77 2.23 67
3 10324.31 2.070 87.5 65.1 12574 59.15 2.17 69
4 10695.25 2.040 88.5 63.2 13031 56.12 2.39 79
5 10827.81 2.350 91.0 66.3 13812 57.42 2.60 104
6 10872.48 2.330 90.6 61.9 14544 63.52 2.67 117
7 10971.19 2.370 91.2 62.1 14931 61.71 2.63 73
8 11145.65 2.590 93.2 66.3 14886 63.01 2.85 97
9 11234.68 2.620 90.1 72.0 16005 68.18 3.10 124
10 11333.88 2.600 95.0 65.3 17064 72.03 3.20 129
11 10997.97 2.830 95.4 67.6 15168 69.75 3.21 122
12 11036.89 2.780 93.7 70.5 16050 74.41 3.36 113
13 11257.35 3.010 93.9 74.2 15839 74.33 3.41 131
14 11533.59 3.060 92.5 77.8 15137 64.24 3.32 155
15 11963.12 3.330 89.2 78.5 14954 60.03 3.24 161
16 12185.15 3.320 93.3 77.8 15648 59.44 3.24 141
17 12377.62 3.600 93.0 81.4 15305 62.50 3.26 116
18 12512.89 3.570 96.1 84.5 15579 55.04 3.48 197
19 12631.48 3.570 96.7 88.0 16348 58.34 3.61 163
20 12268.53 3.830 97.6 93.9 15928 61.92 3.68 154
21 12754.80 3.840 102.6 98.9 16171 67.65 3.67 143
22 13407.75 3.800 107.6 96.7 15937 67.68 3.71 165
23 13480.21 4.070 103.5 98.9 15713 70.30 4.04 169
24 13673.28 4.050 100.8 102.2 15594 75.26 4.20 194
25 13239.71 4.272 94.5 105.4 15683 71.44 4.16 185
26 13557.69 3.858 100.1 105.1 16438 76.36 4.01 160
27 13901.28 4.067 97.4 116.6 17032 81.71 4.03 135
28 13200.58 3.964 103.0 112.0 17696 92.60 4.00 125
29 13406.97 3.782 100.2 108.8 17745 90.60 3.90 131
30 12538.12 4.114 100.2 106.9 19394 92.23 3.93 73
31 12419.57 4.009 99.0 109.5 20148 94.09 3.68 69
32 12193.88 4.025 102.4 106.7 20108 102.79 3.59 54
33 12656.63 4.082 99.0 118.9 18584 109.65 3.55 73
34 12812.48 4.044 103.7 117.5 18441 124.05 3.88 54
35 12056.67 3.916 103.4 113.7 18391 132.69 4.61 58
36 11322.38 4.289 95.3 119.6 19178 135.81 4.87 45
37 11530.75 4.296 93.6 120.6 18079 116.07 4.90 41
38 11114.08 4.193 102.4 117.5 18483 101.42 4.54 74
39 9181.73 3.480 110.5 120.3 19644 75.73 4.47 31
40 8614.55 2.934 109.1 119.8 19195 55.48 4.04 23
41 8595.56 2.221 100.9 108.0 19650 43.80 4.08 123
42 8396.20 1.211 108.1 98.8 20830 45.29 3.47 46
43 7690.50 1.280 105.0 94.6 23595 44.01 3.27 162
44 7235.47 0.960 111.5 84.6 22937 47.48 2.78 62
45 7992.12 0.500 109.5 84.4 21814 51.07 2.83 166
46 8398.37 0.687 110.5 79.1 21928 57.84 2.99 441
47 8593.00 0.344 114.0 73.3 21777 69.04 2.96 313
48 8679.75 0.346 108.2 74.3 21383 65.61 2.99 224
49 9374.63 0.334 110.3 67.8 21467 72.87 2.98 175
50 9634.97 0.340 111.8 64.8 22052 68.41 2.98 305
51 9857.34 0.328 107.5 66.5 22680 73.25 2.92 400
52 10238.83 0.344 114.1 57.7 24320 77.43 2.62 299
53 10433.44 0.341 113.8 53.8 24977 75.28 2.68 204
54 10471.24 0.320 114.5 51.8 25204 77.33 2.59 164
55 10214.51 0.314 114.8 50.9 25739 74.31 2.53 133
56 10677.52 0.325 117.8 49.0 26434 79.70 2.40 83
57 11052.15 0.339 116.7 48.1 27525 85.47 2.16 71
58 10500.19 0.329 122.8 42.6 30695 77.98 2.03 46
59 10159.27 0.480 122.3 40.9 32436 75.69 2.05 68
60 10222.24 0.399 115.0 43.3 30160 75.20 1.91 51
61 10350.40 0.370 118.5 43.7 30236 77.21 1.90 38
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Eonia deposits `2JAAR` Goudkoers
3294.69304 1902.00849 75.84196 -40.72785 -0.03146
Brent gewrentevoet kasbons
19.71174 -859.83688 5.39966
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1622.50 -283.83 -29.23 314.66 1894.40
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3294.69304 1695.74963 1.943 0.057347 .
Eonia 1902.00849 156.88106 12.124 < 2e-16 ***
deposits 75.84196 25.39959 2.986 0.004272 **
`2JAAR` -40.72785 10.00270 -4.072 0.000156 ***
Goudkoers -0.03146 0.05603 -0.562 0.576772
Brent 19.71174 5.88631 3.349 0.001500 **
gewrentevoet -859.83688 345.90073 -2.486 0.016115 *
kasbons 5.39966 1.28034 4.217 9.69e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 676.7 on 53 degrees of freedom
Multiple R-squared: 0.85, Adjusted R-squared: 0.8302
F-statistic: 42.9 on 7 and 53 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,] 8.695336e-03 1.739067e-02 0.9913046641
[2,] 2.509473e-03 5.018947e-03 0.9974905267
[3,] 5.704009e-04 1.140802e-03 0.9994295991
[4,] 1.316016e-04 2.632032e-04 0.9998683984
[5,] 2.615103e-05 5.230205e-05 0.9999738490
[6,] 4.606730e-05 9.213461e-05 0.9999539327
[7,] 2.383379e-05 4.766758e-05 0.9999761662
[8,] 2.408895e-05 4.817790e-05 0.9999759110
[9,] 2.495382e-05 4.990764e-05 0.9999750462
[10,] 3.773143e-05 7.546287e-05 0.9999622686
[11,] 7.841736e-05 1.568347e-04 0.9999215826
[12,] 5.213395e-04 1.042679e-03 0.9994786605
[13,] 1.546104e-03 3.092207e-03 0.9984538964
[14,] 2.937896e-03 5.875792e-03 0.9970621042
[15,] 3.098319e-03 6.196637e-03 0.9969016814
[16,] 3.083565e-03 6.167129e-03 0.9969164353
[17,] 8.214432e-03 1.642886e-02 0.9917855676
[18,] 1.341154e-02 2.682307e-02 0.9865884641
[19,] 2.903977e-02 5.807954e-02 0.9709602292
[20,] 2.889196e-02 5.778391e-02 0.9711080443
[21,] 2.498254e-02 4.996509e-02 0.9750174568
[22,] 2.035456e-02 4.070911e-02 0.9796454430
[23,] 1.622011e-02 3.244022e-02 0.9837798895
[24,] 7.896736e-02 1.579347e-01 0.9210326433
[25,] 1.093972e-01 2.187945e-01 0.8906027642
[26,] 1.869169e-01 3.738337e-01 0.8130831419
[27,] 1.782266e-01 3.564532e-01 0.8217734233
[28,] 5.111406e-01 9.777187e-01 0.4888593518
[29,] 8.329500e-01 3.341001e-01 0.1670500362
[30,] 7.914776e-01 4.170448e-01 0.2085224220
[31,] 7.385714e-01 5.228573e-01 0.2614286251
[32,] 9.432108e-01 1.135785e-01 0.0567892266
[33,] 9.857287e-01 2.854254e-02 0.0142712710
[34,] 9.749714e-01 5.005719e-02 0.0250285926
[35,] 9.978777e-01 4.244585e-03 0.0021222926
[36,] 9.996274e-01 7.452609e-04 0.0003726304
[37,] 9.996611e-01 6.778602e-04 0.0003389301
[38,] 9.987466e-01 2.506841e-03 0.0012534203
[39,] 9.994272e-01 1.145630e-03 0.0005728148
[40,] 9.959729e-01 8.054166e-03 0.0040270832
> postscript(file="/var/www/rcomp/tmp/1247n1293037994.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/www/rcomp/tmp/2247n1293037994.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/www/rcomp/tmp/3247n1293037994.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/www/rcomp/tmp/4vd7q1293037994.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/www/rcomp/tmp/5vd7q1293037994.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 = 61
Frequency = 1
1 2 3 4 5 6
25.596514 -29.227149 -169.389009 314.659177 -161.230247 -334.603338
7 8 9 10 11 12
-98.286988 -290.364728 211.329924 -279.517252 -957.867660 -463.338676
13 14 15 16 17 18
-604.079901 -200.314361 -29.494987 13.552339 -76.088720 -85.150284
19 20 21 22 23 24
384.995945 -275.436795 -38.262173 129.604201 292.530170 763.883466
25 26 27 28 29 30
608.503632 1209.770944 1894.404482 611.980599 1169.148024 -49.823189
31 32 33 34 35 36
21.939427 -775.224405 13.693770 -73.664670 -283.831183 -615.941930
37 38 39 40 41 42
150.698517 -1049.972140 -1411.551406 -795.916769 421.516788 1121.456898
43 44 45 46 47 48
-337.529750 -1054.624610 95.815972 -1622.504453 -837.314570 287.822507
49 50 51 52 53 54
697.048287 114.393905 114.741671 -136.542702 555.297216 603.828657
55 56 57 58 59 60
491.264769 702.269010 876.066857 -72.919616 -734.037544 43.853567
61
2.383971
> postscript(file="/var/www/rcomp/tmp/6vd7q1293037994.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 = 61
Frequency = 1
lag(myerror, k = 1) myerror
0 25.596514 NA
1 -29.227149 25.596514
2 -169.389009 -29.227149
3 314.659177 -169.389009
4 -161.230247 314.659177
5 -334.603338 -161.230247
6 -98.286988 -334.603338
7 -290.364728 -98.286988
8 211.329924 -290.364728
9 -279.517252 211.329924
10 -957.867660 -279.517252
11 -463.338676 -957.867660
12 -604.079901 -463.338676
13 -200.314361 -604.079901
14 -29.494987 -200.314361
15 13.552339 -29.494987
16 -76.088720 13.552339
17 -85.150284 -76.088720
18 384.995945 -85.150284
19 -275.436795 384.995945
20 -38.262173 -275.436795
21 129.604201 -38.262173
22 292.530170 129.604201
23 763.883466 292.530170
24 608.503632 763.883466
25 1209.770944 608.503632
26 1894.404482 1209.770944
27 611.980599 1894.404482
28 1169.148024 611.980599
29 -49.823189 1169.148024
30 21.939427 -49.823189
31 -775.224405 21.939427
32 13.693770 -775.224405
33 -73.664670 13.693770
34 -283.831183 -73.664670
35 -615.941930 -283.831183
36 150.698517 -615.941930
37 -1049.972140 150.698517
38 -1411.551406 -1049.972140
39 -795.916769 -1411.551406
40 421.516788 -795.916769
41 1121.456898 421.516788
42 -337.529750 1121.456898
43 -1054.624610 -337.529750
44 95.815972 -1054.624610
45 -1622.504453 95.815972
46 -837.314570 -1622.504453
47 287.822507 -837.314570
48 697.048287 287.822507
49 114.393905 697.048287
50 114.741671 114.393905
51 -136.542702 114.741671
52 555.297216 -136.542702
53 603.828657 555.297216
54 491.264769 603.828657
55 702.269010 491.264769
56 876.066857 702.269010
57 -72.919616 876.066857
58 -734.037544 -72.919616
59 43.853567 -734.037544
60 2.383971 43.853567
61 NA 2.383971
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -29.227149 25.59651
[2,] -169.389009 -29.22715
[3,] 314.659177 -169.38901
[4,] -161.230247 314.65918
[5,] -334.603338 -161.23025
[6,] -98.286988 -334.60334
[7,] -290.364728 -98.28699
[8,] 211.329924 -290.36473
[9,] -279.517252 211.32992
[10,] -957.867660 -279.51725
[11,] -463.338676 -957.86766
[12,] -604.079901 -463.33868
[13,] -200.314361 -604.07990
[14,] -29.494987 -200.31436
[15,] 13.552339 -29.49499
[16,] -76.088720 13.55234
[17,] -85.150284 -76.08872
[18,] 384.995945 -85.15028
[19,] -275.436795 384.99594
[20,] -38.262173 -275.43679
[21,] 129.604201 -38.26217
[22,] 292.530170 129.60420
[23,] 763.883466 292.53017
[24,] 608.503632 763.88347
[25,] 1209.770944 608.50363
[26,] 1894.404482 1209.77094
[27,] 611.980599 1894.40448
[28,] 1169.148024 611.98060
[29,] -49.823189 1169.14802
[30,] 21.939427 -49.82319
[31,] -775.224405 21.93943
[32,] 13.693770 -775.22440
[33,] -73.664670 13.69377
[34,] -283.831183 -73.66467
[35,] -615.941930 -283.83118
[36,] 150.698517 -615.94193
[37,] -1049.972140 150.69852
[38,] -1411.551406 -1049.97214
[39,] -795.916769 -1411.55141
[40,] 421.516788 -795.91677
[41,] 1121.456898 421.51679
[42,] -337.529750 1121.45690
[43,] -1054.624610 -337.52975
[44,] 95.815972 -1054.62461
[45,] -1622.504453 95.81597
[46,] -837.314570 -1622.50445
[47,] 287.822507 -837.31457
[48,] 697.048287 287.82251
[49,] 114.393905 697.04829
[50,] 114.741671 114.39390
[51,] -136.542702 114.74167
[52,] 555.297216 -136.54270
[53,] 603.828657 555.29722
[54,] 491.264769 603.82866
[55,] 702.269010 491.26477
[56,] 876.066857 702.26901
[57,] -72.919616 876.06686
[58,] -734.037544 -72.91962
[59,] 43.853567 -734.03754
[60,] 2.383971 43.85357
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -29.227149 25.59651
2 -169.389009 -29.22715
3 314.659177 -169.38901
4 -161.230247 314.65918
5 -334.603338 -161.23025
6 -98.286988 -334.60334
7 -290.364728 -98.28699
8 211.329924 -290.36473
9 -279.517252 211.32992
10 -957.867660 -279.51725
11 -463.338676 -957.86766
12 -604.079901 -463.33868
13 -200.314361 -604.07990
14 -29.494987 -200.31436
15 13.552339 -29.49499
16 -76.088720 13.55234
17 -85.150284 -76.08872
18 384.995945 -85.15028
19 -275.436795 384.99594
20 -38.262173 -275.43679
21 129.604201 -38.26217
22 292.530170 129.60420
23 763.883466 292.53017
24 608.503632 763.88347
25 1209.770944 608.50363
26 1894.404482 1209.77094
27 611.980599 1894.40448
28 1169.148024 611.98060
29 -49.823189 1169.14802
30 21.939427 -49.82319
31 -775.224405 21.93943
32 13.693770 -775.22440
33 -73.664670 13.69377
34 -283.831183 -73.66467
35 -615.941930 -283.83118
36 150.698517 -615.94193
37 -1049.972140 150.69852
38 -1411.551406 -1049.97214
39 -795.916769 -1411.55141
40 421.516788 -795.91677
41 1121.456898 421.51679
42 -337.529750 1121.45690
43 -1054.624610 -337.52975
44 95.815972 -1054.62461
45 -1622.504453 95.81597
46 -837.314570 -1622.50445
47 287.822507 -837.31457
48 697.048287 287.82251
49 114.393905 697.04829
50 114.741671 114.39390
51 -136.542702 114.74167
52 555.297216 -136.54270
53 603.828657 555.29722
54 491.264769 603.82866
55 702.269010 491.26477
56 876.066857 702.26901
57 -72.919616 876.06686
58 -734.037544 -72.91962
59 43.853567 -734.03754
60 2.383971 43.85357
> 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/rcomp/tmp/7n4ot1293037994.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/www/rcomp/tmp/8ge5w1293037994.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/www/rcomp/tmp/9ge5w1293037994.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/www/rcomp/tmp/10ge5w1293037994.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/www/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/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/rcomp/tmp/111wm21293037994.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/rcomp/tmp/125x2q1293037994.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/rcomp/tmp/13txh11293037994.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/rcomp/tmp/14m7y41293037994.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/rcomp/tmp/1587xa1293037994.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/rcomp/tmp/16mhu11293037994.tab")
+ }
>
> try(system("convert tmp/1247n1293037994.ps tmp/1247n1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/2247n1293037994.ps tmp/2247n1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/3247n1293037994.ps tmp/3247n1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/4vd7q1293037994.ps tmp/4vd7q1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/5vd7q1293037994.ps tmp/5vd7q1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/6vd7q1293037994.ps tmp/6vd7q1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/7n4ot1293037994.ps tmp/7n4ot1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/8ge5w1293037994.ps tmp/8ge5w1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/9ge5w1293037994.ps tmp/9ge5w1293037994.png",intern=TRUE))
character(0)
> try(system("convert tmp/10ge5w1293037994.ps tmp/10ge5w1293037994.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.210 1.660 4.879