R version 2.6.2 (2008-02-08)
Copyright (C) 2008 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.
Natural language support but running in an English locale
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(56421,53152,53536,52408,41454,38271,35306,26414,31917,38030,27534,18387,50556,43901,48572,43899,37532,40357,35489,29027,34485,42598,30306,26451,47460,50104,61465,53726,39477,43895,31481,29896,33842,39120,33702,25094,51442,45594,52518,48564,41745,49585,32747,33379,35645,37034,35681,20972,58552,54955,65540,51570,51145,46641,35704,33253,35193,41668,34865,21210,56126,49231,59723,48103,47472,50497,40059,34149,36860,46356,36577),dim=c(1,71),dimnames=list(c('geg'),1:71))
> y <- array(NA,dim=c(1,71),dimnames=list(c('geg'),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 = 'No Linear Trend'
> par2 = 'Include Quarterly 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
> 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
geg Q1 Q2 Q3
1 56421 1 0 0
2 53152 0 1 0
3 53536 0 0 1
4 52408 0 0 0
5 41454 1 0 0
6 38271 0 1 0
7 35306 0 0 1
8 26414 0 0 0
9 31917 1 0 0
10 38030 0 1 0
11 27534 0 0 1
12 18387 0 0 0
13 50556 1 0 0
14 43901 0 1 0
15 48572 0 0 1
16 43899 0 0 0
17 37532 1 0 0
18 40357 0 1 0
19 35489 0 0 1
20 29027 0 0 0
21 34485 1 0 0
22 42598 0 1 0
23 30306 0 0 1
24 26451 0 0 0
25 47460 1 0 0
26 50104 0 1 0
27 61465 0 0 1
28 53726 0 0 0
29 39477 1 0 0
30 43895 0 1 0
31 31481 0 0 1
32 29896 0 0 0
33 33842 1 0 0
34 39120 0 1 0
35 33702 0 0 1
36 25094 0 0 0
37 51442 1 0 0
38 45594 0 1 0
39 52518 0 0 1
40 48564 0 0 0
41 41745 1 0 0
42 49585 0 1 0
43 32747 0 0 1
44 33379 0 0 0
45 35645 1 0 0
46 37034 0 1 0
47 35681 0 0 1
48 20972 0 0 0
49 58552 1 0 0
50 54955 0 1 0
51 65540 0 0 1
52 51570 0 0 0
53 51145 1 0 0
54 46641 0 1 0
55 35704 0 0 1
56 33253 0 0 0
57 35193 1 0 0
58 41668 0 1 0
59 34865 0 0 1
60 21210 0 0 0
61 56126 1 0 0
62 49231 0 1 0
63 59723 0 0 1
64 48103 0 0 0
65 47472 1 0 0
66 50497 0 1 0
67 40059 0 0 1
68 34149 0 0 0
69 36860 1 0 0
70 46356 0 1 0
71 36577 0 0 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Q1 Q2 Q3
35088 8652 9967 6623
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-16701 -6953 -1995 7133 23829
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 35088 2394 14.655 < 2e-16 ***
Q1 8652 3339 2.591 0.01172 *
Q2 9967 3339 2.985 0.00395 **
Q3 6623 3339 1.984 0.05138 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 9872 on 67 degrees of freedom
Multiple R-squared: 0.1348, Adjusted R-squared: 0.09609
F-statistic: 3.48 on 3 and 67 DF, p-value: 0.02059
> 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.7345714 0.5308572 0.26542861
[2,] 0.8752745 0.2494511 0.12472554
[3,] 0.8956968 0.2086063 0.10430317
[4,] 0.8517387 0.2965226 0.14826129
[5,] 0.8810441 0.2379118 0.11895589
[6,] 0.9337085 0.1325829 0.06629146
[7,] 0.9090343 0.1819314 0.09096570
[8,] 0.8631987 0.2736027 0.13680133
[9,] 0.8380651 0.3238697 0.16193486
[10,] 0.8256166 0.3487668 0.17438339
[11,] 0.7875410 0.4249180 0.21245899
[12,] 0.7276823 0.5446354 0.27231771
[13,] 0.6745831 0.6508338 0.32541690
[14,] 0.6208953 0.7582093 0.37910465
[15,] 0.5978382 0.8043236 0.40216181
[16,] 0.5209265 0.9581470 0.47907351
[17,] 0.5196396 0.9607208 0.48036038
[18,] 0.4876766 0.9753533 0.51232336
[19,] 0.4259939 0.8519877 0.57400615
[20,] 0.3790407 0.7580813 0.62095934
[21,] 0.6061062 0.7877876 0.39389382
[22,] 0.7599183 0.4801633 0.24008166
[23,] 0.7092035 0.5815930 0.29079652
[24,] 0.6447371 0.7105258 0.35526289
[25,] 0.6426045 0.7147910 0.35739551
[26,] 0.5915115 0.8169771 0.40848855
[27,] 0.5867820 0.8264361 0.41321804
[28,] 0.5402884 0.9194231 0.45971155
[29,] 0.5138411 0.9723179 0.48615893
[30,] 0.5129263 0.9741473 0.48707367
[31,] 0.4834747 0.9669494 0.51652529
[32,] 0.4142603 0.8285205 0.58573973
[33,] 0.4224966 0.8449931 0.57750344
[34,] 0.4746082 0.9492163 0.52539183
[35,] 0.4085069 0.8170139 0.59149307
[36,] 0.3512348 0.7024697 0.64876517
[37,] 0.3415842 0.6831683 0.65841584
[38,] 0.2778843 0.5557686 0.72211570
[39,] 0.2682301 0.5364603 0.73176987
[40,] 0.2553582 0.5107163 0.74464185
[41,] 0.2285587 0.4571173 0.77144135
[42,] 0.2959378 0.5918757 0.70406217
[43,] 0.3491054 0.6982108 0.65089462
[44,] 0.3280058 0.6560117 0.67199416
[45,] 0.6487131 0.7025739 0.35128693
[46,] 0.7640500 0.4719001 0.23595005
[47,] 0.7261394 0.5477212 0.27386060
[48,] 0.6464217 0.7071565 0.35357825
[49,] 0.5924856 0.8150287 0.40751435
[50,] 0.4994469 0.9988937 0.50055314
[51,] 0.4906673 0.9813347 0.50933265
[52,] 0.4200764 0.8401527 0.57992364
[53,] 0.3904478 0.7808956 0.60955220
[54,] 0.5541025 0.8917951 0.44589754
[55,] 0.5984796 0.8030408 0.40152042
[56,] 0.4665806 0.9331612 0.53341942
[57,] 0.8085511 0.3828977 0.19144887
[58,] 0.8809629 0.2380741 0.11903705
> postscript(file="/var/www/html/rcomp/tmp/1cpk91210256070.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/2punu1210256070.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/3rsp41210256070.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/4ihjr1210256070.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/5mqd11210256070.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
12680.7778 8097.0556 11824.6111 17319.6471 -2286.2222 -6783.9444
7 8 9 10 11 12
-6405.3889 -8674.3529 -11823.2222 -7024.9444 -14177.3889 -16701.3529
13 14 15 16 17 18
6815.7778 -1153.9444 6860.6111 8810.6471 -6208.2222 -4697.9444
19 20 21 22 23 24
-6222.3889 -6061.3529 -9255.2222 -2456.9444 -11405.3889 -8637.3529
25 26 27 28 29 30
3719.7778 5049.0556 19753.6111 18637.6471 -4263.2222 -1159.9444
31 32 33 34 35 36
-10230.3889 -5192.3529 -9898.2222 -5934.9444 -8009.3889 -9994.3529
37 38 39 40 41 42
7701.7778 539.0556 10806.6111 13475.6471 -1995.2222 4530.0556
43 44 45 46 47 48
-8964.3889 -1709.3529 -8095.2222 -8020.9444 -6030.3889 -14116.3529
49 50 51 52 53 54
14811.7778 9900.0556 23828.6111 16481.6471 7404.7778 1586.0556
55 56 57 58 59 60
-6007.3889 -1835.3529 -8547.2222 -3386.9444 -6846.3889 -13878.3529
61 62 63 64 65 66
12385.7778 4176.0556 18011.6111 13014.6471 3731.7778 5442.0556
67 68 69 70 71
-1652.3889 -939.3529 -6880.2222 1301.0556 -5134.3889
> postscript(file="/var/www/html/rcomp/tmp/633ak1210256070.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 12680.7778 NA
1 8097.0556 12680.7778
2 11824.6111 8097.0556
3 17319.6471 11824.6111
4 -2286.2222 17319.6471
5 -6783.9444 -2286.2222
6 -6405.3889 -6783.9444
7 -8674.3529 -6405.3889
8 -11823.2222 -8674.3529
9 -7024.9444 -11823.2222
10 -14177.3889 -7024.9444
11 -16701.3529 -14177.3889
12 6815.7778 -16701.3529
13 -1153.9444 6815.7778
14 6860.6111 -1153.9444
15 8810.6471 6860.6111
16 -6208.2222 8810.6471
17 -4697.9444 -6208.2222
18 -6222.3889 -4697.9444
19 -6061.3529 -6222.3889
20 -9255.2222 -6061.3529
21 -2456.9444 -9255.2222
22 -11405.3889 -2456.9444
23 -8637.3529 -11405.3889
24 3719.7778 -8637.3529
25 5049.0556 3719.7778
26 19753.6111 5049.0556
27 18637.6471 19753.6111
28 -4263.2222 18637.6471
29 -1159.9444 -4263.2222
30 -10230.3889 -1159.9444
31 -5192.3529 -10230.3889
32 -9898.2222 -5192.3529
33 -5934.9444 -9898.2222
34 -8009.3889 -5934.9444
35 -9994.3529 -8009.3889
36 7701.7778 -9994.3529
37 539.0556 7701.7778
38 10806.6111 539.0556
39 13475.6471 10806.6111
40 -1995.2222 13475.6471
41 4530.0556 -1995.2222
42 -8964.3889 4530.0556
43 -1709.3529 -8964.3889
44 -8095.2222 -1709.3529
45 -8020.9444 -8095.2222
46 -6030.3889 -8020.9444
47 -14116.3529 -6030.3889
48 14811.7778 -14116.3529
49 9900.0556 14811.7778
50 23828.6111 9900.0556
51 16481.6471 23828.6111
52 7404.7778 16481.6471
53 1586.0556 7404.7778
54 -6007.3889 1586.0556
55 -1835.3529 -6007.3889
56 -8547.2222 -1835.3529
57 -3386.9444 -8547.2222
58 -6846.3889 -3386.9444
59 -13878.3529 -6846.3889
60 12385.7778 -13878.3529
61 4176.0556 12385.7778
62 18011.6111 4176.0556
63 13014.6471 18011.6111
64 3731.7778 13014.6471
65 5442.0556 3731.7778
66 -1652.3889 5442.0556
67 -939.3529 -1652.3889
68 -6880.2222 -939.3529
69 1301.0556 -6880.2222
70 -5134.3889 1301.0556
71 NA -5134.3889
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 8097.0556 12680.7778
[2,] 11824.6111 8097.0556
[3,] 17319.6471 11824.6111
[4,] -2286.2222 17319.6471
[5,] -6783.9444 -2286.2222
[6,] -6405.3889 -6783.9444
[7,] -8674.3529 -6405.3889
[8,] -11823.2222 -8674.3529
[9,] -7024.9444 -11823.2222
[10,] -14177.3889 -7024.9444
[11,] -16701.3529 -14177.3889
[12,] 6815.7778 -16701.3529
[13,] -1153.9444 6815.7778
[14,] 6860.6111 -1153.9444
[15,] 8810.6471 6860.6111
[16,] -6208.2222 8810.6471
[17,] -4697.9444 -6208.2222
[18,] -6222.3889 -4697.9444
[19,] -6061.3529 -6222.3889
[20,] -9255.2222 -6061.3529
[21,] -2456.9444 -9255.2222
[22,] -11405.3889 -2456.9444
[23,] -8637.3529 -11405.3889
[24,] 3719.7778 -8637.3529
[25,] 5049.0556 3719.7778
[26,] 19753.6111 5049.0556
[27,] 18637.6471 19753.6111
[28,] -4263.2222 18637.6471
[29,] -1159.9444 -4263.2222
[30,] -10230.3889 -1159.9444
[31,] -5192.3529 -10230.3889
[32,] -9898.2222 -5192.3529
[33,] -5934.9444 -9898.2222
[34,] -8009.3889 -5934.9444
[35,] -9994.3529 -8009.3889
[36,] 7701.7778 -9994.3529
[37,] 539.0556 7701.7778
[38,] 10806.6111 539.0556
[39,] 13475.6471 10806.6111
[40,] -1995.2222 13475.6471
[41,] 4530.0556 -1995.2222
[42,] -8964.3889 4530.0556
[43,] -1709.3529 -8964.3889
[44,] -8095.2222 -1709.3529
[45,] -8020.9444 -8095.2222
[46,] -6030.3889 -8020.9444
[47,] -14116.3529 -6030.3889
[48,] 14811.7778 -14116.3529
[49,] 9900.0556 14811.7778
[50,] 23828.6111 9900.0556
[51,] 16481.6471 23828.6111
[52,] 7404.7778 16481.6471
[53,] 1586.0556 7404.7778
[54,] -6007.3889 1586.0556
[55,] -1835.3529 -6007.3889
[56,] -8547.2222 -1835.3529
[57,] -3386.9444 -8547.2222
[58,] -6846.3889 -3386.9444
[59,] -13878.3529 -6846.3889
[60,] 12385.7778 -13878.3529
[61,] 4176.0556 12385.7778
[62,] 18011.6111 4176.0556
[63,] 13014.6471 18011.6111
[64,] 3731.7778 13014.6471
[65,] 5442.0556 3731.7778
[66,] -1652.3889 5442.0556
[67,] -939.3529 -1652.3889
[68,] -6880.2222 -939.3529
[69,] 1301.0556 -6880.2222
[70,] -5134.3889 1301.0556
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 8097.0556 12680.7778
2 11824.6111 8097.0556
3 17319.6471 11824.6111
4 -2286.2222 17319.6471
5 -6783.9444 -2286.2222
6 -6405.3889 -6783.9444
7 -8674.3529 -6405.3889
8 -11823.2222 -8674.3529
9 -7024.9444 -11823.2222
10 -14177.3889 -7024.9444
11 -16701.3529 -14177.3889
12 6815.7778 -16701.3529
13 -1153.9444 6815.7778
14 6860.6111 -1153.9444
15 8810.6471 6860.6111
16 -6208.2222 8810.6471
17 -4697.9444 -6208.2222
18 -6222.3889 -4697.9444
19 -6061.3529 -6222.3889
20 -9255.2222 -6061.3529
21 -2456.9444 -9255.2222
22 -11405.3889 -2456.9444
23 -8637.3529 -11405.3889
24 3719.7778 -8637.3529
25 5049.0556 3719.7778
26 19753.6111 5049.0556
27 18637.6471 19753.6111
28 -4263.2222 18637.6471
29 -1159.9444 -4263.2222
30 -10230.3889 -1159.9444
31 -5192.3529 -10230.3889
32 -9898.2222 -5192.3529
33 -5934.9444 -9898.2222
34 -8009.3889 -5934.9444
35 -9994.3529 -8009.3889
36 7701.7778 -9994.3529
37 539.0556 7701.7778
38 10806.6111 539.0556
39 13475.6471 10806.6111
40 -1995.2222 13475.6471
41 4530.0556 -1995.2222
42 -8964.3889 4530.0556
43 -1709.3529 -8964.3889
44 -8095.2222 -1709.3529
45 -8020.9444 -8095.2222
46 -6030.3889 -8020.9444
47 -14116.3529 -6030.3889
48 14811.7778 -14116.3529
49 9900.0556 14811.7778
50 23828.6111 9900.0556
51 16481.6471 23828.6111
52 7404.7778 16481.6471
53 1586.0556 7404.7778
54 -6007.3889 1586.0556
55 -1835.3529 -6007.3889
56 -8547.2222 -1835.3529
57 -3386.9444 -8547.2222
58 -6846.3889 -3386.9444
59 -13878.3529 -6846.3889
60 12385.7778 -13878.3529
61 4176.0556 12385.7778
62 18011.6111 4176.0556
63 13014.6471 18011.6111
64 3731.7778 13014.6471
65 5442.0556 3731.7778
66 -1652.3889 5442.0556
67 -939.3529 -1652.3889
68 -6880.2222 -939.3529
69 1301.0556 -6880.2222
70 -5134.3889 1301.0556
> 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/7rvqt1210256070.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/8a6cv1210256070.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/9i2ln1210256070.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/10v03d1210256070.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/11yn781210256070.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/12mtjk1210256071.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/13plh51210256071.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/148bq91210256071.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/153w161210256071.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/1614ay1210256071.tab")
+ }
>
> system("convert tmp/1cpk91210256070.ps tmp/1cpk91210256070.png")
> system("convert tmp/2punu1210256070.ps tmp/2punu1210256070.png")
> system("convert tmp/3rsp41210256070.ps tmp/3rsp41210256070.png")
> system("convert tmp/4ihjr1210256070.ps tmp/4ihjr1210256070.png")
> system("convert tmp/5mqd11210256070.ps tmp/5mqd11210256070.png")
> system("convert tmp/633ak1210256070.ps tmp/633ak1210256070.png")
> system("convert tmp/7rvqt1210256070.ps tmp/7rvqt1210256070.png")
> system("convert tmp/8a6cv1210256070.ps tmp/8a6cv1210256070.png")
> system("convert tmp/9i2ln1210256070.ps tmp/9i2ln1210256070.png")
> system("convert tmp/10v03d1210256070.ps tmp/10v03d1210256070.png")
>
>
> proc.time()
user system elapsed
5.396 2.716 5.798