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(5560,543,3922,594,3759,611,4138,613,4634,611,3996,594,4308,595,4143,591,4429,589,5219,584,4929,573,5755,567,5592,569,4163,621,4962,629,5208,628,4755,612,4491,595,5732,597,5731,593,5040,590,6102,580,4904,574,5369,573,5578,573,4619,620,4731,626,5011,620,5299,588,4146,566,4625,557,4736,561,4219,549,5116,532,4205,526,4121,511,5103,499,4300,555,4578,565,3809,542,5526,527,4247,510,3830,514,4394,517,4826,508,4409,493,4569,490,4106,469,4794,478,3914,528,3793,534,4405,518,4022,506,4100,502,4788,516,3163,528,3585,533,3903,536,4178,537,3863,524,4187,536),dim=c(2,61),dimnames=list(c('Y','X'),1:61))
> y <- array(NA,dim=c(2,61),dimnames=list(c('Y','X'),1:61))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = 'No Linear Trend'
> par2 = 'Include Monthly Dummies'
> par1 = '1'
> 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 M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 5560 543 1 0 0 0 0 0 0 0 0 0 0
2 3922 594 0 1 0 0 0 0 0 0 0 0 0
3 3759 611 0 0 1 0 0 0 0 0 0 0 0
4 4138 613 0 0 0 1 0 0 0 0 0 0 0
5 4634 611 0 0 0 0 1 0 0 0 0 0 0
6 3996 594 0 0 0 0 0 1 0 0 0 0 0
7 4308 595 0 0 0 0 0 0 1 0 0 0 0
8 4143 591 0 0 0 0 0 0 0 1 0 0 0
9 4429 589 0 0 0 0 0 0 0 0 1 0 0
10 5219 584 0 0 0 0 0 0 0 0 0 1 0
11 4929 573 0 0 0 0 0 0 0 0 0 0 1
12 5755 567 0 0 0 0 0 0 0 0 0 0 0
13 5592 569 1 0 0 0 0 0 0 0 0 0 0
14 4163 621 0 1 0 0 0 0 0 0 0 0 0
15 4962 629 0 0 1 0 0 0 0 0 0 0 0
16 5208 628 0 0 0 1 0 0 0 0 0 0 0
17 4755 612 0 0 0 0 1 0 0 0 0 0 0
18 4491 595 0 0 0 0 0 1 0 0 0 0 0
19 5732 597 0 0 0 0 0 0 1 0 0 0 0
20 5731 593 0 0 0 0 0 0 0 1 0 0 0
21 5040 590 0 0 0 0 0 0 0 0 1 0 0
22 6102 580 0 0 0 0 0 0 0 0 0 1 0
23 4904 574 0 0 0 0 0 0 0 0 0 0 1
24 5369 573 0 0 0 0 0 0 0 0 0 0 0
25 5578 573 1 0 0 0 0 0 0 0 0 0 0
26 4619 620 0 1 0 0 0 0 0 0 0 0 0
27 4731 626 0 0 1 0 0 0 0 0 0 0 0
28 5011 620 0 0 0 1 0 0 0 0 0 0 0
29 5299 588 0 0 0 0 1 0 0 0 0 0 0
30 4146 566 0 0 0 0 0 1 0 0 0 0 0
31 4625 557 0 0 0 0 0 0 1 0 0 0 0
32 4736 561 0 0 0 0 0 0 0 1 0 0 0
33 4219 549 0 0 0 0 0 0 0 0 1 0 0
34 5116 532 0 0 0 0 0 0 0 0 0 1 0
35 4205 526 0 0 0 0 0 0 0 0 0 0 1
36 4121 511 0 0 0 0 0 0 0 0 0 0 0
37 5103 499 1 0 0 0 0 0 0 0 0 0 0
38 4300 555 0 1 0 0 0 0 0 0 0 0 0
39 4578 565 0 0 1 0 0 0 0 0 0 0 0
40 3809 542 0 0 0 1 0 0 0 0 0 0 0
41 5526 527 0 0 0 0 1 0 0 0 0 0 0
42 4247 510 0 0 0 0 0 1 0 0 0 0 0
43 3830 514 0 0 0 0 0 0 1 0 0 0 0
44 4394 517 0 0 0 0 0 0 0 1 0 0 0
45 4826 508 0 0 0 0 0 0 0 0 1 0 0
46 4409 493 0 0 0 0 0 0 0 0 0 1 0
47 4569 490 0 0 0 0 0 0 0 0 0 0 1
48 4106 469 0 0 0 0 0 0 0 0 0 0 0
49 4794 478 1 0 0 0 0 0 0 0 0 0 0
50 3914 528 0 1 0 0 0 0 0 0 0 0 0
51 3793 534 0 0 1 0 0 0 0 0 0 0 0
52 4405 518 0 0 0 1 0 0 0 0 0 0 0
53 4022 506 0 0 0 0 1 0 0 0 0 0 0
54 4100 502 0 0 0 0 0 1 0 0 0 0 0
55 4788 516 0 0 0 0 0 0 1 0 0 0 0
56 3163 528 0 0 0 0 0 0 0 1 0 0 0
57 3585 533 0 0 0 0 0 0 0 0 1 0 0
58 3903 536 0 0 0 0 0 0 0 0 0 1 0
59 4178 537 0 0 0 0 0 0 0 0 0 0 1
60 3863 524 0 0 0 0 0 0 0 0 0 0 0
61 4187 536 1 0 0 0 0 0 0 0 0 0 0
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X M1 M2 M3 M4
687.38 7.48 461.45 -869.10 -758.42 -542.99
M5 M6 M7 M8 M9 M10
-94.80 -630.81 -188.16 -427.82 -410.00 185.82
M11
-169.58
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1046.00 -356.56 92.68 330.33 1035.80
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 687.378 1033.969 0.665 0.509363
X 7.480 1.899 3.938 0.000265 ***
M1 461.451 332.291 1.389 0.171337
M2 -869.104 362.244 -2.399 0.020361 *
M3 -758.416 367.773 -2.062 0.044625 *
M4 -542.992 362.573 -1.498 0.140784
M5 -94.800 355.188 -0.267 0.790689
M6 -630.808 350.099 -1.802 0.077859 .
M7 -188.160 350.737 -0.536 0.594111
M8 -427.816 351.372 -1.218 0.229344
M9 -410.000 350.201 -1.171 0.247477
M10 185.824 348.329 0.533 0.596168
M11 -169.576 347.618 -0.488 0.627897
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 548.6 on 48 degrees of freedom
Multiple R-squared: 0.3998, Adjusted R-squared: 0.2498
F-statistic: 2.665 on 12 and 48 DF, p-value: 0.008006
> 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.5284286 0.94314284 0.47157142
[2,] 0.3808550 0.76171002 0.61914499
[3,] 0.2959467 0.59189343 0.70405328
[4,] 0.5577067 0.88458663 0.44229331
[5,] 0.7943569 0.41128615 0.20564307
[6,] 0.7398966 0.52020676 0.26010338
[7,] 0.8098569 0.38028620 0.19014310
[8,] 0.7307241 0.53855175 0.26927587
[9,] 0.7252882 0.54942367 0.27471184
[10,] 0.6741883 0.65162333 0.32581167
[11,] 0.5938317 0.81233661 0.40616831
[12,] 0.5068898 0.98622043 0.49311022
[13,] 0.4622785 0.92455698 0.53772151
[14,] 0.5048913 0.99021739 0.49510869
[15,] 0.4126872 0.82537444 0.58731278
[16,] 0.3300203 0.66004056 0.66997972
[17,] 0.3458614 0.69172282 0.65413859
[18,] 0.2632569 0.52651371 0.73674315
[19,] 0.2924783 0.58495661 0.70752170
[20,] 0.2173937 0.43478740 0.78260630
[21,] 0.1873726 0.37474512 0.81262744
[22,] 0.1590669 0.31813383 0.84093308
[23,] 0.1495099 0.29901989 0.85049005
[24,] 0.1829405 0.36588100 0.81705950
[25,] 0.1307690 0.26153810 0.86923095
[26,] 0.5030073 0.99398548 0.49699274
[27,] 0.4033273 0.80665451 0.59667275
[28,] 0.4937034 0.98740680 0.50629660
[29,] 0.7034519 0.59309614 0.29654807
[30,] 0.9826392 0.03472161 0.01736080
> postscript(file="/var/www/html/rcomp/tmp/1s7fo1258981457.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/24eka1258981457.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/3eyvb1258981457.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/4ecvq1258981457.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/5muwd1258981457.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 = 61
Frequency = 1
1 2 3 4 5 6
349.53337 -339.39197 -740.23994 -591.62391 -528.85586 -503.68787
7 8 9 10 11 12
-641.81587 -537.23989 -254.09589 -22.51987 125.16011 826.46412
13 14 15 16 17 18
187.05345 -300.35188 328.12012 366.17614 -415.33586 -16.16786
19 20 21 22 23 24
767.22413 1035.80011 349.42412 890.40011 92.68011 395.58414
25 26 27 28 29 30
143.13346 163.12812 119.56011 229.01612 308.18406 -144.24796
31 32 33 34 35 36
-40.57600 280.16001 -164.89602 263.43996 -247.28005 -388.65606
37 38 39 40 41 42
221.65322 330.32791 422.83991 -389.54414 991.46386 375.63186
43 44 45 46 47 48
-513.93614 267.27987 748.78385 -151.84017 385.99984 -89.49619
49 50 51 52 53 54
69.73315 146.28782 -130.28019 385.97578 -355.45620 288.47183
55 56 57 58 59 60
429.10387 -1046.00010 -679.21607 -979.48003 -356.56001 -743.89602
61
-971.10666
> postscript(file="/var/www/html/rcomp/tmp/632ne1258981457.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 = 61
Frequency = 1
lag(myerror, k = 1) myerror
0 349.53337 NA
1 -339.39197 349.53337
2 -740.23994 -339.39197
3 -591.62391 -740.23994
4 -528.85586 -591.62391
5 -503.68787 -528.85586
6 -641.81587 -503.68787
7 -537.23989 -641.81587
8 -254.09589 -537.23989
9 -22.51987 -254.09589
10 125.16011 -22.51987
11 826.46412 125.16011
12 187.05345 826.46412
13 -300.35188 187.05345
14 328.12012 -300.35188
15 366.17614 328.12012
16 -415.33586 366.17614
17 -16.16786 -415.33586
18 767.22413 -16.16786
19 1035.80011 767.22413
20 349.42412 1035.80011
21 890.40011 349.42412
22 92.68011 890.40011
23 395.58414 92.68011
24 143.13346 395.58414
25 163.12812 143.13346
26 119.56011 163.12812
27 229.01612 119.56011
28 308.18406 229.01612
29 -144.24796 308.18406
30 -40.57600 -144.24796
31 280.16001 -40.57600
32 -164.89602 280.16001
33 263.43996 -164.89602
34 -247.28005 263.43996
35 -388.65606 -247.28005
36 221.65322 -388.65606
37 330.32791 221.65322
38 422.83991 330.32791
39 -389.54414 422.83991
40 991.46386 -389.54414
41 375.63186 991.46386
42 -513.93614 375.63186
43 267.27987 -513.93614
44 748.78385 267.27987
45 -151.84017 748.78385
46 385.99984 -151.84017
47 -89.49619 385.99984
48 69.73315 -89.49619
49 146.28782 69.73315
50 -130.28019 146.28782
51 385.97578 -130.28019
52 -355.45620 385.97578
53 288.47183 -355.45620
54 429.10387 288.47183
55 -1046.00010 429.10387
56 -679.21607 -1046.00010
57 -979.48003 -679.21607
58 -356.56001 -979.48003
59 -743.89602 -356.56001
60 -971.10666 -743.89602
61 NA -971.10666
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -339.39197 349.53337
[2,] -740.23994 -339.39197
[3,] -591.62391 -740.23994
[4,] -528.85586 -591.62391
[5,] -503.68787 -528.85586
[6,] -641.81587 -503.68787
[7,] -537.23989 -641.81587
[8,] -254.09589 -537.23989
[9,] -22.51987 -254.09589
[10,] 125.16011 -22.51987
[11,] 826.46412 125.16011
[12,] 187.05345 826.46412
[13,] -300.35188 187.05345
[14,] 328.12012 -300.35188
[15,] 366.17614 328.12012
[16,] -415.33586 366.17614
[17,] -16.16786 -415.33586
[18,] 767.22413 -16.16786
[19,] 1035.80011 767.22413
[20,] 349.42412 1035.80011
[21,] 890.40011 349.42412
[22,] 92.68011 890.40011
[23,] 395.58414 92.68011
[24,] 143.13346 395.58414
[25,] 163.12812 143.13346
[26,] 119.56011 163.12812
[27,] 229.01612 119.56011
[28,] 308.18406 229.01612
[29,] -144.24796 308.18406
[30,] -40.57600 -144.24796
[31,] 280.16001 -40.57600
[32,] -164.89602 280.16001
[33,] 263.43996 -164.89602
[34,] -247.28005 263.43996
[35,] -388.65606 -247.28005
[36,] 221.65322 -388.65606
[37,] 330.32791 221.65322
[38,] 422.83991 330.32791
[39,] -389.54414 422.83991
[40,] 991.46386 -389.54414
[41,] 375.63186 991.46386
[42,] -513.93614 375.63186
[43,] 267.27987 -513.93614
[44,] 748.78385 267.27987
[45,] -151.84017 748.78385
[46,] 385.99984 -151.84017
[47,] -89.49619 385.99984
[48,] 69.73315 -89.49619
[49,] 146.28782 69.73315
[50,] -130.28019 146.28782
[51,] 385.97578 -130.28019
[52,] -355.45620 385.97578
[53,] 288.47183 -355.45620
[54,] 429.10387 288.47183
[55,] -1046.00010 429.10387
[56,] -679.21607 -1046.00010
[57,] -979.48003 -679.21607
[58,] -356.56001 -979.48003
[59,] -743.89602 -356.56001
[60,] -971.10666 -743.89602
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -339.39197 349.53337
2 -740.23994 -339.39197
3 -591.62391 -740.23994
4 -528.85586 -591.62391
5 -503.68787 -528.85586
6 -641.81587 -503.68787
7 -537.23989 -641.81587
8 -254.09589 -537.23989
9 -22.51987 -254.09589
10 125.16011 -22.51987
11 826.46412 125.16011
12 187.05345 826.46412
13 -300.35188 187.05345
14 328.12012 -300.35188
15 366.17614 328.12012
16 -415.33586 366.17614
17 -16.16786 -415.33586
18 767.22413 -16.16786
19 1035.80011 767.22413
20 349.42412 1035.80011
21 890.40011 349.42412
22 92.68011 890.40011
23 395.58414 92.68011
24 143.13346 395.58414
25 163.12812 143.13346
26 119.56011 163.12812
27 229.01612 119.56011
28 308.18406 229.01612
29 -144.24796 308.18406
30 -40.57600 -144.24796
31 280.16001 -40.57600
32 -164.89602 280.16001
33 263.43996 -164.89602
34 -247.28005 263.43996
35 -388.65606 -247.28005
36 221.65322 -388.65606
37 330.32791 221.65322
38 422.83991 330.32791
39 -389.54414 422.83991
40 991.46386 -389.54414
41 375.63186 991.46386
42 -513.93614 375.63186
43 267.27987 -513.93614
44 748.78385 267.27987
45 -151.84017 748.78385
46 385.99984 -151.84017
47 -89.49619 385.99984
48 69.73315 -89.49619
49 146.28782 69.73315
50 -130.28019 146.28782
51 385.97578 -130.28019
52 -355.45620 385.97578
53 288.47183 -355.45620
54 429.10387 288.47183
55 -1046.00010 429.10387
56 -679.21607 -1046.00010
57 -979.48003 -679.21607
58 -356.56001 -979.48003
59 -743.89602 -356.56001
60 -971.10666 -743.89602
> 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/7ppyp1258981457.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/83axo1258981457.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/9e7531258981458.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/10b9bk1258981458.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/11t4dn1258981458.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/12jcy41258981458.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/1333jf1258981458.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/140hpg1258981458.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/15u2eu1258981458.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/160hx31258981458.tab")
+ }
>
> system("convert tmp/1s7fo1258981457.ps tmp/1s7fo1258981457.png")
> system("convert tmp/24eka1258981457.ps tmp/24eka1258981457.png")
> system("convert tmp/3eyvb1258981457.ps tmp/3eyvb1258981457.png")
> system("convert tmp/4ecvq1258981457.ps tmp/4ecvq1258981457.png")
> system("convert tmp/5muwd1258981457.ps tmp/5muwd1258981457.png")
> system("convert tmp/632ne1258981457.ps tmp/632ne1258981457.png")
> system("convert tmp/7ppyp1258981457.ps tmp/7ppyp1258981457.png")
> system("convert tmp/83axo1258981457.ps tmp/83axo1258981457.png")
> system("convert tmp/9e7531258981458.ps tmp/9e7531258981458.png")
> system("convert tmp/10b9bk1258981458.ps tmp/10b9bk1258981458.png")
>
>
> proc.time()
user system elapsed
2.407 1.542 3.730