R version 2.8.0 (2008-10-20)
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.
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(147768,0,137507,0,136919,0,136151,0,133001,0,125554,0,119647,0,114158,0,116193,0,152803,0,161761,0,160942,0,149470,0,139208,0,134588,0,130322,0,126611,0,122401,0,117352,0,112135,0,112879,0,148729,0,157230,0,157221,0,146681,0,136524,0,132111,0,125326,1,122716,1,116615,1,113719,1,110737,1,112093,1,143565,1,149946,1,149147,1,134339,1,122683,1,115614,1,116566,1,111272,1,104609,1,101802,1,94542,1,93051,1,124129,1,130374,1,123946,1,114971,1,105531,1,104919,1,104782,0,101281,0,94545,0,93248,0,84031,0,87486,0,115867,0,120327,0,117008,0,108811,0),dim=c(2,61),dimnames=list(c('jonger_dan_25','plan'),1:61))
> y <- array(NA,dim=c(2,61),dimnames=list(c('jonger_dan_25','plan'),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'
> #'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
jonger_dan_25 plan M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 147768 0 1 0 0 0 0 0 0 0 0 0 0
2 137507 0 0 1 0 0 0 0 0 0 0 0 0
3 136919 0 0 0 1 0 0 0 0 0 0 0 0
4 136151 0 0 0 0 1 0 0 0 0 0 0 0
5 133001 0 0 0 0 0 1 0 0 0 0 0 0
6 125554 0 0 0 0 0 0 1 0 0 0 0 0
7 119647 0 0 0 0 0 0 0 1 0 0 0 0
8 114158 0 0 0 0 0 0 0 0 1 0 0 0
9 116193 0 0 0 0 0 0 0 0 0 1 0 0
10 152803 0 0 0 0 0 0 0 0 0 0 1 0
11 161761 0 0 0 0 0 0 0 0 0 0 0 1
12 160942 0 0 0 0 0 0 0 0 0 0 0 0
13 149470 0 1 0 0 0 0 0 0 0 0 0 0
14 139208 0 0 1 0 0 0 0 0 0 0 0 0
15 134588 0 0 0 1 0 0 0 0 0 0 0 0
16 130322 0 0 0 0 1 0 0 0 0 0 0 0
17 126611 0 0 0 0 0 1 0 0 0 0 0 0
18 122401 0 0 0 0 0 0 1 0 0 0 0 0
19 117352 0 0 0 0 0 0 0 1 0 0 0 0
20 112135 0 0 0 0 0 0 0 0 1 0 0 0
21 112879 0 0 0 0 0 0 0 0 0 1 0 0
22 148729 0 0 0 0 0 0 0 0 0 0 1 0
23 157230 0 0 0 0 0 0 0 0 0 0 0 1
24 157221 0 0 0 0 0 0 0 0 0 0 0 0
25 146681 0 1 0 0 0 0 0 0 0 0 0 0
26 136524 0 0 1 0 0 0 0 0 0 0 0 0
27 132111 0 0 0 1 0 0 0 0 0 0 0 0
28 125326 1 0 0 0 1 0 0 0 0 0 0 0
29 122716 1 0 0 0 0 1 0 0 0 0 0 0
30 116615 1 0 0 0 0 0 1 0 0 0 0 0
31 113719 1 0 0 0 0 0 0 1 0 0 0 0
32 110737 1 0 0 0 0 0 0 0 1 0 0 0
33 112093 1 0 0 0 0 0 0 0 0 1 0 0
34 143565 1 0 0 0 0 0 0 0 0 0 1 0
35 149946 1 0 0 0 0 0 0 0 0 0 0 1
36 149147 1 0 0 0 0 0 0 0 0 0 0 0
37 134339 1 1 0 0 0 0 0 0 0 0 0 0
38 122683 1 0 1 0 0 0 0 0 0 0 0 0
39 115614 1 0 0 1 0 0 0 0 0 0 0 0
40 116566 1 0 0 0 1 0 0 0 0 0 0 0
41 111272 1 0 0 0 0 1 0 0 0 0 0 0
42 104609 1 0 0 0 0 0 1 0 0 0 0 0
43 101802 1 0 0 0 0 0 0 1 0 0 0 0
44 94542 1 0 0 0 0 0 0 0 1 0 0 0
45 93051 1 0 0 0 0 0 0 0 0 1 0 0
46 124129 1 0 0 0 0 0 0 0 0 0 1 0
47 130374 1 0 0 0 0 0 0 0 0 0 0 1
48 123946 1 0 0 0 0 0 0 0 0 0 0 0
49 114971 1 1 0 0 0 0 0 0 0 0 0 0
50 105531 1 0 1 0 0 0 0 0 0 0 0 0
51 104919 1 0 0 1 0 0 0 0 0 0 0 0
52 104782 0 0 0 0 1 0 0 0 0 0 0 0
53 101281 0 0 0 0 0 1 0 0 0 0 0 0
54 94545 0 0 0 0 0 0 1 0 0 0 0 0
55 93248 0 0 0 0 0 0 0 1 0 0 0 0
56 84031 0 0 0 0 0 0 0 0 1 0 0 0
57 87486 0 0 0 0 0 0 0 0 0 1 0 0
58 115867 0 0 0 0 0 0 0 0 0 0 1 0
59 120327 0 0 0 0 0 0 0 0 0 0 0 1
60 117008 0 0 0 0 0 0 0 0 0 0 0 0
61 108811 0 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) plan M1 M2 M3 M4
144914 -8154 -8523 -13362 -16823 -19023
M5 M6 M7 M8 M9 M10
-22677 -28908 -32499 -38532 -37312 -4634
M11
2275
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-27906 -7997 5753 9548 16028
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 144914 6606 21.936 < 2e-16 ***
plan -8154 3772 -2.161 0.035677 *
M1 -8523 8712 -0.978 0.332828
M2 -13362 9096 -1.469 0.148336
M3 -16823 9096 -1.850 0.070545 .
M4 -19023 9096 -2.091 0.041801 *
M5 -22677 9096 -2.493 0.016164 *
M6 -28908 9096 -3.178 0.002593 **
M7 -32499 9096 -3.573 0.000816 ***
M8 -38532 9096 -4.236 0.000102 ***
M9 -37312 9096 -4.102 0.000158 ***
M10 -4634 9096 -0.509 0.612739
M11 2275 9096 0.250 0.803579
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 14380 on 48 degrees of freedom
Multiple R-squared: 0.5469, Adjusted R-squared: 0.4337
F-statistic: 4.829 on 12 and 48 DF, p-value: 3.868e-05
> 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,] 8.741158e-03 1.748232e-02 0.991258842
[2,] 4.710099e-03 9.420198e-03 0.995289901
[3,] 1.216384e-03 2.432769e-03 0.998783616
[4,] 2.728375e-04 5.456749e-04 0.999727163
[5,] 5.947328e-05 1.189466e-04 0.999940527
[6,] 1.636503e-05 3.273006e-05 0.999983635
[7,] 6.472210e-06 1.294442e-05 0.999993528
[8,] 3.603237e-06 7.206475e-06 0.999996397
[9,] 2.427283e-06 4.854566e-06 0.999997573
[10,] 1.922986e-06 3.845971e-06 0.999998077
[11,] 2.256184e-06 4.512367e-06 0.999997744
[12,] 1.430612e-05 2.861225e-05 0.999985694
[13,] 4.265548e-06 8.531097e-06 0.999995734
[14,] 1.352252e-06 2.704504e-06 0.999998648
[15,] 4.323423e-07 8.646846e-07 0.999999568
[16,] 1.752382e-07 3.504764e-07 0.999999825
[17,] 1.961536e-07 3.923072e-07 0.999999804
[18,] 2.084436e-07 4.168872e-07 0.999999792
[19,] 2.785284e-07 5.570567e-07 0.999999721
[20,] 1.369053e-06 2.738105e-06 0.999998631
[21,] 1.181254e-04 2.362507e-04 0.999881875
[22,] 1.578772e-02 3.157544e-02 0.984212282
[23,] 5.502431e-01 8.995138e-01 0.449756910
[24,] 9.938283e-01 1.234341e-02 0.006171705
[25,] 9.966475e-01 6.705095e-03 0.003352548
[26,] 9.960203e-01 7.959426e-03 0.003979713
[27,] 9.949977e-01 1.000462e-02 0.005002311
[28,] 9.872051e-01 2.558973e-02 0.012794864
[29,] 9.873293e-01 2.534134e-02 0.012670672
[30,] 9.787491e-01 4.250186e-02 0.021250928
> postscript(file="/var/www/html/freestat/rcomp/tmp/1a6tc1229961638.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/freestat/rcomp/tmp/20xy91229961638.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/freestat/rcomp/tmp/3hoyz1229961638.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/freestat/rcomp/tmp/4azt91229961638.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/freestat/rcomp/tmp/59gld1229961638.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 7
11376.679 5954.815 8827.215 10260.015 10763.215 9547.615 7231.815
8 9 10 11 12 13 14
7775.815 8591.015 12522.815 14571.815 16027.615 13078.679 7655.815
15 16 17 18 19 20 21
6496.215 4431.015 4373.215 6394.615 4936.815 5752.815 5277.015
22 23 24 25 26 27 28
8448.815 10040.815 12306.615 10289.679 4971.815 4019.215 7588.978
29 30 31 32 33 34 35
8632.178 8762.578 9457.778 12508.778 12644.978 11438.778 10910.778
36 37 38 39 40 41 42
12386.578 6101.642 -715.222 -4323.822 -1171.022 -2811.822 -3243.422
43 44 45 46 47 48 49
-2459.222 -3686.222 -6397.022 -7997.222 -8661.222 -12814.422 -13266.358
50 51 52 53 54 55 56
-17867.222 -15018.822 -21108.985 -20956.785 -21461.385 -19167.185 -22351.185
57 58 59 60 61
-20115.985 -24413.185 -26862.185 -27906.385 -27580.321
> postscript(file="/var/www/html/freestat/rcomp/tmp/6qp921229961638.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 11376.679 NA
1 5954.815 11376.679
2 8827.215 5954.815
3 10260.015 8827.215
4 10763.215 10260.015
5 9547.615 10763.215
6 7231.815 9547.615
7 7775.815 7231.815
8 8591.015 7775.815
9 12522.815 8591.015
10 14571.815 12522.815
11 16027.615 14571.815
12 13078.679 16027.615
13 7655.815 13078.679
14 6496.215 7655.815
15 4431.015 6496.215
16 4373.215 4431.015
17 6394.615 4373.215
18 4936.815 6394.615
19 5752.815 4936.815
20 5277.015 5752.815
21 8448.815 5277.015
22 10040.815 8448.815
23 12306.615 10040.815
24 10289.679 12306.615
25 4971.815 10289.679
26 4019.215 4971.815
27 7588.978 4019.215
28 8632.178 7588.978
29 8762.578 8632.178
30 9457.778 8762.578
31 12508.778 9457.778
32 12644.978 12508.778
33 11438.778 12644.978
34 10910.778 11438.778
35 12386.578 10910.778
36 6101.642 12386.578
37 -715.222 6101.642
38 -4323.822 -715.222
39 -1171.022 -4323.822
40 -2811.822 -1171.022
41 -3243.422 -2811.822
42 -2459.222 -3243.422
43 -3686.222 -2459.222
44 -6397.022 -3686.222
45 -7997.222 -6397.022
46 -8661.222 -7997.222
47 -12814.422 -8661.222
48 -13266.358 -12814.422
49 -17867.222 -13266.358
50 -15018.822 -17867.222
51 -21108.985 -15018.822
52 -20956.785 -21108.985
53 -21461.385 -20956.785
54 -19167.185 -21461.385
55 -22351.185 -19167.185
56 -20115.985 -22351.185
57 -24413.185 -20115.985
58 -26862.185 -24413.185
59 -27906.385 -26862.185
60 -27580.321 -27906.385
61 NA -27580.321
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 5954.815 11376.679
[2,] 8827.215 5954.815
[3,] 10260.015 8827.215
[4,] 10763.215 10260.015
[5,] 9547.615 10763.215
[6,] 7231.815 9547.615
[7,] 7775.815 7231.815
[8,] 8591.015 7775.815
[9,] 12522.815 8591.015
[10,] 14571.815 12522.815
[11,] 16027.615 14571.815
[12,] 13078.679 16027.615
[13,] 7655.815 13078.679
[14,] 6496.215 7655.815
[15,] 4431.015 6496.215
[16,] 4373.215 4431.015
[17,] 6394.615 4373.215
[18,] 4936.815 6394.615
[19,] 5752.815 4936.815
[20,] 5277.015 5752.815
[21,] 8448.815 5277.015
[22,] 10040.815 8448.815
[23,] 12306.615 10040.815
[24,] 10289.679 12306.615
[25,] 4971.815 10289.679
[26,] 4019.215 4971.815
[27,] 7588.978 4019.215
[28,] 8632.178 7588.978
[29,] 8762.578 8632.178
[30,] 9457.778 8762.578
[31,] 12508.778 9457.778
[32,] 12644.978 12508.778
[33,] 11438.778 12644.978
[34,] 10910.778 11438.778
[35,] 12386.578 10910.778
[36,] 6101.642 12386.578
[37,] -715.222 6101.642
[38,] -4323.822 -715.222
[39,] -1171.022 -4323.822
[40,] -2811.822 -1171.022
[41,] -3243.422 -2811.822
[42,] -2459.222 -3243.422
[43,] -3686.222 -2459.222
[44,] -6397.022 -3686.222
[45,] -7997.222 -6397.022
[46,] -8661.222 -7997.222
[47,] -12814.422 -8661.222
[48,] -13266.358 -12814.422
[49,] -17867.222 -13266.358
[50,] -15018.822 -17867.222
[51,] -21108.985 -15018.822
[52,] -20956.785 -21108.985
[53,] -21461.385 -20956.785
[54,] -19167.185 -21461.385
[55,] -22351.185 -19167.185
[56,] -20115.985 -22351.185
[57,] -24413.185 -20115.985
[58,] -26862.185 -24413.185
[59,] -27906.385 -26862.185
[60,] -27580.321 -27906.385
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 5954.815 11376.679
2 8827.215 5954.815
3 10260.015 8827.215
4 10763.215 10260.015
5 9547.615 10763.215
6 7231.815 9547.615
7 7775.815 7231.815
8 8591.015 7775.815
9 12522.815 8591.015
10 14571.815 12522.815
11 16027.615 14571.815
12 13078.679 16027.615
13 7655.815 13078.679
14 6496.215 7655.815
15 4431.015 6496.215
16 4373.215 4431.015
17 6394.615 4373.215
18 4936.815 6394.615
19 5752.815 4936.815
20 5277.015 5752.815
21 8448.815 5277.015
22 10040.815 8448.815
23 12306.615 10040.815
24 10289.679 12306.615
25 4971.815 10289.679
26 4019.215 4971.815
27 7588.978 4019.215
28 8632.178 7588.978
29 8762.578 8632.178
30 9457.778 8762.578
31 12508.778 9457.778
32 12644.978 12508.778
33 11438.778 12644.978
34 10910.778 11438.778
35 12386.578 10910.778
36 6101.642 12386.578
37 -715.222 6101.642
38 -4323.822 -715.222
39 -1171.022 -4323.822
40 -2811.822 -1171.022
41 -3243.422 -2811.822
42 -2459.222 -3243.422
43 -3686.222 -2459.222
44 -6397.022 -3686.222
45 -7997.222 -6397.022
46 -8661.222 -7997.222
47 -12814.422 -8661.222
48 -13266.358 -12814.422
49 -17867.222 -13266.358
50 -15018.822 -17867.222
51 -21108.985 -15018.822
52 -20956.785 -21108.985
53 -21461.385 -20956.785
54 -19167.185 -21461.385
55 -22351.185 -19167.185
56 -20115.985 -22351.185
57 -24413.185 -20115.985
58 -26862.185 -24413.185
59 -27906.385 -26862.185
60 -27580.321 -27906.385
> 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/freestat/rcomp/tmp/7m6391229961638.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/freestat/rcomp/tmp/86krh1229961638.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/freestat/rcomp/tmp/9zgri1229961638.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/freestat/rcomp/tmp/10wzim1229961638.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/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/html/freestat/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/freestat/rcomp/tmp/11lydx1229961638.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/freestat/rcomp/tmp/12k0ka1229961638.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/freestat/rcomp/tmp/13uk251229961638.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/freestat/rcomp/tmp/14rmip1229961639.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/freestat/rcomp/tmp/15swr21229961639.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/freestat/rcomp/tmp/161a8q1229961639.tab")
+ }
> system("convert tmp/1a6tc1229961638.ps tmp/1a6tc1229961638.png")
> system("convert tmp/20xy91229961638.ps tmp/20xy91229961638.png")
> system("convert tmp/3hoyz1229961638.ps tmp/3hoyz1229961638.png")
> system("convert tmp/4azt91229961638.ps tmp/4azt91229961638.png")
> system("convert tmp/59gld1229961638.ps tmp/59gld1229961638.png")
> system("convert tmp/6qp921229961638.ps tmp/6qp921229961638.png")
> system("convert tmp/7m6391229961638.ps tmp/7m6391229961638.png")
> system("convert tmp/86krh1229961638.ps tmp/86krh1229961638.png")
> system("convert tmp/9zgri1229961638.ps tmp/9zgri1229961638.png")
> system("convert tmp/10wzim1229961638.ps tmp/10wzim1229961638.png")
>
>
> proc.time()
user system elapsed
3.699 2.534 4.141