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,36.68,3922,36.7,3759,36.71,4138,36.72,4634,36.73,3996,36.73,4308,36.87,4429,37.31,5219,37.39,4929,37.42,5755,37.51,5592,37.67,4163,37.67,4962,37.71,5208,37.78,4755,37.79,4491,37.84,5732,37.88,5731,38.34,5040,38.58,6102,38.72,4904,38.83,5369,38.9,5578,38.92,4619,38.94,4731,39.1,5011,39.14,5299,39.16,4146,39.32,4625,39.34,4736,39.44,4219,39.92,5116,40.19,4205,40.2,4121,40.27,5103,40.28,4300,40.3,4578,40.34,3809,40.4,5526,40.43,4247,40.48,3830,40.48,4394,40.63,4826,40.74,4409,40.77,4569,40.91,4106,40.92,4794,41.03,3914,41,3793,41.04,4405,41.33,4022,41.44,4100,41.46,4788,41.55,3163,41.55,3585,41.81,3903,41.78,4178,41.84,3863,41.84,4187,41.86),dim=c(2,60),dimnames=list(c('Y','X'),1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('Y','X'),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 = '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
Y X M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 5560 36.68 1 0 0 0 0 0 0 0 0 0 0
2 3922 36.70 0 1 0 0 0 0 0 0 0 0 0
3 3759 36.71 0 0 1 0 0 0 0 0 0 0 0
4 4138 36.72 0 0 0 1 0 0 0 0 0 0 0
5 4634 36.73 0 0 0 0 1 0 0 0 0 0 0
6 3996 36.73 0 0 0 0 0 1 0 0 0 0 0
7 4308 36.87 0 0 0 0 0 0 1 0 0 0 0
8 4429 37.31 0 0 0 0 0 0 0 1 0 0 0
9 5219 37.39 0 0 0 0 0 0 0 0 1 0 0
10 4929 37.42 0 0 0 0 0 0 0 0 0 1 0
11 5755 37.51 0 0 0 0 0 0 0 0 0 0 1
12 5592 37.67 0 0 0 0 0 0 0 0 0 0 0
13 4163 37.67 1 0 0 0 0 0 0 0 0 0 0
14 4962 37.71 0 1 0 0 0 0 0 0 0 0 0
15 5208 37.78 0 0 1 0 0 0 0 0 0 0 0
16 4755 37.79 0 0 0 1 0 0 0 0 0 0 0
17 4491 37.84 0 0 0 0 1 0 0 0 0 0 0
18 5732 37.88 0 0 0 0 0 1 0 0 0 0 0
19 5731 38.34 0 0 0 0 0 0 1 0 0 0 0
20 5040 38.58 0 0 0 0 0 0 0 1 0 0 0
21 6102 38.72 0 0 0 0 0 0 0 0 1 0 0
22 4904 38.83 0 0 0 0 0 0 0 0 0 1 0
23 5369 38.90 0 0 0 0 0 0 0 0 0 0 1
24 5578 38.92 0 0 0 0 0 0 0 0 0 0 0
25 4619 38.94 1 0 0 0 0 0 0 0 0 0 0
26 4731 39.10 0 1 0 0 0 0 0 0 0 0 0
27 5011 39.14 0 0 1 0 0 0 0 0 0 0 0
28 5299 39.16 0 0 0 1 0 0 0 0 0 0 0
29 4146 39.32 0 0 0 0 1 0 0 0 0 0 0
30 4625 39.34 0 0 0 0 0 1 0 0 0 0 0
31 4736 39.44 0 0 0 0 0 0 1 0 0 0 0
32 4219 39.92 0 0 0 0 0 0 0 1 0 0 0
33 5116 40.19 0 0 0 0 0 0 0 0 1 0 0
34 4205 40.20 0 0 0 0 0 0 0 0 0 1 0
35 4121 40.27 0 0 0 0 0 0 0 0 0 0 1
36 5103 40.28 0 0 0 0 0 0 0 0 0 0 0
37 4300 40.30 1 0 0 0 0 0 0 0 0 0 0
38 4578 40.34 0 1 0 0 0 0 0 0 0 0 0
39 3809 40.40 0 0 1 0 0 0 0 0 0 0 0
40 5526 40.43 0 0 0 1 0 0 0 0 0 0 0
41 4247 40.48 0 0 0 0 1 0 0 0 0 0 0
42 3830 40.48 0 0 0 0 0 1 0 0 0 0 0
43 4394 40.63 0 0 0 0 0 0 1 0 0 0 0
44 4826 40.74 0 0 0 0 0 0 0 1 0 0 0
45 4409 40.77 0 0 0 0 0 0 0 0 1 0 0
46 4569 40.91 0 0 0 0 0 0 0 0 0 1 0
47 4106 40.92 0 0 0 0 0 0 0 0 0 0 1
48 4794 41.03 0 0 0 0 0 0 0 0 0 0 0
49 3914 41.00 1 0 0 0 0 0 0 0 0 0 0
50 3793 41.04 0 1 0 0 0 0 0 0 0 0 0
51 4405 41.33 0 0 1 0 0 0 0 0 0 0 0
52 4022 41.44 0 0 0 1 0 0 0 0 0 0 0
53 4100 41.46 0 0 0 0 1 0 0 0 0 0 0
54 4788 41.55 0 0 0 0 0 1 0 0 0 0 0
55 3163 41.55 0 0 0 0 0 0 1 0 0 0 0
56 3585 41.81 0 0 0 0 0 0 0 1 0 0 0
57 3903 41.78 0 0 0 0 0 0 0 0 1 0 0
58 4178 41.84 0 0 0 0 0 0 0 0 0 1 0
59 3863 41.84 0 0 0 0 0 0 0 0 0 0 1
60 4187 41.86 0 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
12017.8 -174.4 -719.9 -823.4 -765.9 -450.0
M5 M6 M7 M8 M9 M10
-864.3 -588.4 -686.6 -679.8 -132.7 -513.3
M11
-419.2
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1091.2946 -374.2667 -0.2169 356.3958 1085.6822
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 12017.78 1895.38 6.341 8.2e-08 ***
X -174.38 46.98 -3.712 0.000544 ***
M1 -719.91 375.90 -1.915 0.061564 .
M2 -823.45 375.54 -2.193 0.033315 *
M3 -765.86 375.03 -2.042 0.046772 *
M4 -449.98 374.85 -1.200 0.235985
M5 -864.27 374.57 -2.307 0.025487 *
M6 -588.43 374.43 -1.572 0.122766
M7 -686.59 373.76 -1.837 0.072540 .
M8 -679.83 372.98 -1.823 0.074712 .
M9 -132.74 372.84 -0.356 0.723421
M10 -513.33 372.78 -1.377 0.175026
M11 -419.16 372.76 -1.124 0.266515
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 589.4 on 47 degrees of freedom
Multiple R-squared: 0.3161, Adjusted R-squared: 0.1415
F-statistic: 1.811 on 12 and 47 DF, p-value: 0.07387
> 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.9871828 0.02563449 0.01281725
[2,] 0.9825383 0.03492346 0.01746173
[3,] 0.9882491 0.02350190 0.01175095
[4,] 0.9854849 0.02903014 0.01451507
[5,] 0.9714872 0.05702562 0.02851281
[6,] 0.9647894 0.07042116 0.03521058
[7,] 0.9596014 0.08079719 0.04039859
[8,] 0.9693200 0.06136008 0.03068004
[9,] 0.9540595 0.09188108 0.04594054
[10,] 0.9494679 0.10106419 0.05053210
[11,] 0.9198646 0.16027083 0.08013542
[12,] 0.8827968 0.23440647 0.11720324
[13,] 0.8328670 0.33426605 0.16713303
[14,] 0.8571939 0.28561220 0.14280610
[15,] 0.8313361 0.33732789 0.16866394
[16,] 0.7955710 0.40885802 0.20442901
[17,] 0.7889645 0.42207095 0.21103548
[18,] 0.7783655 0.44326902 0.22163451
[19,] 0.7892311 0.42153773 0.21076887
[20,] 0.8235293 0.35294136 0.17647068
[21,] 0.7536150 0.49276991 0.24638496
[22,] 0.6675220 0.66495598 0.33247799
[23,] 0.5861048 0.82779038 0.41389519
[24,] 0.6425654 0.71486913 0.35743456
[25,] 0.7090381 0.58192382 0.29096191
[26,] 0.5995098 0.80098034 0.40049017
[27,] 0.9484821 0.10303576 0.05151788
[28,] 0.9520041 0.09599171 0.04799586
[29,] 0.9890431 0.02191381 0.01095690
> postscript(file="/var/www/html/rcomp/tmp/1f3nu1258198161.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/2l88e1258198161.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/3zp8b1258198161.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/4w9cr1258198161.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/5kgw31258198161.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> qqnorm(mysum$resid, main='Residual Normal Q-Q Plot')
> qqline(mysum$resid)
> grid()
> dev.off()
null device
1
> (myerror <- as.ts(mysum$resid))
Time Series:
Start = 1
End = 60
Frequency = 1
1 2 3 4 5 6
658.52897 -872.44638 -1091.29462 -1026.42860 -114.39903 -1028.23054
7 8 9 10 11 12
-593.66206 -402.69462 -145.83353 -50.00888 697.51524 143.25609
13 14 15 16 17 18
-565.83103 343.68130 544.29608 -222.83790 -63.83297 908.31087
19 20 21 22 23 24
1085.68218 429.77285 969.09697 170.87233 553.90877 347.23588
25 26 27 28 29 30
111.63644 355.07483 584.45810 560.06796 -150.74489 55.91127
31 32 33 34 35 36
282.50440 -157.55281 239.44121 -289.22182 -455.18537 109.39790
37 38 39 40 41 42
29.79846 418.31078 -397.81827 1008.53543 152.54036 -540.29115
43 44 45 46 47 48
148.02117 592.44194 -366.41616 198.59070 -356.83588 -68.81422
49 50 51 52 53 54
-234.13285 -244.62053 360.35870 -319.33689 176.43652 604.29955
55 56 57 58 59 60
-922.54570 -461.96736 -696.28849 -30.23233 -439.40275 -531.07564
> postscript(file="/var/www/html/rcomp/tmp/6uyla1258198161.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> dum <- cbind(lag(myerror,k=1),myerror)
> dum
Time Series:
Start = 0
End = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 658.52897 NA
1 -872.44638 658.52897
2 -1091.29462 -872.44638
3 -1026.42860 -1091.29462
4 -114.39903 -1026.42860
5 -1028.23054 -114.39903
6 -593.66206 -1028.23054
7 -402.69462 -593.66206
8 -145.83353 -402.69462
9 -50.00888 -145.83353
10 697.51524 -50.00888
11 143.25609 697.51524
12 -565.83103 143.25609
13 343.68130 -565.83103
14 544.29608 343.68130
15 -222.83790 544.29608
16 -63.83297 -222.83790
17 908.31087 -63.83297
18 1085.68218 908.31087
19 429.77285 1085.68218
20 969.09697 429.77285
21 170.87233 969.09697
22 553.90877 170.87233
23 347.23588 553.90877
24 111.63644 347.23588
25 355.07483 111.63644
26 584.45810 355.07483
27 560.06796 584.45810
28 -150.74489 560.06796
29 55.91127 -150.74489
30 282.50440 55.91127
31 -157.55281 282.50440
32 239.44121 -157.55281
33 -289.22182 239.44121
34 -455.18537 -289.22182
35 109.39790 -455.18537
36 29.79846 109.39790
37 418.31078 29.79846
38 -397.81827 418.31078
39 1008.53543 -397.81827
40 152.54036 1008.53543
41 -540.29115 152.54036
42 148.02117 -540.29115
43 592.44194 148.02117
44 -366.41616 592.44194
45 198.59070 -366.41616
46 -356.83588 198.59070
47 -68.81422 -356.83588
48 -234.13285 -68.81422
49 -244.62053 -234.13285
50 360.35870 -244.62053
51 -319.33689 360.35870
52 176.43652 -319.33689
53 604.29955 176.43652
54 -922.54570 604.29955
55 -461.96736 -922.54570
56 -696.28849 -461.96736
57 -30.23233 -696.28849
58 -439.40275 -30.23233
59 -531.07564 -439.40275
60 NA -531.07564
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -872.44638 658.52897
[2,] -1091.29462 -872.44638
[3,] -1026.42860 -1091.29462
[4,] -114.39903 -1026.42860
[5,] -1028.23054 -114.39903
[6,] -593.66206 -1028.23054
[7,] -402.69462 -593.66206
[8,] -145.83353 -402.69462
[9,] -50.00888 -145.83353
[10,] 697.51524 -50.00888
[11,] 143.25609 697.51524
[12,] -565.83103 143.25609
[13,] 343.68130 -565.83103
[14,] 544.29608 343.68130
[15,] -222.83790 544.29608
[16,] -63.83297 -222.83790
[17,] 908.31087 -63.83297
[18,] 1085.68218 908.31087
[19,] 429.77285 1085.68218
[20,] 969.09697 429.77285
[21,] 170.87233 969.09697
[22,] 553.90877 170.87233
[23,] 347.23588 553.90877
[24,] 111.63644 347.23588
[25,] 355.07483 111.63644
[26,] 584.45810 355.07483
[27,] 560.06796 584.45810
[28,] -150.74489 560.06796
[29,] 55.91127 -150.74489
[30,] 282.50440 55.91127
[31,] -157.55281 282.50440
[32,] 239.44121 -157.55281
[33,] -289.22182 239.44121
[34,] -455.18537 -289.22182
[35,] 109.39790 -455.18537
[36,] 29.79846 109.39790
[37,] 418.31078 29.79846
[38,] -397.81827 418.31078
[39,] 1008.53543 -397.81827
[40,] 152.54036 1008.53543
[41,] -540.29115 152.54036
[42,] 148.02117 -540.29115
[43,] 592.44194 148.02117
[44,] -366.41616 592.44194
[45,] 198.59070 -366.41616
[46,] -356.83588 198.59070
[47,] -68.81422 -356.83588
[48,] -234.13285 -68.81422
[49,] -244.62053 -234.13285
[50,] 360.35870 -244.62053
[51,] -319.33689 360.35870
[52,] 176.43652 -319.33689
[53,] 604.29955 176.43652
[54,] -922.54570 604.29955
[55,] -461.96736 -922.54570
[56,] -696.28849 -461.96736
[57,] -30.23233 -696.28849
[58,] -439.40275 -30.23233
[59,] -531.07564 -439.40275
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -872.44638 658.52897
2 -1091.29462 -872.44638
3 -1026.42860 -1091.29462
4 -114.39903 -1026.42860
5 -1028.23054 -114.39903
6 -593.66206 -1028.23054
7 -402.69462 -593.66206
8 -145.83353 -402.69462
9 -50.00888 -145.83353
10 697.51524 -50.00888
11 143.25609 697.51524
12 -565.83103 143.25609
13 343.68130 -565.83103
14 544.29608 343.68130
15 -222.83790 544.29608
16 -63.83297 -222.83790
17 908.31087 -63.83297
18 1085.68218 908.31087
19 429.77285 1085.68218
20 969.09697 429.77285
21 170.87233 969.09697
22 553.90877 170.87233
23 347.23588 553.90877
24 111.63644 347.23588
25 355.07483 111.63644
26 584.45810 355.07483
27 560.06796 584.45810
28 -150.74489 560.06796
29 55.91127 -150.74489
30 282.50440 55.91127
31 -157.55281 282.50440
32 239.44121 -157.55281
33 -289.22182 239.44121
34 -455.18537 -289.22182
35 109.39790 -455.18537
36 29.79846 109.39790
37 418.31078 29.79846
38 -397.81827 418.31078
39 1008.53543 -397.81827
40 152.54036 1008.53543
41 -540.29115 152.54036
42 148.02117 -540.29115
43 592.44194 148.02117
44 -366.41616 592.44194
45 198.59070 -366.41616
46 -356.83588 198.59070
47 -68.81422 -356.83588
48 -234.13285 -68.81422
49 -244.62053 -234.13285
50 360.35870 -244.62053
51 -319.33689 360.35870
52 176.43652 -319.33689
53 604.29955 176.43652
54 -922.54570 604.29955
55 -461.96736 -922.54570
56 -696.28849 -461.96736
57 -30.23233 -696.28849
58 -439.40275 -30.23233
59 -531.07564 -439.40275
> 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/7saau1258198161.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/80ib41258198161.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/9jz7r1258198161.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/107jm91258198161.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/11j24b1258198161.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/12flf51258198161.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/13verl1258198161.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/1401vy1258198161.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/153ej41258198161.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/16e8u61258198161.tab")
+ }
>
> system("convert tmp/1f3nu1258198161.ps tmp/1f3nu1258198161.png")
> system("convert tmp/2l88e1258198161.ps tmp/2l88e1258198161.png")
> system("convert tmp/3zp8b1258198161.ps tmp/3zp8b1258198161.png")
> system("convert tmp/4w9cr1258198161.ps tmp/4w9cr1258198161.png")
> system("convert tmp/5kgw31258198161.ps tmp/5kgw31258198161.png")
> system("convert tmp/6uyla1258198161.ps tmp/6uyla1258198161.png")
> system("convert tmp/7saau1258198161.ps tmp/7saau1258198161.png")
> system("convert tmp/80ib41258198161.ps tmp/80ib41258198161.png")
> system("convert tmp/9jz7r1258198161.ps tmp/9jz7r1258198161.png")
> system("convert tmp/107jm91258198161.ps tmp/107jm91258198161.png")
>
>
> proc.time()
user system elapsed
2.356 1.520 3.458