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(8.4,420,8.4,418,8.4,410,8.6,418,8.9,426,8.8,428,8.3,430,7.5,424,7.2,423,7.4,427,8.8,441,9.3,449,9.3,452,8.7,462,8.2,455,8.3,461,8.5,461,8.6,463,8.5,462,8.2,456,8.1,455,7.9,456,8.6,472,8.7,472,8.7,471,8.5,465,8.4,459,8.5,465,8.7,468,8.7,467,8.6,463,8.5,460,8.3,462,8.00,461,8.2,476,8.1,476,8.1,471,8.00,453,7.9,443,7.9,442,8.00,444,8.00,438,7.9,427,8.00,424,7.7,416,7.2,406,7.5,431,7.3,434,7.00,418,7.00,412,7.00,404,7.2,409,7.3,412,7.1,406,6.8,398,6.4,397,6.1,385,6.5,390,7.7,413,7.9,413,7.5,401,6.9,397,6.6,397,6.9,409,7.7,419,8.00,424,8.00,428,7.7,430,7.3,424,7.4,433,8.1,456,8.3,459,8.2,446),dim=c(2,73),dimnames=list(c('wgb','nwwz'),1:73))
> y <- array(NA,dim=c(2,73),dimnames=list(c('wgb','nwwz'),1:73))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = 'No Linear Trend'
> par2 = 'Do not include Seasonal Dummies'
> par1 = '1'
> #'GNU S' R Code compiled by R2WASP v. 1.0.44 ()
> #Author: Prof. Dr. P. Wessa
> #To cite this work: AUTHOR(S), (YEAR), YOUR SOFTWARE TITLE (vNUMBER) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_YOURPAGE.wasp/
> #Source of accompanying publication: Office for Research, Development, and Education
> #Technical description: Write here your technical program description (don't use hard returns!)
> library(lattice)
> library(lmtest)
Loading required package: zoo
Attaching package: 'zoo'
The following object(s) are masked from package:base :
as.Date.numeric
> n25 <- 25 #minimum number of obs. for Goldfeld-Quandt test
> par1 <- as.numeric(par1)
> x <- t(y)
> k <- length(x[1,])
> n <- length(x[,1])
> x1 <- cbind(x[,par1], x[,1:k!=par1])
> mycolnames <- c(colnames(x)[par1], colnames(x)[1:k!=par1])
> colnames(x1) <- mycolnames #colnames(x)[par1]
> x <- x1
> if (par3 == 'First Differences'){
+ x2 <- array(0, dim=c(n-1,k), dimnames=list(1:(n-1), paste('(1-B)',colnames(x),sep='')))
+ for (i in 1:n-1) {
+ for (j in 1:k) {
+ x2[i,j] <- x[i+1,j] - x[i,j]
+ }
+ }
+ x <- x2
+ }
> if (par2 == 'Include Monthly Dummies'){
+ x2 <- array(0, dim=c(n,11), dimnames=list(1:n, paste('M', seq(1:11), sep ='')))
+ for (i in 1:11){
+ x2[seq(i,n,12),i] <- 1
+ }
+ x <- cbind(x, x2)
+ }
> if (par2 == 'Include Quarterly Dummies'){
+ x2 <- array(0, dim=c(n,3), dimnames=list(1:n, paste('Q', seq(1:3), sep ='')))
+ for (i in 1:3){
+ x2[seq(i,n,4),i] <- 1
+ }
+ x <- cbind(x, x2)
+ }
> k <- length(x[1,])
> if (par3 == 'Linear Trend'){
+ x <- cbind(x, c(1:n))
+ colnames(x)[k+1] <- 't'
+ }
> x
wgb nwwz
1 8.4 420
2 8.4 418
3 8.4 410
4 8.6 418
5 8.9 426
6 8.8 428
7 8.3 430
8 7.5 424
9 7.2 423
10 7.4 427
11 8.8 441
12 9.3 449
13 9.3 452
14 8.7 462
15 8.2 455
16 8.3 461
17 8.5 461
18 8.6 463
19 8.5 462
20 8.2 456
21 8.1 455
22 7.9 456
23 8.6 472
24 8.7 472
25 8.7 471
26 8.5 465
27 8.4 459
28 8.5 465
29 8.7 468
30 8.7 467
31 8.6 463
32 8.5 460
33 8.3 462
34 8.0 461
35 8.2 476
36 8.1 476
37 8.1 471
38 8.0 453
39 7.9 443
40 7.9 442
41 8.0 444
42 8.0 438
43 7.9 427
44 8.0 424
45 7.7 416
46 7.2 406
47 7.5 431
48 7.3 434
49 7.0 418
50 7.0 412
51 7.0 404
52 7.2 409
53 7.3 412
54 7.1 406
55 6.8 398
56 6.4 397
57 6.1 385
58 6.5 390
59 7.7 413
60 7.9 413
61 7.5 401
62 6.9 397
63 6.6 397
64 6.9 409
65 7.7 419
66 8.0 424
67 8.0 428
68 7.7 430
69 7.3 424
70 7.4 433
71 8.1 456
72 8.3 459
73 8.2 446
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) nwwz
-1.18542 0.02090
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.76008 -0.32136 -0.09892 0.16223 1.18313
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.185421 0.972694 -1.219 0.227
nwwz 0.020897 0.002226 9.386 4.56e-14 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.4675 on 71 degrees of freedom
Multiple R-squared: 0.5537, Adjusted R-squared: 0.5475
F-statistic: 88.1 on 1 and 71 DF, p-value: 4.556e-14
> 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.09681849 0.1936369721 9.031815e-01
[2,] 0.04277190 0.0855437940 9.572281e-01
[3,] 0.17130989 0.3426197756 8.286901e-01
[4,] 0.83498654 0.3300269197 1.650135e-01
[5,] 0.98550784 0.0289843126 1.449216e-02
[6,] 0.99309520 0.0138096075 6.904804e-03
[7,] 0.99659805 0.0068038946 3.401947e-03
[8,] 0.99943170 0.0011366011 5.683006e-04
[9,] 0.99990969 0.0001806217 9.031084e-05
[10,] 0.99992711 0.0001457750 7.288751e-05
[11,] 0.99994778 0.0001044331 5.221655e-05
[12,] 0.99993918 0.0001216451 6.082253e-05
[13,] 0.99989119 0.0002176140 1.088070e-04
[14,] 0.99980650 0.0003870004 1.935002e-04
[15,] 0.99965826 0.0006834767 3.417383e-04
[16,] 0.99950798 0.0009840480 4.920240e-04
[17,] 0.99936418 0.0012716479 6.358239e-04
[18,] 0.99950780 0.0009843916 4.921958e-04
[19,] 0.99908390 0.0018321997 9.160999e-04
[20,] 0.99843318 0.0031336392 1.566820e-03
[21,] 0.99743672 0.0051265655 2.563283e-03
[22,] 0.99574779 0.0085044235 4.252212e-03
[23,] 0.99333347 0.0133330624 6.666531e-03
[24,] 0.98955452 0.0208909581 1.044548e-02
[25,] 0.98574705 0.0285059091 1.425295e-02
[26,] 0.98162488 0.0367502393 1.837512e-02
[27,] 0.97641955 0.0471609083 2.358045e-02
[28,] 0.96957417 0.0608516533 3.042583e-02
[29,] 0.95848021 0.0830395896 4.151979e-02
[30,] 0.95759148 0.0848170308 4.240852e-02
[31,] 0.95695785 0.0860842987 4.304215e-02
[32,] 0.96882257 0.0623548645 3.117743e-02
[33,] 0.97714117 0.0457176644 2.285883e-02
[34,] 0.97448397 0.0510320633 2.551603e-02
[35,] 0.96934269 0.0613146260 3.065731e-02
[36,] 0.96230770 0.0753846093 3.769230e-02
[37,] 0.95018326 0.0996334820 4.981674e-02
[38,] 0.93399400 0.1320120015 6.600600e-02
[39,] 0.92401115 0.1519777001 7.598885e-02
[40,] 0.92912100 0.1417580010 7.087900e-02
[41,] 0.93213848 0.1357230493 6.786152e-02
[42,] 0.93624855 0.1275029010 6.375145e-02
[43,] 0.93233743 0.1353251321 6.766257e-02
[44,] 0.95908959 0.0818208248 4.091041e-02
[45,] 0.97276569 0.0544686162 2.723431e-02
[46,] 0.97336001 0.0532799765 2.663999e-02
[47,] 0.96545885 0.0690823072 3.454115e-02
[48,] 0.95130186 0.0973962883 4.869814e-02
[49,] 0.93077706 0.1384458788 6.922294e-02
[50,] 0.90510161 0.1897967864 9.489839e-02
[51,] 0.87796237 0.2440752632 1.220376e-01
[52,] 0.90680910 0.1863818056 9.319090e-02
[53,] 0.94273733 0.1145253364 5.726267e-02
[54,] 0.94155904 0.1168819199 5.844096e-02
[55,] 0.92707379 0.1458524248 7.292621e-02
[56,] 0.95151016 0.0969796812 4.848984e-02
[57,] 0.96492475 0.0701505099 3.507525e-02
[58,] 0.93772818 0.1245436314 6.227182e-02
[59,] 0.92601373 0.1479725353 7.398627e-02
[60,] 0.93946662 0.1210667679 6.053338e-02
[61,] 0.89232014 0.2153597140 1.076799e-01
[62,] 0.90725429 0.1854914209 9.274571e-02
[63,] 0.95569868 0.0886026478 4.430132e-02
[64,] 0.92197154 0.1560569278 7.802846e-02
> postscript(file="/var/www/html/rcomp/tmp/165801258620892.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/2f8yj1258620892.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/35je01258620892.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/4kdle1258620892.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/5a88r1258620892.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 = 73
Frequency = 1
1 2 3 4 5 6
0.808515292 0.850310079 1.017489226 1.050310079 1.183130932 1.041336146
7 8 9 10 11 12
0.499541359 -0.175074281 -0.454176888 -0.337766461 0.769670032 1.102490885
13 14 15 16 17 18
1.039798705 0.230824771 -0.122893475 -0.148277835 0.051722165 0.109927378
19 20 21 22 23 24
0.030824771 -0.143790869 -0.222893475 -0.443790869 -0.078149162 0.021850838
25 26 27 28 29 30
0.042748231 -0.031867409 -0.006483049 -0.031867409 0.105440411 0.126337804
31 32 33 34 35 36
0.109927378 0.072619558 -0.169175229 -0.448277835 -0.561738736 -0.661738736
37 38 39 40 41 42
-0.557251769 -0.281098688 -0.172124755 -0.151227362 -0.093022148 0.032362212
43 44 45 46 47 48
0.162233539 0.324925719 0.192104866 -0.098921201 -0.321356035 -0.584048215
49 50 51 52 53 54
-0.549689921 -0.424305561 -0.257126414 -0.161613381 -0.124305561 -0.198921201
55 56 57 58 59 60
-0.331742054 -0.710844660 -0.760075940 -0.464562907 0.254797046 0.454797046
61 62 63 64 65 66
0.305565766 -0.210844660 -0.510844660 -0.461613381 0.129412686 0.324925719
67 68 69 70 71 72
0.241336146 -0.100458641 -0.375074281 -0.463150821 -0.243790869 -0.106483049
73
0.065183065
> postscript(file="/var/www/html/rcomp/tmp/6lvmk1258620892.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 = 73
Frequency = 1
lag(myerror, k = 1) myerror
0 0.808515292 NA
1 0.850310079 0.808515292
2 1.017489226 0.850310079
3 1.050310079 1.017489226
4 1.183130932 1.050310079
5 1.041336146 1.183130932
6 0.499541359 1.041336146
7 -0.175074281 0.499541359
8 -0.454176888 -0.175074281
9 -0.337766461 -0.454176888
10 0.769670032 -0.337766461
11 1.102490885 0.769670032
12 1.039798705 1.102490885
13 0.230824771 1.039798705
14 -0.122893475 0.230824771
15 -0.148277835 -0.122893475
16 0.051722165 -0.148277835
17 0.109927378 0.051722165
18 0.030824771 0.109927378
19 -0.143790869 0.030824771
20 -0.222893475 -0.143790869
21 -0.443790869 -0.222893475
22 -0.078149162 -0.443790869
23 0.021850838 -0.078149162
24 0.042748231 0.021850838
25 -0.031867409 0.042748231
26 -0.006483049 -0.031867409
27 -0.031867409 -0.006483049
28 0.105440411 -0.031867409
29 0.126337804 0.105440411
30 0.109927378 0.126337804
31 0.072619558 0.109927378
32 -0.169175229 0.072619558
33 -0.448277835 -0.169175229
34 -0.561738736 -0.448277835
35 -0.661738736 -0.561738736
36 -0.557251769 -0.661738736
37 -0.281098688 -0.557251769
38 -0.172124755 -0.281098688
39 -0.151227362 -0.172124755
40 -0.093022148 -0.151227362
41 0.032362212 -0.093022148
42 0.162233539 0.032362212
43 0.324925719 0.162233539
44 0.192104866 0.324925719
45 -0.098921201 0.192104866
46 -0.321356035 -0.098921201
47 -0.584048215 -0.321356035
48 -0.549689921 -0.584048215
49 -0.424305561 -0.549689921
50 -0.257126414 -0.424305561
51 -0.161613381 -0.257126414
52 -0.124305561 -0.161613381
53 -0.198921201 -0.124305561
54 -0.331742054 -0.198921201
55 -0.710844660 -0.331742054
56 -0.760075940 -0.710844660
57 -0.464562907 -0.760075940
58 0.254797046 -0.464562907
59 0.454797046 0.254797046
60 0.305565766 0.454797046
61 -0.210844660 0.305565766
62 -0.510844660 -0.210844660
63 -0.461613381 -0.510844660
64 0.129412686 -0.461613381
65 0.324925719 0.129412686
66 0.241336146 0.324925719
67 -0.100458641 0.241336146
68 -0.375074281 -0.100458641
69 -0.463150821 -0.375074281
70 -0.243790869 -0.463150821
71 -0.106483049 -0.243790869
72 0.065183065 -0.106483049
73 NA 0.065183065
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 0.850310079 0.808515292
[2,] 1.017489226 0.850310079
[3,] 1.050310079 1.017489226
[4,] 1.183130932 1.050310079
[5,] 1.041336146 1.183130932
[6,] 0.499541359 1.041336146
[7,] -0.175074281 0.499541359
[8,] -0.454176888 -0.175074281
[9,] -0.337766461 -0.454176888
[10,] 0.769670032 -0.337766461
[11,] 1.102490885 0.769670032
[12,] 1.039798705 1.102490885
[13,] 0.230824771 1.039798705
[14,] -0.122893475 0.230824771
[15,] -0.148277835 -0.122893475
[16,] 0.051722165 -0.148277835
[17,] 0.109927378 0.051722165
[18,] 0.030824771 0.109927378
[19,] -0.143790869 0.030824771
[20,] -0.222893475 -0.143790869
[21,] -0.443790869 -0.222893475
[22,] -0.078149162 -0.443790869
[23,] 0.021850838 -0.078149162
[24,] 0.042748231 0.021850838
[25,] -0.031867409 0.042748231
[26,] -0.006483049 -0.031867409
[27,] -0.031867409 -0.006483049
[28,] 0.105440411 -0.031867409
[29,] 0.126337804 0.105440411
[30,] 0.109927378 0.126337804
[31,] 0.072619558 0.109927378
[32,] -0.169175229 0.072619558
[33,] -0.448277835 -0.169175229
[34,] -0.561738736 -0.448277835
[35,] -0.661738736 -0.561738736
[36,] -0.557251769 -0.661738736
[37,] -0.281098688 -0.557251769
[38,] -0.172124755 -0.281098688
[39,] -0.151227362 -0.172124755
[40,] -0.093022148 -0.151227362
[41,] 0.032362212 -0.093022148
[42,] 0.162233539 0.032362212
[43,] 0.324925719 0.162233539
[44,] 0.192104866 0.324925719
[45,] -0.098921201 0.192104866
[46,] -0.321356035 -0.098921201
[47,] -0.584048215 -0.321356035
[48,] -0.549689921 -0.584048215
[49,] -0.424305561 -0.549689921
[50,] -0.257126414 -0.424305561
[51,] -0.161613381 -0.257126414
[52,] -0.124305561 -0.161613381
[53,] -0.198921201 -0.124305561
[54,] -0.331742054 -0.198921201
[55,] -0.710844660 -0.331742054
[56,] -0.760075940 -0.710844660
[57,] -0.464562907 -0.760075940
[58,] 0.254797046 -0.464562907
[59,] 0.454797046 0.254797046
[60,] 0.305565766 0.454797046
[61,] -0.210844660 0.305565766
[62,] -0.510844660 -0.210844660
[63,] -0.461613381 -0.510844660
[64,] 0.129412686 -0.461613381
[65,] 0.324925719 0.129412686
[66,] 0.241336146 0.324925719
[67,] -0.100458641 0.241336146
[68,] -0.375074281 -0.100458641
[69,] -0.463150821 -0.375074281
[70,] -0.243790869 -0.463150821
[71,] -0.106483049 -0.243790869
[72,] 0.065183065 -0.106483049
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 0.850310079 0.808515292
2 1.017489226 0.850310079
3 1.050310079 1.017489226
4 1.183130932 1.050310079
5 1.041336146 1.183130932
6 0.499541359 1.041336146
7 -0.175074281 0.499541359
8 -0.454176888 -0.175074281
9 -0.337766461 -0.454176888
10 0.769670032 -0.337766461
11 1.102490885 0.769670032
12 1.039798705 1.102490885
13 0.230824771 1.039798705
14 -0.122893475 0.230824771
15 -0.148277835 -0.122893475
16 0.051722165 -0.148277835
17 0.109927378 0.051722165
18 0.030824771 0.109927378
19 -0.143790869 0.030824771
20 -0.222893475 -0.143790869
21 -0.443790869 -0.222893475
22 -0.078149162 -0.443790869
23 0.021850838 -0.078149162
24 0.042748231 0.021850838
25 -0.031867409 0.042748231
26 -0.006483049 -0.031867409
27 -0.031867409 -0.006483049
28 0.105440411 -0.031867409
29 0.126337804 0.105440411
30 0.109927378 0.126337804
31 0.072619558 0.109927378
32 -0.169175229 0.072619558
33 -0.448277835 -0.169175229
34 -0.561738736 -0.448277835
35 -0.661738736 -0.561738736
36 -0.557251769 -0.661738736
37 -0.281098688 -0.557251769
38 -0.172124755 -0.281098688
39 -0.151227362 -0.172124755
40 -0.093022148 -0.151227362
41 0.032362212 -0.093022148
42 0.162233539 0.032362212
43 0.324925719 0.162233539
44 0.192104866 0.324925719
45 -0.098921201 0.192104866
46 -0.321356035 -0.098921201
47 -0.584048215 -0.321356035
48 -0.549689921 -0.584048215
49 -0.424305561 -0.549689921
50 -0.257126414 -0.424305561
51 -0.161613381 -0.257126414
52 -0.124305561 -0.161613381
53 -0.198921201 -0.124305561
54 -0.331742054 -0.198921201
55 -0.710844660 -0.331742054
56 -0.760075940 -0.710844660
57 -0.464562907 -0.760075940
58 0.254797046 -0.464562907
59 0.454797046 0.254797046
60 0.305565766 0.454797046
61 -0.210844660 0.305565766
62 -0.510844660 -0.210844660
63 -0.461613381 -0.510844660
64 0.129412686 -0.461613381
65 0.324925719 0.129412686
66 0.241336146 0.324925719
67 -0.100458641 0.241336146
68 -0.375074281 -0.100458641
69 -0.463150821 -0.375074281
70 -0.243790869 -0.463150821
71 -0.106483049 -0.243790869
72 0.065183065 -0.106483049
> 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/7z1pe1258620892.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/8zcqb1258620892.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/9d7x81258620892.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/10dtc61258620892.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/11loau1258620892.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/12faut1258620892.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/13u7nc1258620892.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/14gzpa1258620892.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/15cyag1258620892.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/16uhtb1258620892.tab")
+ }
>
> system("convert tmp/165801258620892.ps tmp/165801258620892.png")
> system("convert tmp/2f8yj1258620892.ps tmp/2f8yj1258620892.png")
> system("convert tmp/35je01258620892.ps tmp/35je01258620892.png")
> system("convert tmp/4kdle1258620892.ps tmp/4kdle1258620892.png")
> system("convert tmp/5a88r1258620892.ps tmp/5a88r1258620892.png")
> system("convert tmp/6lvmk1258620892.ps tmp/6lvmk1258620892.png")
> system("convert tmp/7z1pe1258620892.ps tmp/7z1pe1258620892.png")
> system("convert tmp/8zcqb1258620892.ps tmp/8zcqb1258620892.png")
> system("convert tmp/9d7x81258620892.ps tmp/9d7x81258620892.png")
> system("convert tmp/10dtc61258620892.ps tmp/10dtc61258620892.png")
>
>
> proc.time()
user system elapsed
2.578 1.530 4.646