R version 2.13.0 (2011-04-13)
Copyright (C) 2011 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: i486-pc-linux-gnu (32-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> x <- array(list(10
+ ,539.51
+ ,10
+ ,407
+ ,10
+ ,723.78
+ ,10
+ ,463
+ ,10
+ ,682.06
+ ,10
+ ,556
+ ,10
+ ,283.19
+ ,10
+ ,646
+ ,10
+ ,377.18
+ ,10
+ ,702
+ ,10
+ ,486.64
+ ,11
+ ,353
+ ,10
+ ,545.38
+ ,11
+ ,346
+ ,10
+ ,554.27
+ ,11
+ ,451
+ ,10
+ ,532.54
+ ,11
+ ,964
+ ,10
+ ,324.31
+ ,12
+ ,574
+ ,10
+ ,695.25
+ ,13
+ ,031
+ ,10
+ ,827.81
+ ,13
+ ,812
+ ,10
+ ,872.48
+ ,14
+ ,544
+ ,10
+ ,971.19
+ ,14
+ ,931
+ ,11
+ ,145.65
+ ,14
+ ,886
+ ,11
+ ,234.68
+ ,16
+ ,005
+ ,11
+ ,333.88
+ ,17
+ ,064
+ ,10
+ ,997.97
+ ,15
+ ,168
+ ,11
+ ,036.89
+ ,16
+ ,050
+ ,11
+ ,257.35
+ ,15
+ ,839
+ ,11
+ ,533.59
+ ,15
+ ,137
+ ,11
+ ,963.12
+ ,14
+ ,954
+ ,12
+ ,185.15
+ ,15
+ ,648
+ ,12
+ ,377.62
+ ,15
+ ,305
+ ,12
+ ,512.89
+ ,15
+ ,579
+ ,12
+ ,631.48
+ ,16
+ ,348
+ ,12
+ ,268.53
+ ,15
+ ,928
+ ,12
+ ,754.80
+ ,16
+ ,171
+ ,13
+ ,407.75
+ ,15
+ ,937
+ ,13
+ ,480.21
+ ,15
+ ,713
+ ,13
+ ,673.28
+ ,15
+ ,594
+ ,13
+ ,239.71
+ ,15
+ ,683
+ ,13
+ ,557.69
+ ,16
+ ,438
+ ,13
+ ,901.28
+ ,17
+ ,032
+ ,13
+ ,200.58
+ ,17
+ ,696
+ ,13
+ ,406.97
+ ,17
+ ,745
+ ,12
+ ,538.12
+ ,19
+ ,394
+ ,12
+ ,419.57
+ ,20
+ ,148
+ ,12
+ ,193.88
+ ,20
+ ,108
+ ,12
+ ,656.63
+ ,18
+ ,584
+ ,12
+ ,812.48
+ ,18
+ ,441
+ ,12
+ ,056.67
+ ,18
+ ,391
+ ,11
+ ,322.38
+ ,19
+ ,178
+ ,11
+ ,530.75
+ ,18
+ ,079
+ ,11
+ ,114.08
+ ,18
+ ,483
+ ,9
+ ,181.73
+ ,19
+ ,644
+ ,8
+ ,614.55
+ ,19
+ ,195
+ ,8
+ ,595.56
+ ,19
+ ,650
+ ,8
+ ,396.20
+ ,20
+ ,830
+ ,7
+ ,690.50
+ ,23
+ ,595
+ ,7
+ ,235.47
+ ,22
+ ,937
+ ,7
+ ,992.12
+ ,21
+ ,814
+ ,8
+ ,398.37
+ ,21
+ ,928
+ ,8
+ ,593.00
+ ,21
+ ,777
+ ,8
+ ,679.75
+ ,21
+ ,383
+ ,9
+ ,374.63
+ ,21
+ ,467
+ ,9
+ ,634.97
+ ,22
+ ,052
+ ,9
+ ,857.34
+ ,22
+ ,680
+ ,10
+ ,238.83
+ ,24
+ ,320
+ ,10
+ ,433.44
+ ,24
+ ,977
+ ,10
+ ,471.24
+ ,25
+ ,204)
+ ,dim=c(4
+ ,61)
+ ,dimnames=list(c('Dow'
+ ,'Jones'
+ ,''
+ ,'Gold')
+ ,1:61))
> y <- array(NA,dim=c(4,61),dimnames=list(c('Dow','Jones','','Gold'),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
Dow Jones Gold
1 10 539.51 10 407
2 10 723.78 10 463
3 10 682.06 10 556
4 10 283.19 10 646
5 10 377.18 10 702
6 10 486.64 11 353
7 10 545.38 11 346
8 10 554.27 11 451
9 10 532.54 11 964
10 10 324.31 12 574
11 10 695.25 13 31
12 10 827.81 13 812
13 10 872.48 14 544
14 10 971.19 14 931
15 11 145.65 14 886
16 11 234.68 16 5
17 11 333.88 17 64
18 10 997.97 15 168
19 11 36.89 16 50
20 11 257.35 15 839
21 11 533.59 15 137
22 11 963.12 14 954
23 12 185.15 15 648
24 12 377.62 15 305
25 12 512.89 15 579
26 12 631.48 16 348
27 12 268.53 15 928
28 12 754.80 16 171
29 13 407.75 15 937
30 13 480.21 15 713
31 13 673.28 15 594
32 13 239.71 15 683
33 13 557.69 16 438
34 13 901.28 17 32
35 13 200.58 17 696
36 13 406.97 17 745
37 12 538.12 19 394
38 12 419.57 20 148
39 12 193.88 20 108
40 12 656.63 18 584
41 12 812.48 18 441
42 12 56.67 18 391
43 11 322.38 19 178
44 11 530.75 18 79
45 11 114.08 18 483
46 9 181.73 19 644
47 8 614.55 19 195
48 8 595.56 19 650
49 8 396.20 20 830
50 7 690.50 23 595
51 7 235.47 22 937
52 7 992.12 21 814
53 8 398.37 21 928
54 8 593.00 21 777
55 8 679.75 21 383
56 9 374.63 21 467
57 9 634.97 22 52
58 9 857.34 22 680
59 10 238.83 24 320
60 10 433.44 24 977
61 10 471.24 25 204
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Jones V3 Gold
14.205161 -0.001543 -0.141364 -0.001071
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-2.728 -1.393 -0.101 1.324 2.624
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 14.2051610 1.0712525 13.260 < 2e-16 ***
Jones -0.0015435 0.0008161 -1.891 0.06366 .
V3 -0.1413641 0.0512763 -2.757 0.00782 **
Gold -0.0010711 0.0006846 -1.565 0.12321
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1.553 on 57 degrees of freedom
Multiple R-squared: 0.178, Adjusted R-squared: 0.1348
F-statistic: 4.115 on 3 and 57 DF, p-value: 0.01035
> 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,] 1.413822e-47 2.827644e-47 1.000000000
[2,] 2.226591e-61 4.453183e-61 1.000000000
[3,] 4.933145e-77 9.866290e-77 1.000000000
[4,] 5.661980e-102 1.132396e-101 1.000000000
[5,] 1.190371e-105 2.380742e-105 1.000000000
[6,] 1.266610e-118 2.533219e-118 1.000000000
[7,] 6.066834e-141 1.213367e-140 1.000000000
[8,] 3.312279e-150 6.624558e-150 1.000000000
[9,] 7.408202e-09 1.481640e-08 0.999999993
[10,] 1.977207e-09 3.954413e-09 0.999999998
[11,] 2.702687e-10 5.405374e-10 1.000000000
[12,] 5.084489e-11 1.016898e-10 1.000000000
[13,] 1.028162e-11 2.056325e-11 1.000000000
[14,] 2.672682e-12 5.345363e-12 1.000000000
[15,] 1.723403e-11 3.446806e-11 1.000000000
[16,] 9.355655e-10 1.871131e-09 0.999999999
[17,] 2.212152e-08 4.424304e-08 0.999999978
[18,] 4.026929e-07 8.053857e-07 0.999999597
[19,] 1.596444e-06 3.192887e-06 0.999998404
[20,] 2.795153e-06 5.590305e-06 0.999997205
[21,] 1.530835e-06 3.061669e-06 0.999998469
[22,] 2.613544e-06 5.227088e-06 0.999997386
[23,] 1.536024e-05 3.072048e-05 0.999984640
[24,] 5.727320e-05 1.145464e-04 0.999942727
[25,] 2.107188e-04 4.214376e-04 0.999789281
[26,] 2.287339e-04 4.574677e-04 0.999771266
[27,] 3.034044e-04 6.068088e-04 0.999696596
[28,] 5.172087e-04 1.034417e-03 0.999482791
[29,] 4.432868e-04 8.865736e-04 0.999556713
[30,] 8.968391e-04 1.793678e-03 0.999103161
[31,] 1.471053e-03 2.942106e-03 0.998528947
[32,] 1.813363e-03 3.626726e-03 0.998186637
[33,] 1.728766e-03 3.457531e-03 0.998271234
[34,] 5.312383e-03 1.062477e-02 0.994687617
[35,] 5.433111e-02 1.086622e-01 0.945668886
[36,] 6.764316e-02 1.352863e-01 0.932356843
[37,] 8.613089e-02 1.722618e-01 0.913869106
[38,] 1.516104e-01 3.032207e-01 0.848389646
[39,] 4.088999e-01 8.177998e-01 0.591100106
[40,] 6.755761e-01 6.488478e-01 0.324423924
[41,] 7.809435e-01 4.381130e-01 0.219056502
[42,] 8.416634e-01 3.166731e-01 0.158336561
[43,] 8.516016e-01 2.967968e-01 0.148398375
[44,] 9.635978e-01 7.280437e-02 0.036402183
[45,] 9.956598e-01 8.680433e-03 0.004340217
[46,] 9.960579e-01 7.884140e-03 0.003942070
[47,] 9.912922e-01 1.741559e-02 0.008707796
[48,] 9.858402e-01 2.831968e-02 0.014159839
> postscript(file="/var/wessaorg/rcomp/tmp/1s4v01321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index')
> points(x[,1]-mysum$resid)
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/2rbch1321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/3asm81321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/4l2k51321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals')
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/5ttj71321567248.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
-1.52284162 -1.17843845 -1.14321920 -1.66247394 -1.45741784 -1.52092270
7 8 9 10 11 12
-1.43775540 -1.31156626 -0.79562271 -1.39339758 -1.26110533 -0.21995560
13 14 15 16 17 18
-0.29670315 0.27017847 -0.05224192 -0.57575137 -0.21807619 -0.36438591
19 20 21 22 23 24
-0.83283968 0.21118844 -0.11436059 1.28235816 0.89516420 0.82484779
25 26 27 28 29 30
1.32712334 1.40410262 1.32377430 1.40485893 2.54830028 2.42021167
31 32 33 34 35 36
2.59075191 2.01686664 2.38660849 2.62342928 2.25312228 2.62416982
37 38 39 40 41 42
1.73336532 1.42825250 1.03705550 1.97843378 2.06581814 0.84567033
43 44 45 46 47 48
0.16900929 0.24322332 0.03282550 -1.54894262 -2.36181763 -1.90376977
49 50 51 52 53 54
-1.87731628 -2.25068508 -2.72806497 -1.83328822 -1.62763317 -1.48896087
55 56 57 58 59 60
-1.77708302 -1.15806145 -1.05937701 -0.04348677 -0.10103069 0.90307445
61
0.27480824
> postscript(file="/var/wessaorg/rcomp/tmp/6wu021321567248.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 -1.52284162 NA
1 -1.17843845 -1.52284162
2 -1.14321920 -1.17843845
3 -1.66247394 -1.14321920
4 -1.45741784 -1.66247394
5 -1.52092270 -1.45741784
6 -1.43775540 -1.52092270
7 -1.31156626 -1.43775540
8 -0.79562271 -1.31156626
9 -1.39339758 -0.79562271
10 -1.26110533 -1.39339758
11 -0.21995560 -1.26110533
12 -0.29670315 -0.21995560
13 0.27017847 -0.29670315
14 -0.05224192 0.27017847
15 -0.57575137 -0.05224192
16 -0.21807619 -0.57575137
17 -0.36438591 -0.21807619
18 -0.83283968 -0.36438591
19 0.21118844 -0.83283968
20 -0.11436059 0.21118844
21 1.28235816 -0.11436059
22 0.89516420 1.28235816
23 0.82484779 0.89516420
24 1.32712334 0.82484779
25 1.40410262 1.32712334
26 1.32377430 1.40410262
27 1.40485893 1.32377430
28 2.54830028 1.40485893
29 2.42021167 2.54830028
30 2.59075191 2.42021167
31 2.01686664 2.59075191
32 2.38660849 2.01686664
33 2.62342928 2.38660849
34 2.25312228 2.62342928
35 2.62416982 2.25312228
36 1.73336532 2.62416982
37 1.42825250 1.73336532
38 1.03705550 1.42825250
39 1.97843378 1.03705550
40 2.06581814 1.97843378
41 0.84567033 2.06581814
42 0.16900929 0.84567033
43 0.24322332 0.16900929
44 0.03282550 0.24322332
45 -1.54894262 0.03282550
46 -2.36181763 -1.54894262
47 -1.90376977 -2.36181763
48 -1.87731628 -1.90376977
49 -2.25068508 -1.87731628
50 -2.72806497 -2.25068508
51 -1.83328822 -2.72806497
52 -1.62763317 -1.83328822
53 -1.48896087 -1.62763317
54 -1.77708302 -1.48896087
55 -1.15806145 -1.77708302
56 -1.05937701 -1.15806145
57 -0.04348677 -1.05937701
58 -0.10103069 -0.04348677
59 0.90307445 -0.10103069
60 0.27480824 0.90307445
61 NA 0.27480824
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -1.17843845 -1.52284162
[2,] -1.14321920 -1.17843845
[3,] -1.66247394 -1.14321920
[4,] -1.45741784 -1.66247394
[5,] -1.52092270 -1.45741784
[6,] -1.43775540 -1.52092270
[7,] -1.31156626 -1.43775540
[8,] -0.79562271 -1.31156626
[9,] -1.39339758 -0.79562271
[10,] -1.26110533 -1.39339758
[11,] -0.21995560 -1.26110533
[12,] -0.29670315 -0.21995560
[13,] 0.27017847 -0.29670315
[14,] -0.05224192 0.27017847
[15,] -0.57575137 -0.05224192
[16,] -0.21807619 -0.57575137
[17,] -0.36438591 -0.21807619
[18,] -0.83283968 -0.36438591
[19,] 0.21118844 -0.83283968
[20,] -0.11436059 0.21118844
[21,] 1.28235816 -0.11436059
[22,] 0.89516420 1.28235816
[23,] 0.82484779 0.89516420
[24,] 1.32712334 0.82484779
[25,] 1.40410262 1.32712334
[26,] 1.32377430 1.40410262
[27,] 1.40485893 1.32377430
[28,] 2.54830028 1.40485893
[29,] 2.42021167 2.54830028
[30,] 2.59075191 2.42021167
[31,] 2.01686664 2.59075191
[32,] 2.38660849 2.01686664
[33,] 2.62342928 2.38660849
[34,] 2.25312228 2.62342928
[35,] 2.62416982 2.25312228
[36,] 1.73336532 2.62416982
[37,] 1.42825250 1.73336532
[38,] 1.03705550 1.42825250
[39,] 1.97843378 1.03705550
[40,] 2.06581814 1.97843378
[41,] 0.84567033 2.06581814
[42,] 0.16900929 0.84567033
[43,] 0.24322332 0.16900929
[44,] 0.03282550 0.24322332
[45,] -1.54894262 0.03282550
[46,] -2.36181763 -1.54894262
[47,] -1.90376977 -2.36181763
[48,] -1.87731628 -1.90376977
[49,] -2.25068508 -1.87731628
[50,] -2.72806497 -2.25068508
[51,] -1.83328822 -2.72806497
[52,] -1.62763317 -1.83328822
[53,] -1.48896087 -1.62763317
[54,] -1.77708302 -1.48896087
[55,] -1.15806145 -1.77708302
[56,] -1.05937701 -1.15806145
[57,] -0.04348677 -1.05937701
[58,] -0.10103069 -0.04348677
[59,] 0.90307445 -0.10103069
[60,] 0.27480824 0.90307445
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -1.17843845 -1.52284162
2 -1.14321920 -1.17843845
3 -1.66247394 -1.14321920
4 -1.45741784 -1.66247394
5 -1.52092270 -1.45741784
6 -1.43775540 -1.52092270
7 -1.31156626 -1.43775540
8 -0.79562271 -1.31156626
9 -1.39339758 -0.79562271
10 -1.26110533 -1.39339758
11 -0.21995560 -1.26110533
12 -0.29670315 -0.21995560
13 0.27017847 -0.29670315
14 -0.05224192 0.27017847
15 -0.57575137 -0.05224192
16 -0.21807619 -0.57575137
17 -0.36438591 -0.21807619
18 -0.83283968 -0.36438591
19 0.21118844 -0.83283968
20 -0.11436059 0.21118844
21 1.28235816 -0.11436059
22 0.89516420 1.28235816
23 0.82484779 0.89516420
24 1.32712334 0.82484779
25 1.40410262 1.32712334
26 1.32377430 1.40410262
27 1.40485893 1.32377430
28 2.54830028 1.40485893
29 2.42021167 2.54830028
30 2.59075191 2.42021167
31 2.01686664 2.59075191
32 2.38660849 2.01686664
33 2.62342928 2.38660849
34 2.25312228 2.62342928
35 2.62416982 2.25312228
36 1.73336532 2.62416982
37 1.42825250 1.73336532
38 1.03705550 1.42825250
39 1.97843378 1.03705550
40 2.06581814 1.97843378
41 0.84567033 2.06581814
42 0.16900929 0.84567033
43 0.24322332 0.16900929
44 0.03282550 0.24322332
45 -1.54894262 0.03282550
46 -2.36181763 -1.54894262
47 -1.90376977 -2.36181763
48 -1.87731628 -1.90376977
49 -2.25068508 -1.87731628
50 -2.72806497 -2.25068508
51 -1.83328822 -2.72806497
52 -1.62763317 -1.83328822
53 -1.48896087 -1.62763317
54 -1.77708302 -1.48896087
55 -1.15806145 -1.77708302
56 -1.05937701 -1.15806145
57 -0.04348677 -1.05937701
58 -0.10103069 -0.04348677
59 0.90307445 -0.10103069
60 0.27480824 0.90307445
> plot(z,main=paste('Residual Lag plot, lowess, and regression line'), ylab='values of Residuals', xlab='lagged values of Residuals')
> lines(lowess(z))
> abline(lm(z))
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/773ft1321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/855my1321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/99iz01321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
> plot(mylm, las = 1, sub='Residual Diagnostics')
> par(opar)
> dev.off()
null device
1
> if (n > n25) {
+ postscript(file="/var/wessaorg/rcomp/tmp/10tbh41321567248.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
+ plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint')
+ grid()
+ dev.off()
+ }
null device
1
>
> #Note: the /var/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/wessaorg/rcomp/createtable")
>
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Estimated Regression Equation', 1, TRUE)
> a<-table.row.end(a)
> myeq <- colnames(x)[1]
> myeq <- paste(myeq, '[t] = ', sep='')
> for (i in 1:k){
+ if (mysum$coefficients[i,1] > 0) myeq <- paste(myeq, '+', '')
+ myeq <- paste(myeq, mysum$coefficients[i,1], sep=' ')
+ if (rownames(mysum$coefficients)[i] != '(Intercept)') {
+ myeq <- paste(myeq, rownames(mysum$coefficients)[i], sep='')
+ if (rownames(mysum$coefficients)[i] != 't') myeq <- paste(myeq, '[t]', sep='')
+ }
+ }
> myeq <- paste(myeq, ' + e[t]')
> a<-table.row.start(a)
> a<-table.element(a, myeq)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/1143cq1321567248.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,hyperlink('http://www.xycoon.com/ols1.htm','Multiple Linear Regression - Ordinary Least Squares',''), 6, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Variable',header=TRUE)
> a<-table.element(a,'Parameter',header=TRUE)
> a<-table.element(a,'S.D.',header=TRUE)
> a<-table.element(a,'T-STAT
H0: parameter = 0',header=TRUE)
> a<-table.element(a,'2-tail p-value',header=TRUE)
> a<-table.element(a,'1-tail p-value',header=TRUE)
> a<-table.row.end(a)
> for (i in 1:k){
+ a<-table.row.start(a)
+ a<-table.element(a,rownames(mysum$coefficients)[i],header=TRUE)
+ a<-table.element(a,mysum$coefficients[i,1])
+ a<-table.element(a, round(mysum$coefficients[i,2],6))
+ a<-table.element(a, round(mysum$coefficients[i,3],4))
+ a<-table.element(a, round(mysum$coefficients[i,4],6))
+ a<-table.element(a, round(mysum$coefficients[i,4]/2,6))
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/126wvz1321567248.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Regression Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple R',1,TRUE)
> a<-table.element(a, sqrt(mysum$r.squared))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'R-squared',1,TRUE)
> a<-table.element(a, mysum$r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Adjusted R-squared',1,TRUE)
> a<-table.element(a, mysum$adj.r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (value)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[3])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'p-value',1,TRUE)
> a<-table.element(a, 1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3]))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Residual Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Residual Standard Deviation',1,TRUE)
> a<-table.element(a, mysum$sigma)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Sum Squared Residuals',1,TRUE)
> a<-table.element(a, sum(myerror*myerror))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/13qe7m1321567248.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Actuals, Interpolation, and Residuals', 4, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Time or Index', 1, TRUE)
> a<-table.element(a, 'Actuals', 1, TRUE)
> a<-table.element(a, 'Interpolation
Forecast', 1, TRUE)
> a<-table.element(a, 'Residuals
Prediction Error', 1, TRUE)
> a<-table.row.end(a)
> for (i in 1:n) {
+ a<-table.row.start(a)
+ a<-table.element(a,i, 1, TRUE)
+ a<-table.element(a,x[i])
+ a<-table.element(a,x[i]-mysum$resid[i])
+ a<-table.element(a,mysum$resid[i])
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/1498pg1321567248.tab")
> if (n > n25) {
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'p-values',header=TRUE)
+ a<-table.element(a,'Alternative Hypothesis',3,header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'breakpoint index',header=TRUE)
+ a<-table.element(a,'greater',header=TRUE)
+ a<-table.element(a,'2-sided',header=TRUE)
+ a<-table.element(a,'less',header=TRUE)
+ a<-table.row.end(a)
+ for (mypoint in kp3:nmkm3) {
+ a<-table.row.start(a)
+ a<-table.element(a,mypoint,header=TRUE)
+ a<-table.element(a,gqarr[mypoint-kp3+1,1])
+ a<-table.element(a,gqarr[mypoint-kp3+1,2])
+ a<-table.element(a,gqarr[mypoint-kp3+1,3])
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/15b4ag1321567248.tab")
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Meta Analysis of Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'Description',header=TRUE)
+ a<-table.element(a,'# significant tests',header=TRUE)
+ a<-table.element(a,'% significant tests',header=TRUE)
+ a<-table.element(a,'OK/NOK',header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'1% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant1)
+ a<-table.element(a,numsignificant1/numgqtests)
+ if (numsignificant1/numgqtests < 0.01) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'5% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant5)
+ a<-table.element(a,numsignificant5/numgqtests)
+ if (numsignificant5/numgqtests < 0.05) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'10% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant10)
+ a<-table.element(a,numsignificant10/numgqtests)
+ if (numsignificant10/numgqtests < 0.1) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/169wfi1321567248.tab")
+ }
>
> try(system("convert tmp/1s4v01321567248.ps tmp/1s4v01321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/2rbch1321567248.ps tmp/2rbch1321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/3asm81321567248.ps tmp/3asm81321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/4l2k51321567248.ps tmp/4l2k51321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/5ttj71321567248.ps tmp/5ttj71321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/6wu021321567248.ps tmp/6wu021321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/773ft1321567248.ps tmp/773ft1321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/855my1321567248.ps tmp/855my1321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/99iz01321567248.ps tmp/99iz01321567248.png",intern=TRUE))
character(0)
> try(system("convert tmp/10tbh41321567248.ps tmp/10tbh41321567248.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.371 0.514 3.901