R version 2.13.0 (2011-04-13)
Copyright (C) 2011 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: i486-pc-linux-gnu (32-bit)
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(4290.89,-1,2.1,4443.91,1,1.7,4502.64,-1,1.8,4356.98,2,1.8,4591.27,2,1.8,4696.96,1,1.3,4621.4,-1,1.3,4562.84,-2,1.3,4202.52,-2,1.2,4296.49,-1,1.4,4435.23,-8,2.2,4105.18,-4,2.9,4116.68,-6,3.1,3844.49,-3,3.5,3720.98,-3,3.6,3674.4,-7,4.4,3857.62,-9,4.1,3801.06,-11,5.1,3504.37,-13,5.8,3032.6,-11,5.9,3047.03,-9,5.4,2962.34,-17,5.5,2197.82,-22,4.8,2014.45,-25,3.2,1862.83,-20,2.7,1905.41,-24,2.1,1810.99,-24,1.9,1670.07,-22,0.6,1864.44,-19,0.7,2052.02,-18,-0.2,2029.6,-17,-1,2070.83,-11,-1.7,2293.41,-11,-0.7,2443.27,-12,-1,2513.17,-10,-0.9,2466.92,-15,0,2502.66,-15,0.3,2539.91,-15,0.8,2482.6,-13,0.8,2626.15,-8,1.9,2656.32,-13,2.1,2446.66,-9,2.5,2467.38,-7,2.7,2462.32,-4,2.4,2504.58,-4,2.4,2579.39,-2,2.9,2649.24,0,3.1,2636.87,-2,3,2613.94,-3,3.4,2634.01,1,3.7,2711.94,-2,3.5,2646.43,-1,3.5,2717.79,1,3.3,2701.54,-3,3.1,2572.98,-4,3.4,2488.92,-9,4,2204.91,-9,3.4,2123.99,-7,3.4,2149.1,-14,3.4,2036.71,-12,3.7),dim=c(3,60),dimnames=list(c('BEL','CON','INF
'),1:60))
> y <- array(NA,dim=c(3,60),dimnames=list(c('BEL','CON','INF
'),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 = '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
> 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
BEL CON INF\r
1 4290.89 -1 2.1
2 4443.91 1 1.7
3 4502.64 -1 1.8
4 4356.98 2 1.8
5 4591.27 2 1.8
6 4696.96 1 1.3
7 4621.40 -1 1.3
8 4562.84 -2 1.3
9 4202.52 -2 1.2
10 4296.49 -1 1.4
11 4435.23 -8 2.2
12 4105.18 -4 2.9
13 4116.68 -6 3.1
14 3844.49 -3 3.5
15 3720.98 -3 3.6
16 3674.40 -7 4.4
17 3857.62 -9 4.1
18 3801.06 -11 5.1
19 3504.37 -13 5.8
20 3032.60 -11 5.9
21 3047.03 -9 5.4
22 2962.34 -17 5.5
23 2197.82 -22 4.8
24 2014.45 -25 3.2
25 1862.83 -20 2.7
26 1905.41 -24 2.1
27 1810.99 -24 1.9
28 1670.07 -22 0.6
29 1864.44 -19 0.7
30 2052.02 -18 -0.2
31 2029.60 -17 -1.0
32 2070.83 -11 -1.7
33 2293.41 -11 -0.7
34 2443.27 -12 -1.0
35 2513.17 -10 -0.9
36 2466.92 -15 0.0
37 2502.66 -15 0.3
38 2539.91 -15 0.8
39 2482.60 -13 0.8
40 2626.15 -8 1.9
41 2656.32 -13 2.1
42 2446.66 -9 2.5
43 2467.38 -7 2.7
44 2462.32 -4 2.4
45 2504.58 -4 2.4
46 2579.39 -2 2.9
47 2649.24 0 3.1
48 2636.87 -2 3.0
49 2613.94 -3 3.4
50 2634.01 1 3.7
51 2711.94 -2 3.5
52 2646.43 -1 3.5
53 2717.79 1 3.3
54 2701.54 -3 3.1
55 2572.98 -4 3.4
56 2488.92 -9 4.0
57 2204.91 -9 3.4
58 2123.99 -7 3.4
59 2149.10 -14 3.4
60 2036.71 -12 3.7
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) CON `INF\r`
3575.03 78.67 23.56
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1106.87 -635.81 -17.91 632.57 1437.71
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3575.03 196.47 18.197 < 2e-16 ***
CON 78.67 12.46 6.315 4.36e-08 ***
`INF\r` 23.56 52.83 0.446 0.657
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 702.8 on 57 degrees of freedom
Multiple R-squared: 0.4201, Adjusted R-squared: 0.3997
F-statistic: 20.65 on 2 and 57 DF, p-value: 1.801e-07
> 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.0056985526 1.139711e-02 9.943014e-01
[2,] 0.0010186423 2.037285e-03 9.989814e-01
[3,] 0.0002131949 4.263897e-04 9.997868e-01
[4,] 0.0020324124 4.064825e-03 9.979676e-01
[5,] 0.0012537245 2.507449e-03 9.987463e-01
[6,] 0.0015499161 3.099832e-03 9.984501e-01
[7,] 0.0012399710 2.479942e-03 9.987600e-01
[8,] 0.0009506558 1.901312e-03 9.990493e-01
[9,] 0.0010172723 2.034545e-03 9.989827e-01
[10,] 0.0011553366 2.310673e-03 9.988447e-01
[11,] 0.0008590323 1.718065e-03 9.991410e-01
[12,] 0.0017704869 3.540974e-03 9.982295e-01
[13,] 0.0114275035 2.285501e-02 9.885725e-01
[14,] 0.0308564831 6.171297e-02 9.691435e-01
[15,] 0.0896526278 1.793053e-01 9.103474e-01
[16,] 0.2345450428 4.690901e-01 7.654550e-01
[17,] 0.7329429395 5.341141e-01 2.670571e-01
[18,] 0.9866121279 2.677574e-02 1.338787e-02
[19,] 0.9950453094 9.909381e-03 4.954691e-03
[20,] 0.9992380079 1.523984e-03 7.619921e-04
[21,] 0.9987590875 2.481825e-03 1.240913e-03
[22,] 0.9980310117 3.937977e-03 1.968988e-03
[23,] 0.9991804248 1.639150e-03 8.195752e-04
[24,] 0.9996027229 7.945542e-04 3.972771e-04
[25,] 0.9995082623 9.834755e-04 4.917377e-04
[26,] 0.9996092761 7.814478e-04 3.907239e-04
[27,] 0.9999892534 2.149310e-05 1.074655e-05
[28,] 0.9999967543 6.491378e-06 3.245689e-06
[29,] 0.9999956468 8.706439e-06 4.353220e-06
[30,] 0.9999972163 5.567339e-06 2.783669e-06
[31,] 0.9999933669 1.326611e-05 6.633055e-06
[32,] 0.9999822210 3.555794e-05 1.777897e-05
[33,] 0.9999575613 8.487740e-05 4.243870e-05
[34,] 0.9999045379 1.909243e-04 9.546214e-05
[35,] 0.9999145908 1.708185e-04 8.540925e-05
[36,] 0.9999934757 1.304863e-05 6.524314e-06
[37,] 0.9999959632 8.073678e-06 4.036839e-06
[38,] 0.9999970712 5.857631e-06 2.928815e-06
[39,] 0.9999974314 5.137202e-06 2.568601e-06
[40,] 0.9999960103 7.979471e-06 3.989735e-06
[41,] 0.9999932619 1.347623e-05 6.738116e-06
[42,] 0.9999872451 2.550985e-05 1.275492e-05
[43,] 0.9999632703 7.345948e-05 3.672974e-05
[44,] 0.9998852245 2.295510e-04 1.147755e-04
[45,] 0.9997781933 4.436134e-04 2.218067e-04
[46,] 0.9992097696 1.580461e-03 7.902304e-04
[47,] 0.9970125354 5.974929e-03 2.987465e-03
[48,] 0.9894010672 2.119787e-02 1.059893e-02
[49,] 0.9758042304 4.839154e-02 2.419577e-02
> postscript(file="/var/wessaorg/rcomp/tmp/1qnyf1324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/tmp/2o8fm1324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/tmp/33cgc1324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/tmp/4gtm01324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/tmp/5lnrx1324479296.ps",horizontal=F,onefile=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
745.04574 750.15392 963.86439 582.19936 816.48936 1012.62879
7 8 9 10 11 12
1094.40548 1114.51382 756.55004 767.13926 1437.70793 776.49103
13 14 15 16 17 18
940.61529 422.99539 297.12917 546.37281 893.99815 971.21266
19 20 21 22 23 24
815.36583 183.90292 52.77732 595.07786 240.39309 330.72760
25 26 27 28 29 30
-202.45303 168.93765 79.23008 -188.39578 -232.38703 -102.26942
31 32 33 34 35 36
-184.50803 -598.79457 -399.77674 -164.17975 -253.97265 71.91311
37 38 39 40 41 42
100.58446 126.05337 -88.59331 -364.30342 54.49586 -479.26238
43 44 45 46 47 48
-620.59150 -854.58788 -812.32788 -906.63566 -998.83478 -851.51187
49 50 51 52 53 54
-805.19840 -1106.87042 -788.22296 -932.40130 -1013.66556 -710.52975
55 56 57 58 59 60
-767.49005 -472.34564 -742.21833 -980.47502 -404.68662 -681.48195
> postscript(file="/var/wessaorg/rcomp/tmp/6lt0i1324479296.ps",horizontal=F,onefile=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 745.04574 NA
1 750.15392 745.04574
2 963.86439 750.15392
3 582.19936 963.86439
4 816.48936 582.19936
5 1012.62879 816.48936
6 1094.40548 1012.62879
7 1114.51382 1094.40548
8 756.55004 1114.51382
9 767.13926 756.55004
10 1437.70793 767.13926
11 776.49103 1437.70793
12 940.61529 776.49103
13 422.99539 940.61529
14 297.12917 422.99539
15 546.37281 297.12917
16 893.99815 546.37281
17 971.21266 893.99815
18 815.36583 971.21266
19 183.90292 815.36583
20 52.77732 183.90292
21 595.07786 52.77732
22 240.39309 595.07786
23 330.72760 240.39309
24 -202.45303 330.72760
25 168.93765 -202.45303
26 79.23008 168.93765
27 -188.39578 79.23008
28 -232.38703 -188.39578
29 -102.26942 -232.38703
30 -184.50803 -102.26942
31 -598.79457 -184.50803
32 -399.77674 -598.79457
33 -164.17975 -399.77674
34 -253.97265 -164.17975
35 71.91311 -253.97265
36 100.58446 71.91311
37 126.05337 100.58446
38 -88.59331 126.05337
39 -364.30342 -88.59331
40 54.49586 -364.30342
41 -479.26238 54.49586
42 -620.59150 -479.26238
43 -854.58788 -620.59150
44 -812.32788 -854.58788
45 -906.63566 -812.32788
46 -998.83478 -906.63566
47 -851.51187 -998.83478
48 -805.19840 -851.51187
49 -1106.87042 -805.19840
50 -788.22296 -1106.87042
51 -932.40130 -788.22296
52 -1013.66556 -932.40130
53 -710.52975 -1013.66556
54 -767.49005 -710.52975
55 -472.34564 -767.49005
56 -742.21833 -472.34564
57 -980.47502 -742.21833
58 -404.68662 -980.47502
59 -681.48195 -404.68662
60 NA -681.48195
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 750.15392 745.04574
[2,] 963.86439 750.15392
[3,] 582.19936 963.86439
[4,] 816.48936 582.19936
[5,] 1012.62879 816.48936
[6,] 1094.40548 1012.62879
[7,] 1114.51382 1094.40548
[8,] 756.55004 1114.51382
[9,] 767.13926 756.55004
[10,] 1437.70793 767.13926
[11,] 776.49103 1437.70793
[12,] 940.61529 776.49103
[13,] 422.99539 940.61529
[14,] 297.12917 422.99539
[15,] 546.37281 297.12917
[16,] 893.99815 546.37281
[17,] 971.21266 893.99815
[18,] 815.36583 971.21266
[19,] 183.90292 815.36583
[20,] 52.77732 183.90292
[21,] 595.07786 52.77732
[22,] 240.39309 595.07786
[23,] 330.72760 240.39309
[24,] -202.45303 330.72760
[25,] 168.93765 -202.45303
[26,] 79.23008 168.93765
[27,] -188.39578 79.23008
[28,] -232.38703 -188.39578
[29,] -102.26942 -232.38703
[30,] -184.50803 -102.26942
[31,] -598.79457 -184.50803
[32,] -399.77674 -598.79457
[33,] -164.17975 -399.77674
[34,] -253.97265 -164.17975
[35,] 71.91311 -253.97265
[36,] 100.58446 71.91311
[37,] 126.05337 100.58446
[38,] -88.59331 126.05337
[39,] -364.30342 -88.59331
[40,] 54.49586 -364.30342
[41,] -479.26238 54.49586
[42,] -620.59150 -479.26238
[43,] -854.58788 -620.59150
[44,] -812.32788 -854.58788
[45,] -906.63566 -812.32788
[46,] -998.83478 -906.63566
[47,] -851.51187 -998.83478
[48,] -805.19840 -851.51187
[49,] -1106.87042 -805.19840
[50,] -788.22296 -1106.87042
[51,] -932.40130 -788.22296
[52,] -1013.66556 -932.40130
[53,] -710.52975 -1013.66556
[54,] -767.49005 -710.52975
[55,] -472.34564 -767.49005
[56,] -742.21833 -472.34564
[57,] -980.47502 -742.21833
[58,] -404.68662 -980.47502
[59,] -681.48195 -404.68662
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 750.15392 745.04574
2 963.86439 750.15392
3 582.19936 963.86439
4 816.48936 582.19936
5 1012.62879 816.48936
6 1094.40548 1012.62879
7 1114.51382 1094.40548
8 756.55004 1114.51382
9 767.13926 756.55004
10 1437.70793 767.13926
11 776.49103 1437.70793
12 940.61529 776.49103
13 422.99539 940.61529
14 297.12917 422.99539
15 546.37281 297.12917
16 893.99815 546.37281
17 971.21266 893.99815
18 815.36583 971.21266
19 183.90292 815.36583
20 52.77732 183.90292
21 595.07786 52.77732
22 240.39309 595.07786
23 330.72760 240.39309
24 -202.45303 330.72760
25 168.93765 -202.45303
26 79.23008 168.93765
27 -188.39578 79.23008
28 -232.38703 -188.39578
29 -102.26942 -232.38703
30 -184.50803 -102.26942
31 -598.79457 -184.50803
32 -399.77674 -598.79457
33 -164.17975 -399.77674
34 -253.97265 -164.17975
35 71.91311 -253.97265
36 100.58446 71.91311
37 126.05337 100.58446
38 -88.59331 126.05337
39 -364.30342 -88.59331
40 54.49586 -364.30342
41 -479.26238 54.49586
42 -620.59150 -479.26238
43 -854.58788 -620.59150
44 -812.32788 -854.58788
45 -906.63566 -812.32788
46 -998.83478 -906.63566
47 -851.51187 -998.83478
48 -805.19840 -851.51187
49 -1106.87042 -805.19840
50 -788.22296 -1106.87042
51 -932.40130 -788.22296
52 -1013.66556 -932.40130
53 -710.52975 -1013.66556
54 -767.49005 -710.52975
55 -472.34564 -767.49005
56 -742.21833 -472.34564
57 -980.47502 -742.21833
58 -404.68662 -980.47502
59 -681.48195 -404.68662
> 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/wessaorg/rcomp/tmp/7ttse1324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/tmp/80syr1324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/tmp/9fx731324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/tmp/105mq91324479296.ps",horizontal=F,onefile=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/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/wessaorg/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/wessaorg/rcomp/tmp/11z7mf1324479296.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/wessaorg/rcomp/tmp/1218gn1324479296.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/wessaorg/rcomp/tmp/13v7kh1324479296.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/wessaorg/rcomp/tmp/141lzv1324479296.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/wessaorg/rcomp/tmp/15fs9g1324479296.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/wessaorg/rcomp/tmp/16k8701324479296.tab")
+ }
>
> try(system("convert tmp/1qnyf1324479296.ps tmp/1qnyf1324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/2o8fm1324479296.ps tmp/2o8fm1324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/33cgc1324479296.ps tmp/33cgc1324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/4gtm01324479296.ps tmp/4gtm01324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/5lnrx1324479296.ps tmp/5lnrx1324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/6lt0i1324479296.ps tmp/6lt0i1324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/7ttse1324479296.ps tmp/7ttse1324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/80syr1324479296.ps tmp/80syr1324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/9fx731324479296.ps tmp/9fx731324479296.png",intern=TRUE))
character(0)
> try(system("convert tmp/105mq91324479296.ps tmp/105mq91324479296.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.433 0.603 4.045