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(2350.44,0,2440.25,0,2408.64,0,2472.81,0,2407.6,0,2454.62,0,2448.05,0,2497.84,0,2645.64,0,2756.76,0,2849.27,0,2921.44,0,2981.85,0,3080.58,0,3106.22,0,3119.31,0,3061.26,0,3097.31,0,3161.69,0,3257.16,0,3277.01,0,3295.32,0,3363.99,0,3494.17,0,3667.03,0,3813.06,0,3917.96,0,3895.51,0,3801.06,0,3570.12,0,3701.61,0,3862.27,0,3970.1,0,4138.52,0,4199.75,0,4290.89,0,4443.91,0,4502.64,0,4356.98,0,4591.27,0,4696.96,0,4621.4,0,4562.84,0,4202.52,0,4296.49,0,4435.23,0,4105.18,0,4116.68,0,3844.49,0,3720.98,0,3674.4,0,3857.62,0,3801.06,0,3504.37,0,3032.6,0,3047.03,0,2962.34,1,2197.82,1,2014.45,1,1862.83,1,1905.41,1,1810.99,1,1670.07,1,1864.44,1,2052.02,1,2029.6,1,2070.83,1,2293.41,1,2443.27,1,2513.17,1,2466.92,1),dim=c(2,71),dimnames=list(c('Y','X'),1:71))
> y <- array(NA,dim=c(2,71),dimnames=list(c('Y','X'),1:71))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = '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 t
1 2350.44 0 1
2 2440.25 0 2
3 2408.64 0 3
4 2472.81 0 4
5 2407.60 0 5
6 2454.62 0 6
7 2448.05 0 7
8 2497.84 0 8
9 2645.64 0 9
10 2756.76 0 10
11 2849.27 0 11
12 2921.44 0 12
13 2981.85 0 13
14 3080.58 0 14
15 3106.22 0 15
16 3119.31 0 16
17 3061.26 0 17
18 3097.31 0 18
19 3161.69 0 19
20 3257.16 0 20
21 3277.01 0 21
22 3295.32 0 22
23 3363.99 0 23
24 3494.17 0 24
25 3667.03 0 25
26 3813.06 0 26
27 3917.96 0 27
28 3895.51 0 28
29 3801.06 0 29
30 3570.12 0 30
31 3701.61 0 31
32 3862.27 0 32
33 3970.10 0 33
34 4138.52 0 34
35 4199.75 0 35
36 4290.89 0 36
37 4443.91 0 37
38 4502.64 0 38
39 4356.98 0 39
40 4591.27 0 40
41 4696.96 0 41
42 4621.40 0 42
43 4562.84 0 43
44 4202.52 0 44
45 4296.49 0 45
46 4435.23 0 46
47 4105.18 0 47
48 4116.68 0 48
49 3844.49 0 49
50 3720.98 0 50
51 3674.40 0 51
52 3857.62 0 52
53 3801.06 0 53
54 3504.37 0 54
55 3032.60 0 55
56 3047.03 0 56
57 2962.34 1 57
58 2197.82 1 58
59 2014.45 1 59
60 1862.83 1 60
61 1905.41 1 61
62 1810.99 1 62
63 1670.07 1 63
64 1864.44 1 64
65 2052.02 1 65
66 2029.60 1 66
67 2070.83 1 67
68 2293.41 1 68
69 2443.27 1 69
70 2513.17 1 70
71 2466.92 1 71
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X t
2616.26 -2504.75 31.76
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1347.516 -262.878 -5.033 254.071 1040.788
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2616.262 116.540 22.449 < 2e-16 ***
X -2504.749 178.160 -14.059 < 2e-16 ***
t 31.755 3.549 8.948 4.18e-13 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 433.3 on 68 degrees of freedom
Multiple R-squared: 0.7459, Adjusted R-squared: 0.7384
F-statistic: 99.81 on 2 and 68 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,] 1.510118e-03 3.020236e-03 0.9984898819
[2,] 1.443521e-04 2.887043e-04 0.9998556479
[3,] 1.391486e-05 2.782972e-05 0.9999860851
[4,] 3.015549e-05 6.031097e-05 0.9999698445
[5,] 4.401884e-05 8.803768e-05 0.9999559812
[6,] 3.729577e-05 7.459154e-05 0.9999627042
[7,] 1.992880e-05 3.985760e-05 0.9999800712
[8,] 7.697870e-06 1.539574e-05 0.9999923021
[9,] 3.439626e-06 6.879252e-06 0.9999965604
[10,] 9.488566e-07 1.897713e-06 0.9999990511
[11,] 2.159010e-07 4.318020e-07 0.9999997841
[12,] 8.912461e-08 1.782492e-07 0.9999999109
[13,] 3.759412e-08 7.518823e-08 0.9999999624
[14,] 1.276089e-08 2.552179e-08 0.9999999872
[15,] 3.407503e-09 6.815007e-09 0.9999999966
[16,] 1.062223e-09 2.124447e-09 0.9999999989
[17,] 4.324286e-10 8.648572e-10 0.9999999996
[18,] 1.571623e-10 3.143246e-10 0.9999999998
[19,] 5.701021e-11 1.140204e-10 0.9999999999
[20,] 7.412656e-11 1.482531e-10 0.9999999999
[21,] 2.823949e-10 5.647899e-10 0.9999999997
[22,] 9.271237e-10 1.854247e-09 0.9999999991
[23,] 5.530417e-10 1.106083e-09 0.9999999994
[24,] 2.181975e-10 4.363949e-10 0.9999999998
[25,] 4.767361e-09 9.534722e-09 0.9999999952
[26,] 1.280472e-08 2.560944e-08 0.9999999872
[27,] 1.125880e-08 2.251761e-08 0.9999999887
[28,] 8.043159e-09 1.608632e-08 0.9999999920
[29,] 5.843554e-09 1.168711e-08 0.9999999942
[30,] 3.814046e-09 7.628091e-09 0.9999999962
[31,] 2.483848e-09 4.967695e-09 0.9999999975
[32,] 2.791948e-09 5.583896e-09 0.9999999972
[33,] 2.334281e-09 4.668562e-09 0.9999999977
[34,] 7.816796e-10 1.563359e-09 0.9999999992
[35,] 4.996574e-10 9.993148e-10 0.9999999995
[36,] 6.526656e-10 1.305331e-09 0.9999999993
[37,] 4.341599e-10 8.683197e-10 0.9999999996
[38,] 4.316598e-10 8.633196e-10 0.9999999996
[39,] 2.317208e-08 4.634416e-08 0.9999999768
[40,] 1.999474e-07 3.998947e-07 0.9999998001
[41,] 1.353581e-06 2.707163e-06 0.9999986464
[42,] 6.288497e-05 1.257699e-04 0.9999371150
[43,] 9.876255e-04 1.975251e-03 0.9990123745
[44,] 1.431839e-02 2.863678e-02 0.9856816086
[45,] 7.173376e-02 1.434675e-01 0.9282662362
[46,] 1.698860e-01 3.397719e-01 0.8301140408
[47,] 2.860755e-01 5.721511e-01 0.7139244658
[48,] 4.647347e-01 9.294695e-01 0.5352652646
[49,] 6.329270e-01 7.341460e-01 0.3670730119
[50,] 7.639358e-01 4.721285e-01 0.2360642377
[51,] 8.107534e-01 3.784931e-01 0.1892465545
[52,] 9.959024e-01 8.195150e-03 0.0040975752
[53,] 9.991749e-01 1.650189e-03 0.0008250945
[54,] 9.997292e-01 5.416310e-04 0.0002708155
[55,] 9.996958e-01 6.083520e-04 0.0003041760
[56,] 9.998811e-01 2.378295e-04 0.0001189147
[57,] 9.997959e-01 4.081102e-04 0.0002040551
[58,] 9.993689e-01 1.262191e-03 0.0006310955
[59,] 9.965488e-01 6.902354e-03 0.0034511772
[60,] 9.870672e-01 2.586563e-02 0.0129328161
> postscript(file="/var/www/html/rcomp/tmp/1a7031260884501.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/256dp1260884501.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/3165d1260884501.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/43z9w1260884501.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/5jsvp1260884501.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 = 71
Frequency = 1
1 2 3 4 5 6
-297.576715 -239.521796 -302.886876 -270.471957 -367.437037 -352.172118
7 8 9 10 11 12
-390.497198 -372.462279 -256.417359 -177.052440 -116.297520 -75.882601
13 14 15 16 17 18
-47.227681 19.747238 13.632158 -5.032923 -94.838003 -90.543084
19 20 21 22 23 24
-57.918164 5.796756 -6.108325 -19.553405 17.361514 115.786434
25 26 27 28 29 30
256.891353 371.166273 444.311192 390.106112 263.901031 1.205951
31 32 33 34 35 36
100.940870 229.845790 305.920709 442.585629 472.060548 531.445468
37 38 39 40 41 42
652.710387 679.685307 502.270226 704.805146 778.740065 671.424985
43 44 45 46 47 48
581.109904 189.034824 251.249743 358.234663 -3.570418 -23.825498
49 50 51 52 53 54
-327.770578 -483.035659 -561.370739 -409.905820 -498.220900 -826.665981
55 56 57 58 59 60
-1330.191061 -1347.516142 1040.787563 244.512483 29.387402 -153.987678
61 62 63 64 65 66
-143.162759 -269.337839 -442.012920 -279.398000 -123.573080 -177.748161
67 68 69 70 71
-168.273241 22.551678 140.656598 178.801517 100.796437
> postscript(file="/var/www/html/rcomp/tmp/6l3u01260884501.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 = 71
Frequency = 1
lag(myerror, k = 1) myerror
0 -297.576715 NA
1 -239.521796 -297.576715
2 -302.886876 -239.521796
3 -270.471957 -302.886876
4 -367.437037 -270.471957
5 -352.172118 -367.437037
6 -390.497198 -352.172118
7 -372.462279 -390.497198
8 -256.417359 -372.462279
9 -177.052440 -256.417359
10 -116.297520 -177.052440
11 -75.882601 -116.297520
12 -47.227681 -75.882601
13 19.747238 -47.227681
14 13.632158 19.747238
15 -5.032923 13.632158
16 -94.838003 -5.032923
17 -90.543084 -94.838003
18 -57.918164 -90.543084
19 5.796756 -57.918164
20 -6.108325 5.796756
21 -19.553405 -6.108325
22 17.361514 -19.553405
23 115.786434 17.361514
24 256.891353 115.786434
25 371.166273 256.891353
26 444.311192 371.166273
27 390.106112 444.311192
28 263.901031 390.106112
29 1.205951 263.901031
30 100.940870 1.205951
31 229.845790 100.940870
32 305.920709 229.845790
33 442.585629 305.920709
34 472.060548 442.585629
35 531.445468 472.060548
36 652.710387 531.445468
37 679.685307 652.710387
38 502.270226 679.685307
39 704.805146 502.270226
40 778.740065 704.805146
41 671.424985 778.740065
42 581.109904 671.424985
43 189.034824 581.109904
44 251.249743 189.034824
45 358.234663 251.249743
46 -3.570418 358.234663
47 -23.825498 -3.570418
48 -327.770578 -23.825498
49 -483.035659 -327.770578
50 -561.370739 -483.035659
51 -409.905820 -561.370739
52 -498.220900 -409.905820
53 -826.665981 -498.220900
54 -1330.191061 -826.665981
55 -1347.516142 -1330.191061
56 1040.787563 -1347.516142
57 244.512483 1040.787563
58 29.387402 244.512483
59 -153.987678 29.387402
60 -143.162759 -153.987678
61 -269.337839 -143.162759
62 -442.012920 -269.337839
63 -279.398000 -442.012920
64 -123.573080 -279.398000
65 -177.748161 -123.573080
66 -168.273241 -177.748161
67 22.551678 -168.273241
68 140.656598 22.551678
69 178.801517 140.656598
70 100.796437 178.801517
71 NA 100.796437
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -239.521796 -297.576715
[2,] -302.886876 -239.521796
[3,] -270.471957 -302.886876
[4,] -367.437037 -270.471957
[5,] -352.172118 -367.437037
[6,] -390.497198 -352.172118
[7,] -372.462279 -390.497198
[8,] -256.417359 -372.462279
[9,] -177.052440 -256.417359
[10,] -116.297520 -177.052440
[11,] -75.882601 -116.297520
[12,] -47.227681 -75.882601
[13,] 19.747238 -47.227681
[14,] 13.632158 19.747238
[15,] -5.032923 13.632158
[16,] -94.838003 -5.032923
[17,] -90.543084 -94.838003
[18,] -57.918164 -90.543084
[19,] 5.796756 -57.918164
[20,] -6.108325 5.796756
[21,] -19.553405 -6.108325
[22,] 17.361514 -19.553405
[23,] 115.786434 17.361514
[24,] 256.891353 115.786434
[25,] 371.166273 256.891353
[26,] 444.311192 371.166273
[27,] 390.106112 444.311192
[28,] 263.901031 390.106112
[29,] 1.205951 263.901031
[30,] 100.940870 1.205951
[31,] 229.845790 100.940870
[32,] 305.920709 229.845790
[33,] 442.585629 305.920709
[34,] 472.060548 442.585629
[35,] 531.445468 472.060548
[36,] 652.710387 531.445468
[37,] 679.685307 652.710387
[38,] 502.270226 679.685307
[39,] 704.805146 502.270226
[40,] 778.740065 704.805146
[41,] 671.424985 778.740065
[42,] 581.109904 671.424985
[43,] 189.034824 581.109904
[44,] 251.249743 189.034824
[45,] 358.234663 251.249743
[46,] -3.570418 358.234663
[47,] -23.825498 -3.570418
[48,] -327.770578 -23.825498
[49,] -483.035659 -327.770578
[50,] -561.370739 -483.035659
[51,] -409.905820 -561.370739
[52,] -498.220900 -409.905820
[53,] -826.665981 -498.220900
[54,] -1330.191061 -826.665981
[55,] -1347.516142 -1330.191061
[56,] 1040.787563 -1347.516142
[57,] 244.512483 1040.787563
[58,] 29.387402 244.512483
[59,] -153.987678 29.387402
[60,] -143.162759 -153.987678
[61,] -269.337839 -143.162759
[62,] -442.012920 -269.337839
[63,] -279.398000 -442.012920
[64,] -123.573080 -279.398000
[65,] -177.748161 -123.573080
[66,] -168.273241 -177.748161
[67,] 22.551678 -168.273241
[68,] 140.656598 22.551678
[69,] 178.801517 140.656598
[70,] 100.796437 178.801517
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -239.521796 -297.576715
2 -302.886876 -239.521796
3 -270.471957 -302.886876
4 -367.437037 -270.471957
5 -352.172118 -367.437037
6 -390.497198 -352.172118
7 -372.462279 -390.497198
8 -256.417359 -372.462279
9 -177.052440 -256.417359
10 -116.297520 -177.052440
11 -75.882601 -116.297520
12 -47.227681 -75.882601
13 19.747238 -47.227681
14 13.632158 19.747238
15 -5.032923 13.632158
16 -94.838003 -5.032923
17 -90.543084 -94.838003
18 -57.918164 -90.543084
19 5.796756 -57.918164
20 -6.108325 5.796756
21 -19.553405 -6.108325
22 17.361514 -19.553405
23 115.786434 17.361514
24 256.891353 115.786434
25 371.166273 256.891353
26 444.311192 371.166273
27 390.106112 444.311192
28 263.901031 390.106112
29 1.205951 263.901031
30 100.940870 1.205951
31 229.845790 100.940870
32 305.920709 229.845790
33 442.585629 305.920709
34 472.060548 442.585629
35 531.445468 472.060548
36 652.710387 531.445468
37 679.685307 652.710387
38 502.270226 679.685307
39 704.805146 502.270226
40 778.740065 704.805146
41 671.424985 778.740065
42 581.109904 671.424985
43 189.034824 581.109904
44 251.249743 189.034824
45 358.234663 251.249743
46 -3.570418 358.234663
47 -23.825498 -3.570418
48 -327.770578 -23.825498
49 -483.035659 -327.770578
50 -561.370739 -483.035659
51 -409.905820 -561.370739
52 -498.220900 -409.905820
53 -826.665981 -498.220900
54 -1330.191061 -826.665981
55 -1347.516142 -1330.191061
56 1040.787563 -1347.516142
57 244.512483 1040.787563
58 29.387402 244.512483
59 -153.987678 29.387402
60 -143.162759 -153.987678
61 -269.337839 -143.162759
62 -442.012920 -269.337839
63 -279.398000 -442.012920
64 -123.573080 -279.398000
65 -177.748161 -123.573080
66 -168.273241 -177.748161
67 22.551678 -168.273241
68 140.656598 22.551678
69 178.801517 140.656598
70 100.796437 178.801517
> 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/728cw1260884501.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/80y9s1260884501.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/9ysoz1260884501.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/10bvoz1260884501.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/11dcgl1260884501.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/127y2z1260884501.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/131dcb1260884501.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/14lsxl1260884501.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/15ya0x1260884501.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/16seyx1260884501.tab")
+ }
>
> try(system("convert tmp/1a7031260884501.ps tmp/1a7031260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/256dp1260884501.ps tmp/256dp1260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/3165d1260884501.ps tmp/3165d1260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/43z9w1260884501.ps tmp/43z9w1260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/5jsvp1260884501.ps tmp/5jsvp1260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/6l3u01260884501.ps tmp/6l3u01260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/728cw1260884501.ps tmp/728cw1260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/80y9s1260884501.ps tmp/80y9s1260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/9ysoz1260884501.ps tmp/9ysoz1260884501.png",intern=TRUE))
character(0)
> try(system("convert tmp/10bvoz1260884501.ps tmp/10bvoz1260884501.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.497 1.588 4.471