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(2284
+ ,33
+ ,41
+ ,76403
+ ,194493
+ ,3160
+ ,108
+ ,90
+ ,108094
+ ,530670
+ ,4150
+ ,150
+ ,136
+ ,134759
+ ,518365
+ ,7285
+ ,115
+ ,97
+ ,188873
+ ,491303
+ ,1134
+ ,162
+ ,63
+ ,146216
+ ,527021
+ ,4658
+ ,158
+ ,114
+ ,156608
+ ,233773
+ ,2384
+ ,97
+ ,77
+ ,61348
+ ,405972
+ ,3748
+ ,9
+ ,6
+ ,50350
+ ,652925
+ ,5371
+ ,66
+ ,47
+ ,87720
+ ,446211
+ ,1285
+ ,107
+ ,51
+ ,99489
+ ,341340
+ ,9327
+ ,101
+ ,85
+ ,87419
+ ,387699
+ ,5565
+ ,47
+ ,43
+ ,94355
+ ,493408
+ ,1528
+ ,38
+ ,32
+ ,60326
+ ,146494
+ ,3122
+ ,34
+ ,25
+ ,94670
+ ,414462
+ ,7561
+ ,87
+ ,77
+ ,82425
+ ,364304
+ ,2675
+ ,79
+ ,54
+ ,59017
+ ,355178
+ ,13253
+ ,947
+ ,251
+ ,90829
+ ,357760
+ ,880
+ ,74
+ ,15
+ ,80791
+ ,261216
+ ,2053
+ ,53
+ ,44
+ ,100423
+ ,397144
+ ,1424
+ ,94
+ ,73
+ ,131116
+ ,374943
+ ,4036
+ ,63
+ ,85
+ ,100269
+ ,424898
+ ,3045
+ ,58
+ ,49
+ ,27330
+ ,202055
+ ,5119
+ ,49
+ ,38
+ ,39039
+ ,378525
+ ,1431
+ ,34
+ ,35
+ ,106885
+ ,310768
+ ,554
+ ,11
+ ,9
+ ,79285
+ ,325738
+ ,1975
+ ,35
+ ,34
+ ,118881
+ ,394510
+ ,1765
+ ,20
+ ,20
+ ,77623
+ ,247060
+ ,1012
+ ,47
+ ,29
+ ,114768
+ ,368078
+ ,810
+ ,43
+ ,11
+ ,74015
+ ,236761
+ ,1280
+ ,117
+ ,52
+ ,69465
+ ,312378
+ ,666
+ ,171
+ ,13
+ ,117869
+ ,339836
+ ,1380
+ ,26
+ ,29
+ ,60982
+ ,347385
+ ,4677
+ ,75
+ ,66
+ ,90131
+ ,426280
+ ,876
+ ,59
+ ,33
+ ,138971
+ ,352850
+ ,814
+ ,18
+ ,15
+ ,39625
+ ,301881
+ ,514
+ ,15
+ ,15
+ ,102725
+ ,377516
+ ,5692
+ ,72
+ ,68
+ ,64239
+ ,357312
+ ,3642
+ ,86
+ ,100
+ ,90262
+ ,458343
+ ,540
+ ,14
+ ,13
+ ,103960
+ ,354228
+ ,2099
+ ,64
+ ,45
+ ,106611
+ ,308636
+ ,567
+ ,11
+ ,14
+ ,103345
+ ,386212
+ ,2001
+ ,52
+ ,36
+ ,95551
+ ,393343
+ ,2949
+ ,41
+ ,40
+ ,82903
+ ,378509
+ ,2253
+ ,99
+ ,68
+ ,63593
+ ,452469
+ ,6533
+ ,75
+ ,29
+ ,126910
+ ,364839
+ ,1889
+ ,45
+ ,43
+ ,37527
+ ,358649
+ ,3055
+ ,43
+ ,30
+ ,60247
+ ,376641
+ ,272
+ ,8
+ ,9
+ ,112995
+ ,429112
+ ,1414
+ ,198
+ ,22
+ ,70184
+ ,330546
+ ,2564
+ ,22
+ ,19
+ ,130140
+ ,403560
+ ,1383
+ ,11
+ ,9
+ ,73221
+ ,317892)
+ ,dim=c(5
+ ,51)
+ ,dimnames=list(c('Costs'
+ ,'Trades'
+ ,'Orders'
+ ,'Dividends'
+ ,'Wealth')
+ ,1:51))
> y <- array(NA,dim=c(5,51),dimnames=list(c('Costs','Trades','Orders','Dividends','Wealth'),1:51))
> 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 = '5'
> #'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
Wealth Costs Trades Orders Dividends
1 194493 2284 33 41 76403
2 530670 3160 108 90 108094
3 518365 4150 150 136 134759
4 491303 7285 115 97 188873
5 527021 1134 162 63 146216
6 233773 4658 158 114 156608
7 405972 2384 97 77 61348
8 652925 3748 9 6 50350
9 446211 5371 66 47 87720
10 341340 1285 107 51 99489
11 387699 9327 101 85 87419
12 493408 5565 47 43 94355
13 146494 1528 38 32 60326
14 414462 3122 34 25 94670
15 364304 7561 87 77 82425
16 355178 2675 79 54 59017
17 357760 13253 947 251 90829
18 261216 880 74 15 80791
19 397144 2053 53 44 100423
20 374943 1424 94 73 131116
21 424898 4036 63 85 100269
22 202055 3045 58 49 27330
23 378525 5119 49 38 39039
24 310768 1431 34 35 106885
25 325738 554 11 9 79285
26 394510 1975 35 34 118881
27 247060 1765 20 20 77623
28 368078 1012 47 29 114768
29 236761 810 43 11 74015
30 312378 1280 117 52 69465
31 339836 666 171 13 117869
32 347385 1380 26 29 60982
33 426280 4677 75 66 90131
34 352850 876 59 33 138971
35 301881 814 18 15 39625
36 377516 514 15 15 102725
37 357312 5692 72 68 64239
38 458343 3642 86 100 90262
39 354228 540 14 13 103960
40 308636 2099 64 45 106611
41 386212 567 11 14 103345
42 393343 2001 52 36 95551
43 378509 2949 41 40 82903
44 452469 2253 99 68 63593
45 364839 6533 75 29 126910
46 358649 1889 45 43 37527
47 376641 3055 43 30 60247
48 429112 272 8 9 112995
49 330546 1414 198 22 70184
50 403560 2564 22 19 130140
51 317892 1383 11 9 73221
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Costs Trades Orders Dividends
2.774e+05 9.686e+00 -2.158e+02 3.975e+02 6.934e-01
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-208507 -30473 6574 29056 303911
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.774e+05 3.919e+04 7.078 6.97e-09 ***
Costs 9.686e+00 7.267e+00 1.333 0.1891
Trades -2.158e+02 1.584e+02 -1.363 0.1797
Orders 3.975e+02 5.997e+02 0.663 0.5108
Dividends 6.934e-01 4.002e-01 1.733 0.0899 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 87540 on 46 degrees of freedom
Multiple R-squared: 0.1543, Adjusted R-squared: 0.08078
F-statistic: 2.098 on 4 and 46 DF, p-value: 0.09633
> 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.9999936 1.285915e-05 6.429574e-06
[2,] 0.9999941 1.184760e-05 5.923798e-06
[3,] 0.9999870 2.590622e-05 1.295311e-05
[4,] 0.9999887 2.263995e-05 1.131998e-05
[5,] 0.9999861 2.788045e-05 1.394022e-05
[6,] 0.9999999 1.155128e-07 5.775641e-08
[7,] 0.9999999 2.617452e-07 1.308726e-07
[8,] 0.9999997 5.853527e-07 2.926763e-07
[9,] 0.9999991 1.888610e-06 9.443050e-07
[10,] 0.9999978 4.473623e-06 2.236811e-06
[11,] 0.9999976 4.710282e-06 2.355141e-06
[12,] 0.9999935 1.296275e-05 6.481376e-06
[13,] 0.9999872 2.560046e-05 1.280023e-05
[14,] 0.9999666 6.681985e-05 3.340992e-05
[15,] 0.9999957 8.505612e-06 4.252806e-06
[16,] 0.9999908 1.834216e-05 9.171080e-06
[17,] 0.9999905 1.901407e-05 9.507036e-06
[18,] 0.9999736 5.278914e-05 2.639457e-05
[19,] 0.9999292 1.415452e-04 7.077260e-05
[20,] 0.9999751 4.988423e-05 2.494212e-05
[21,] 0.9999315 1.370080e-04 6.850402e-05
[22,] 0.9999815 3.706038e-05 1.853019e-05
[23,] 0.9999802 3.963637e-05 1.981819e-05
[24,] 0.9999391 1.217861e-04 6.089303e-05
[25,] 0.9998316 3.368102e-04 1.684051e-04
[26,] 0.9996078 7.843303e-04 3.921652e-04
[27,] 0.9995450 9.100155e-04 4.550078e-04
[28,] 0.9990365 1.927069e-03 9.635343e-04
[29,] 0.9974091 5.181858e-03 2.590929e-03
[30,] 0.9943300 1.133999e-02 5.669996e-03
[31,] 0.9870434 2.591323e-02 1.295662e-02
[32,] 0.9726464 5.470729e-02 2.735364e-02
[33,] 0.9990585 1.883047e-03 9.415235e-04
[34,] 0.9962262 7.547631e-03 3.773815e-03
[35,] 0.9890358 2.192830e-02 1.096415e-02
[36,] 0.9683557 6.328865e-02 3.164432e-02
> postscript(file="/var/www/html/rcomp/tmp/1cbkt1291122107.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/html/rcomp/tmp/252je1291122107.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/html/rcomp/tmp/352je1291122107.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/html/rcomp/tmp/452je1291122107.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/html/rcomp/tmp/5xtih1291122107.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 = 51
Frequency = 1
1 2 3 4 5
-167137.21111 135288.44161 85684.82334 -1315.51455 147213.89510
6 7 8 9 10
-208507.06549 53312.45146 303910.75092 51567.98734 -14628.54870
11 12 13 14 15
-52605.60173 89775.65332 -192010.28309 38623.99703 -55273.07556
16 17 18 19 20
6573.59715 -6377.21079 -70677.06589 24219.65140 -15849.84815
21 22 23 24 25
18734.68149 -130706.49568 19985.65408 -61134.47304 -13161.60814
26 27 28 29 30
9633.81309 -104847.39875 -42.06220 -94854.69385 -20964.35349
31 32 33 34 35
6031.32893 8462.61546 31077.86392 -29734.85151 -12912.57032
36 37 38 39 40
21229.73543 -31211.31565 61935.38642 -2587.27417 -67048.49355
41 42 43 44 45
28516.86939 27264.25589 8054.11930 103530.34452 -59136.39514
46 47 48 49 50
29594.36912 25274.27931 68923.31899 24807.65445 8328.84088
51
-24832.96880
> postscript(file="/var/www/html/rcomp/tmp/6xtih1291122107.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 = 51
Frequency = 1
lag(myerror, k = 1) myerror
0 -167137.21111 NA
1 135288.44161 -167137.21111
2 85684.82334 135288.44161
3 -1315.51455 85684.82334
4 147213.89510 -1315.51455
5 -208507.06549 147213.89510
6 53312.45146 -208507.06549
7 303910.75092 53312.45146
8 51567.98734 303910.75092
9 -14628.54870 51567.98734
10 -52605.60173 -14628.54870
11 89775.65332 -52605.60173
12 -192010.28309 89775.65332
13 38623.99703 -192010.28309
14 -55273.07556 38623.99703
15 6573.59715 -55273.07556
16 -6377.21079 6573.59715
17 -70677.06589 -6377.21079
18 24219.65140 -70677.06589
19 -15849.84815 24219.65140
20 18734.68149 -15849.84815
21 -130706.49568 18734.68149
22 19985.65408 -130706.49568
23 -61134.47304 19985.65408
24 -13161.60814 -61134.47304
25 9633.81309 -13161.60814
26 -104847.39875 9633.81309
27 -42.06220 -104847.39875
28 -94854.69385 -42.06220
29 -20964.35349 -94854.69385
30 6031.32893 -20964.35349
31 8462.61546 6031.32893
32 31077.86392 8462.61546
33 -29734.85151 31077.86392
34 -12912.57032 -29734.85151
35 21229.73543 -12912.57032
36 -31211.31565 21229.73543
37 61935.38642 -31211.31565
38 -2587.27417 61935.38642
39 -67048.49355 -2587.27417
40 28516.86939 -67048.49355
41 27264.25589 28516.86939
42 8054.11930 27264.25589
43 103530.34452 8054.11930
44 -59136.39514 103530.34452
45 29594.36912 -59136.39514
46 25274.27931 29594.36912
47 68923.31899 25274.27931
48 24807.65445 68923.31899
49 8328.84088 24807.65445
50 -24832.96880 8328.84088
51 NA -24832.96880
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 135288.44161 -167137.21111
[2,] 85684.82334 135288.44161
[3,] -1315.51455 85684.82334
[4,] 147213.89510 -1315.51455
[5,] -208507.06549 147213.89510
[6,] 53312.45146 -208507.06549
[7,] 303910.75092 53312.45146
[8,] 51567.98734 303910.75092
[9,] -14628.54870 51567.98734
[10,] -52605.60173 -14628.54870
[11,] 89775.65332 -52605.60173
[12,] -192010.28309 89775.65332
[13,] 38623.99703 -192010.28309
[14,] -55273.07556 38623.99703
[15,] 6573.59715 -55273.07556
[16,] -6377.21079 6573.59715
[17,] -70677.06589 -6377.21079
[18,] 24219.65140 -70677.06589
[19,] -15849.84815 24219.65140
[20,] 18734.68149 -15849.84815
[21,] -130706.49568 18734.68149
[22,] 19985.65408 -130706.49568
[23,] -61134.47304 19985.65408
[24,] -13161.60814 -61134.47304
[25,] 9633.81309 -13161.60814
[26,] -104847.39875 9633.81309
[27,] -42.06220 -104847.39875
[28,] -94854.69385 -42.06220
[29,] -20964.35349 -94854.69385
[30,] 6031.32893 -20964.35349
[31,] 8462.61546 6031.32893
[32,] 31077.86392 8462.61546
[33,] -29734.85151 31077.86392
[34,] -12912.57032 -29734.85151
[35,] 21229.73543 -12912.57032
[36,] -31211.31565 21229.73543
[37,] 61935.38642 -31211.31565
[38,] -2587.27417 61935.38642
[39,] -67048.49355 -2587.27417
[40,] 28516.86939 -67048.49355
[41,] 27264.25589 28516.86939
[42,] 8054.11930 27264.25589
[43,] 103530.34452 8054.11930
[44,] -59136.39514 103530.34452
[45,] 29594.36912 -59136.39514
[46,] 25274.27931 29594.36912
[47,] 68923.31899 25274.27931
[48,] 24807.65445 68923.31899
[49,] 8328.84088 24807.65445
[50,] -24832.96880 8328.84088
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 135288.44161 -167137.21111
2 85684.82334 135288.44161
3 -1315.51455 85684.82334
4 147213.89510 -1315.51455
5 -208507.06549 147213.89510
6 53312.45146 -208507.06549
7 303910.75092 53312.45146
8 51567.98734 303910.75092
9 -14628.54870 51567.98734
10 -52605.60173 -14628.54870
11 89775.65332 -52605.60173
12 -192010.28309 89775.65332
13 38623.99703 -192010.28309
14 -55273.07556 38623.99703
15 6573.59715 -55273.07556
16 -6377.21079 6573.59715
17 -70677.06589 -6377.21079
18 24219.65140 -70677.06589
19 -15849.84815 24219.65140
20 18734.68149 -15849.84815
21 -130706.49568 18734.68149
22 19985.65408 -130706.49568
23 -61134.47304 19985.65408
24 -13161.60814 -61134.47304
25 9633.81309 -13161.60814
26 -104847.39875 9633.81309
27 -42.06220 -104847.39875
28 -94854.69385 -42.06220
29 -20964.35349 -94854.69385
30 6031.32893 -20964.35349
31 8462.61546 6031.32893
32 31077.86392 8462.61546
33 -29734.85151 31077.86392
34 -12912.57032 -29734.85151
35 21229.73543 -12912.57032
36 -31211.31565 21229.73543
37 61935.38642 -31211.31565
38 -2587.27417 61935.38642
39 -67048.49355 -2587.27417
40 28516.86939 -67048.49355
41 27264.25589 28516.86939
42 8054.11930 27264.25589
43 103530.34452 8054.11930
44 -59136.39514 103530.34452
45 29594.36912 -59136.39514
46 25274.27931 29594.36912
47 68923.31899 25274.27931
48 24807.65445 68923.31899
49 8328.84088 24807.65445
50 -24832.96880 8328.84088
> 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/78l0k1291122107.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/html/rcomp/tmp/88l0k1291122107.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/html/rcomp/tmp/91uh51291122107.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/html/rcomp/tmp/101uh51291122107.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/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/114cxt1291122107.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/127dwy1291122107.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/13webs1291122107.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/147nsv1291122107.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/15aor11291122107.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/16og7a1291122107.tab")
+ }
>
> try(system("convert tmp/1cbkt1291122107.ps tmp/1cbkt1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/252je1291122107.ps tmp/252je1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/352je1291122107.ps tmp/352je1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/452je1291122107.ps tmp/452je1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/5xtih1291122107.ps tmp/5xtih1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/6xtih1291122107.ps tmp/6xtih1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/78l0k1291122107.ps tmp/78l0k1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/88l0k1291122107.ps tmp/88l0k1291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/91uh51291122107.ps tmp/91uh51291122107.png",intern=TRUE))
character(0)
> try(system("convert tmp/101uh51291122107.ps tmp/101uh51291122107.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.443 1.624 8.413