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(363,14.3,364,14.2,363,15.9,358,15.3,357,15.5,357,15.1,380,15,378,12.1,376,15.8,380,16.9,379,15.1,384,13.7,392,14.8,394,14.7,392,16,396,15.4,392,15,396,15.5,419,15.1,421,11.7,420,16.3,418,16.7,410,15,418,14.9,426,14.6,428,15.3,430,17.9,424,16.4,423,15.4,427,17.9,441,15.9,449,13.9,452,17.8,462,17.9,455,17.4,461,16.7,461,16,463,16.6,462,19.1,456,17.8,455,17.2,456,18.6,472,16.3,472,15.1,471,19.2,465,17.7,459,19.1,465,18,468,17.5,467,17.8,463,21.1,460,17.2,462,19.4,461,19.8,476,17.6,476,16.2,471,19.5,453,19.9,443,20,442,17.3),dim=c(2,60),dimnames=list(c('WK>25j','ExpBE'),1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('WK>25j','ExpBE'),1:60))
> 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 = 'Include Monthly 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
WK>25j ExpBE M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 t
1 363 14.3 1 0 0 0 0 0 0 0 0 0 0 1
2 364 14.2 0 1 0 0 0 0 0 0 0 0 0 2
3 363 15.9 0 0 1 0 0 0 0 0 0 0 0 3
4 358 15.3 0 0 0 1 0 0 0 0 0 0 0 4
5 357 15.5 0 0 0 0 1 0 0 0 0 0 0 5
6 357 15.1 0 0 0 0 0 1 0 0 0 0 0 6
7 380 15.0 0 0 0 0 0 0 1 0 0 0 0 7
8 378 12.1 0 0 0 0 0 0 0 1 0 0 0 8
9 376 15.8 0 0 0 0 0 0 0 0 1 0 0 9
10 380 16.9 0 0 0 0 0 0 0 0 0 1 0 10
11 379 15.1 0 0 0 0 0 0 0 0 0 0 1 11
12 384 13.7 0 0 0 0 0 0 0 0 0 0 0 12
13 392 14.8 1 0 0 0 0 0 0 0 0 0 0 13
14 394 14.7 0 1 0 0 0 0 0 0 0 0 0 14
15 392 16.0 0 0 1 0 0 0 0 0 0 0 0 15
16 396 15.4 0 0 0 1 0 0 0 0 0 0 0 16
17 392 15.0 0 0 0 0 1 0 0 0 0 0 0 17
18 396 15.5 0 0 0 0 0 1 0 0 0 0 0 18
19 419 15.1 0 0 0 0 0 0 1 0 0 0 0 19
20 421 11.7 0 0 0 0 0 0 0 1 0 0 0 20
21 420 16.3 0 0 0 0 0 0 0 0 1 0 0 21
22 418 16.7 0 0 0 0 0 0 0 0 0 1 0 22
23 410 15.0 0 0 0 0 0 0 0 0 0 0 1 23
24 418 14.9 0 0 0 0 0 0 0 0 0 0 0 24
25 426 14.6 1 0 0 0 0 0 0 0 0 0 0 25
26 428 15.3 0 1 0 0 0 0 0 0 0 0 0 26
27 430 17.9 0 0 1 0 0 0 0 0 0 0 0 27
28 424 16.4 0 0 0 1 0 0 0 0 0 0 0 28
29 423 15.4 0 0 0 0 1 0 0 0 0 0 0 29
30 427 17.9 0 0 0 0 0 1 0 0 0 0 0 30
31 441 15.9 0 0 0 0 0 0 1 0 0 0 0 31
32 449 13.9 0 0 0 0 0 0 0 1 0 0 0 32
33 452 17.8 0 0 0 0 0 0 0 0 1 0 0 33
34 462 17.9 0 0 0 0 0 0 0 0 0 1 0 34
35 455 17.4 0 0 0 0 0 0 0 0 0 0 1 35
36 461 16.7 0 0 0 0 0 0 0 0 0 0 0 36
37 461 16.0 1 0 0 0 0 0 0 0 0 0 0 37
38 463 16.6 0 1 0 0 0 0 0 0 0 0 0 38
39 462 19.1 0 0 1 0 0 0 0 0 0 0 0 39
40 456 17.8 0 0 0 1 0 0 0 0 0 0 0 40
41 455 17.2 0 0 0 0 1 0 0 0 0 0 0 41
42 456 18.6 0 0 0 0 0 1 0 0 0 0 0 42
43 472 16.3 0 0 0 0 0 0 1 0 0 0 0 43
44 472 15.1 0 0 0 0 0 0 0 1 0 0 0 44
45 471 19.2 0 0 0 0 0 0 0 0 1 0 0 45
46 465 17.7 0 0 0 0 0 0 0 0 0 1 0 46
47 459 19.1 0 0 0 0 0 0 0 0 0 0 1 47
48 465 18.0 0 0 0 0 0 0 0 0 0 0 0 48
49 468 17.5 1 0 0 0 0 0 0 0 0 0 0 49
50 467 17.8 0 1 0 0 0 0 0 0 0 0 0 50
51 463 21.1 0 0 1 0 0 0 0 0 0 0 0 51
52 460 17.2 0 0 0 1 0 0 0 0 0 0 0 52
53 462 19.4 0 0 0 0 1 0 0 0 0 0 0 53
54 461 19.8 0 0 0 0 0 1 0 0 0 0 0 54
55 476 17.6 0 0 0 0 0 0 1 0 0 0 0 55
56 476 16.2 0 0 0 0 0 0 0 1 0 0 0 56
57 471 19.5 0 0 0 0 0 0 0 0 1 0 0 57
58 453 19.9 0 0 0 0 0 0 0 0 0 1 0 58
59 443 20.0 0 0 0 0 0 0 0 0 0 0 1 59
60 442 17.3 0 0 0 0 0 0 0 0 0 0 0 60
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) ExpBE M1 M2 M3 M4
390.49516 -2.19378 10.60697 10.23043 11.84146 2.98449
M5 M6 M7 M8 M9 M10
-0.03080 1.30893 14.24684 8.87360 14.08243 9.71102
M11 t
0.02333 2.19079
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-41.990 -8.304 1.626 8.759 28.272
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 390.49516 41.99352 9.299 3.90e-12 ***
ExpBE -2.19378 3.13535 -0.700 0.488
M1 10.60697 9.52556 1.114 0.271
M2 10.23043 9.57727 1.068 0.291
M3 11.84146 12.53824 0.944 0.350
M4 2.98449 9.92610 0.301 0.765
M5 -0.03080 9.91429 -0.003 0.998
M6 1.30893 10.92569 0.120 0.905
M7 14.24684 9.47269 1.504 0.139
M8 8.87360 11.30781 0.785 0.437
M9 14.08243 11.06017 1.273 0.209
M10 9.71102 11.08656 0.876 0.386
M11 0.02333 10.24028 0.002 0.998
t 2.19079 0.28158 7.780 6.25e-10 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 14.89 on 46 degrees of freedom
Multiple R-squared: 0.8786, Adjusted R-squared: 0.8443
F-statistic: 25.6 on 13 and 46 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.318510e-02 2.637021e-02 0.9868149
[2,] 1.056258e-02 2.112515e-02 0.9894374
[3,] 4.722780e-03 9.445560e-03 0.9952772
[4,] 2.048355e-03 4.096710e-03 0.9979516
[5,] 3.042321e-03 6.084642e-03 0.9969577
[6,] 1.155144e-03 2.310287e-03 0.9988449
[7,] 5.495414e-04 1.099083e-03 0.9994505
[8,] 1.821230e-04 3.642460e-04 0.9998179
[9,] 1.175800e-04 2.351600e-04 0.9998824
[10,] 5.747313e-05 1.149463e-04 0.9999425
[11,] 2.645898e-05 5.291795e-05 0.9999735
[12,] 2.752818e-05 5.505635e-05 0.9999725
[13,] 2.264624e-05 4.529247e-05 0.9999774
[14,] 2.853615e-05 5.707230e-05 0.9999715
[15,] 8.084164e-04 1.616833e-03 0.9991916
[16,] 3.008594e-03 6.017188e-03 0.9969914
[17,] 1.896004e-02 3.792008e-02 0.9810400
[18,] 4.591131e-02 9.182263e-02 0.9540887
[19,] 4.777238e-02 9.554477e-02 0.9522276
[20,] 3.765733e-02 7.531465e-02 0.9623427
[21,] 2.689853e-02 5.379706e-02 0.9731015
[22,] 1.711556e-02 3.423112e-02 0.9828844
[23,] 9.099063e-03 1.819813e-02 0.9909009
[24,] 2.094170e-02 4.188339e-02 0.9790583
[25,] 1.290610e-02 2.581221e-02 0.9870939
[26,] 1.416865e-02 2.833731e-02 0.9858313
[27,] 2.020786e-02 4.041573e-02 0.9797921
> postscript(file="/var/www/html/rcomp/tmp/1zl301258736416.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/2wqmg1258736416.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/3lrsi1258736416.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/422l31258736416.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/5yoq11258736416.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 = 60
Frequency = 1
1 2 3 4 5 6 7
-8.921849 -9.955486 -11.027879 -10.677974 -10.414720 -14.822760 -7.170844
8 9 10 11 12 13 14
-12.350367 -13.632998 -5.039217 -2.491133 -2.729889 -5.114489 -5.148126
15 16 17 18 19 20 21
-8.098031 1.251874 -2.801141 -1.234778 5.759003 3.482590 5.174362
22 23 24 25 26 27 28
6.232496 1.999959 7.613118 2.157224 3.878612 7.780622 5.156124
29 30 31 32 33 34 35
2.786841 8.740766 3.224498 10.019378 14.175502 26.575502 25.975502
36 37 38 39 40 41 42
28.272393 13.938987 15.440997 16.123629 13.937887 12.446116 12.986882
43 44 45 46 47 48 49
8.812479 9.362385 9.957265 2.847216 7.415400 8.834778 -2.059872
50 51 52 53 54 55 56
-4.215997 -4.778340 -9.667912 -2.017096 -5.670111 -10.625136 -10.513987
57 58 59 60
-15.674131 -30.615997 -32.899728 -41.990400
> postscript(file="/var/www/html/rcomp/tmp/6vkug1258736416.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 = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 -8.921849 NA
1 -9.955486 -8.921849
2 -11.027879 -9.955486
3 -10.677974 -11.027879
4 -10.414720 -10.677974
5 -14.822760 -10.414720
6 -7.170844 -14.822760
7 -12.350367 -7.170844
8 -13.632998 -12.350367
9 -5.039217 -13.632998
10 -2.491133 -5.039217
11 -2.729889 -2.491133
12 -5.114489 -2.729889
13 -5.148126 -5.114489
14 -8.098031 -5.148126
15 1.251874 -8.098031
16 -2.801141 1.251874
17 -1.234778 -2.801141
18 5.759003 -1.234778
19 3.482590 5.759003
20 5.174362 3.482590
21 6.232496 5.174362
22 1.999959 6.232496
23 7.613118 1.999959
24 2.157224 7.613118
25 3.878612 2.157224
26 7.780622 3.878612
27 5.156124 7.780622
28 2.786841 5.156124
29 8.740766 2.786841
30 3.224498 8.740766
31 10.019378 3.224498
32 14.175502 10.019378
33 26.575502 14.175502
34 25.975502 26.575502
35 28.272393 25.975502
36 13.938987 28.272393
37 15.440997 13.938987
38 16.123629 15.440997
39 13.937887 16.123629
40 12.446116 13.937887
41 12.986882 12.446116
42 8.812479 12.986882
43 9.362385 8.812479
44 9.957265 9.362385
45 2.847216 9.957265
46 7.415400 2.847216
47 8.834778 7.415400
48 -2.059872 8.834778
49 -4.215997 -2.059872
50 -4.778340 -4.215997
51 -9.667912 -4.778340
52 -2.017096 -9.667912
53 -5.670111 -2.017096
54 -10.625136 -5.670111
55 -10.513987 -10.625136
56 -15.674131 -10.513987
57 -30.615997 -15.674131
58 -32.899728 -30.615997
59 -41.990400 -32.899728
60 NA -41.990400
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -9.955486 -8.921849
[2,] -11.027879 -9.955486
[3,] -10.677974 -11.027879
[4,] -10.414720 -10.677974
[5,] -14.822760 -10.414720
[6,] -7.170844 -14.822760
[7,] -12.350367 -7.170844
[8,] -13.632998 -12.350367
[9,] -5.039217 -13.632998
[10,] -2.491133 -5.039217
[11,] -2.729889 -2.491133
[12,] -5.114489 -2.729889
[13,] -5.148126 -5.114489
[14,] -8.098031 -5.148126
[15,] 1.251874 -8.098031
[16,] -2.801141 1.251874
[17,] -1.234778 -2.801141
[18,] 5.759003 -1.234778
[19,] 3.482590 5.759003
[20,] 5.174362 3.482590
[21,] 6.232496 5.174362
[22,] 1.999959 6.232496
[23,] 7.613118 1.999959
[24,] 2.157224 7.613118
[25,] 3.878612 2.157224
[26,] 7.780622 3.878612
[27,] 5.156124 7.780622
[28,] 2.786841 5.156124
[29,] 8.740766 2.786841
[30,] 3.224498 8.740766
[31,] 10.019378 3.224498
[32,] 14.175502 10.019378
[33,] 26.575502 14.175502
[34,] 25.975502 26.575502
[35,] 28.272393 25.975502
[36,] 13.938987 28.272393
[37,] 15.440997 13.938987
[38,] 16.123629 15.440997
[39,] 13.937887 16.123629
[40,] 12.446116 13.937887
[41,] 12.986882 12.446116
[42,] 8.812479 12.986882
[43,] 9.362385 8.812479
[44,] 9.957265 9.362385
[45,] 2.847216 9.957265
[46,] 7.415400 2.847216
[47,] 8.834778 7.415400
[48,] -2.059872 8.834778
[49,] -4.215997 -2.059872
[50,] -4.778340 -4.215997
[51,] -9.667912 -4.778340
[52,] -2.017096 -9.667912
[53,] -5.670111 -2.017096
[54,] -10.625136 -5.670111
[55,] -10.513987 -10.625136
[56,] -15.674131 -10.513987
[57,] -30.615997 -15.674131
[58,] -32.899728 -30.615997
[59,] -41.990400 -32.899728
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -9.955486 -8.921849
2 -11.027879 -9.955486
3 -10.677974 -11.027879
4 -10.414720 -10.677974
5 -14.822760 -10.414720
6 -7.170844 -14.822760
7 -12.350367 -7.170844
8 -13.632998 -12.350367
9 -5.039217 -13.632998
10 -2.491133 -5.039217
11 -2.729889 -2.491133
12 -5.114489 -2.729889
13 -5.148126 -5.114489
14 -8.098031 -5.148126
15 1.251874 -8.098031
16 -2.801141 1.251874
17 -1.234778 -2.801141
18 5.759003 -1.234778
19 3.482590 5.759003
20 5.174362 3.482590
21 6.232496 5.174362
22 1.999959 6.232496
23 7.613118 1.999959
24 2.157224 7.613118
25 3.878612 2.157224
26 7.780622 3.878612
27 5.156124 7.780622
28 2.786841 5.156124
29 8.740766 2.786841
30 3.224498 8.740766
31 10.019378 3.224498
32 14.175502 10.019378
33 26.575502 14.175502
34 25.975502 26.575502
35 28.272393 25.975502
36 13.938987 28.272393
37 15.440997 13.938987
38 16.123629 15.440997
39 13.937887 16.123629
40 12.446116 13.937887
41 12.986882 12.446116
42 8.812479 12.986882
43 9.362385 8.812479
44 9.957265 9.362385
45 2.847216 9.957265
46 7.415400 2.847216
47 8.834778 7.415400
48 -2.059872 8.834778
49 -4.215997 -2.059872
50 -4.778340 -4.215997
51 -9.667912 -4.778340
52 -2.017096 -9.667912
53 -5.670111 -2.017096
54 -10.625136 -5.670111
55 -10.513987 -10.625136
56 -15.674131 -10.513987
57 -30.615997 -15.674131
58 -32.899728 -30.615997
59 -41.990400 -32.899728
> 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/7733z1258736416.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/8xzw81258736416.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/9kzgs1258736416.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/104mok1258736416.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/11msy11258736416.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/12qmnf1258736416.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/13fenr1258736416.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/14v6s31258736416.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/15ufeb1258736416.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/164xxw1258736416.tab")
+ }
>
> system("convert tmp/1zl301258736416.ps tmp/1zl301258736416.png")
> system("convert tmp/2wqmg1258736416.ps tmp/2wqmg1258736416.png")
> system("convert tmp/3lrsi1258736416.ps tmp/3lrsi1258736416.png")
> system("convert tmp/422l31258736416.ps tmp/422l31258736416.png")
> system("convert tmp/5yoq11258736416.ps tmp/5yoq11258736416.png")
> system("convert tmp/6vkug1258736416.ps tmp/6vkug1258736416.png")
> system("convert tmp/7733z1258736416.ps tmp/7733z1258736416.png")
> system("convert tmp/8xzw81258736416.ps tmp/8xzw81258736416.png")
> system("convert tmp/9kzgs1258736416.ps tmp/9kzgs1258736416.png")
> system("convert tmp/104mok1258736416.ps tmp/104mok1258736416.png")
>
>
> proc.time()
user system elapsed
2.561 1.643 5.883