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(299.63,0,305.945,0,382.252,0,348.846,0,335.367,0,373.617,0,312.612,0,312.232,0,337.161,0,331.476,0,350.103,0,345.127,0,297.256,0,295.979,0,361.007,0,321.803,0,354.937,0,349.432,0,290.979,0,349.576,0,327.625,0,349.377,0,336.777,0,339.134,0,323.321,0,318.86,0,373.583,0,333.03,0,408.556,0,414.646,0,291.514,0,348.857,0,349.368,0,375.765,0,364.136,0,349.53,0,348.167,1,332.856,1,360.551,1,346.969,1,392.815,1,372.02,1,371.027,1,342.672,1,367.343,1,390.786,1,343.785,1,362.6,1,349.468,1,340.624,1,369.536,1,407.782,1,392.239,1,404.824,1,373.669,1,344.902,1,396.7,1,398.911,1,366.009,1,392.484,1),dim=c(2,60),dimnames=list(c('x','y'),1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('x','y'),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
x y M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 299.630 0 1 0 0 0 0 0 0 0 0 0 0
2 305.945 0 0 1 0 0 0 0 0 0 0 0 0
3 382.252 0 0 0 1 0 0 0 0 0 0 0 0
4 348.846 0 0 0 0 1 0 0 0 0 0 0 0
5 335.367 0 0 0 0 0 1 0 0 0 0 0 0
6 373.617 0 0 0 0 0 0 1 0 0 0 0 0
7 312.612 0 0 0 0 0 0 0 1 0 0 0 0
8 312.232 0 0 0 0 0 0 0 0 1 0 0 0
9 337.161 0 0 0 0 0 0 0 0 0 1 0 0
10 331.476 0 0 0 0 0 0 0 0 0 0 1 0
11 350.103 0 0 0 0 0 0 0 0 0 0 0 1
12 345.127 0 0 0 0 0 0 0 0 0 0 0 0
13 297.256 0 1 0 0 0 0 0 0 0 0 0 0
14 295.979 0 0 1 0 0 0 0 0 0 0 0 0
15 361.007 0 0 0 1 0 0 0 0 0 0 0 0
16 321.803 0 0 0 0 1 0 0 0 0 0 0 0
17 354.937 0 0 0 0 0 1 0 0 0 0 0 0
18 349.432 0 0 0 0 0 0 1 0 0 0 0 0
19 290.979 0 0 0 0 0 0 0 1 0 0 0 0
20 349.576 0 0 0 0 0 0 0 0 1 0 0 0
21 327.625 0 0 0 0 0 0 0 0 0 1 0 0
22 349.377 0 0 0 0 0 0 0 0 0 0 1 0
23 336.777 0 0 0 0 0 0 0 0 0 0 0 1
24 339.134 0 0 0 0 0 0 0 0 0 0 0 0
25 323.321 0 1 0 0 0 0 0 0 0 0 0 0
26 318.860 0 0 1 0 0 0 0 0 0 0 0 0
27 373.583 0 0 0 1 0 0 0 0 0 0 0 0
28 333.030 0 0 0 0 1 0 0 0 0 0 0 0
29 408.556 0 0 0 0 0 1 0 0 0 0 0 0
30 414.646 0 0 0 0 0 0 1 0 0 0 0 0
31 291.514 0 0 0 0 0 0 0 1 0 0 0 0
32 348.857 0 0 0 0 0 0 0 0 1 0 0 0
33 349.368 0 0 0 0 0 0 0 0 0 1 0 0
34 375.765 0 0 0 0 0 0 0 0 0 0 1 0
35 364.136 0 0 0 0 0 0 0 0 0 0 0 1
36 349.530 0 0 0 0 0 0 0 0 0 0 0 0
37 348.167 1 1 0 0 0 0 0 0 0 0 0 0
38 332.856 1 0 1 0 0 0 0 0 0 0 0 0
39 360.551 1 0 0 1 0 0 0 0 0 0 0 0
40 346.969 1 0 0 0 1 0 0 0 0 0 0 0
41 392.815 1 0 0 0 0 1 0 0 0 0 0 0
42 372.020 1 0 0 0 0 0 1 0 0 0 0 0
43 371.027 1 0 0 0 0 0 0 1 0 0 0 0
44 342.672 1 0 0 0 0 0 0 0 1 0 0 0
45 367.343 1 0 0 0 0 0 0 0 0 1 0 0
46 390.786 1 0 0 0 0 0 0 0 0 0 1 0
47 343.785 1 0 0 0 0 0 0 0 0 0 0 1
48 362.600 1 0 0 0 0 0 0 0 0 0 0 0
49 349.468 1 1 0 0 0 0 0 0 0 0 0 0
50 340.624 1 0 1 0 0 0 0 0 0 0 0 0
51 369.536 1 0 0 1 0 0 0 0 0 0 0 0
52 407.782 1 0 0 0 1 0 0 0 0 0 0 0
53 392.239 1 0 0 0 0 1 0 0 0 0 0 0
54 404.824 1 0 0 0 0 0 1 0 0 0 0 0
55 373.669 1 0 0 0 0 0 0 1 0 0 0 0
56 344.902 1 0 0 0 0 0 0 0 1 0 0 0
57 396.700 1 0 0 0 0 0 0 0 0 1 0 0
58 398.911 1 0 0 0 0 0 0 0 0 0 1 0
59 366.009 1 0 0 0 0 0 0 0 0 0 0 1
60 392.484 1 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) y M1 M2 M3 M4
346.178 28.991 -34.207 -38.922 11.611 -6.089
M5 M6 M7 M8 M9 M10
19.008 25.133 -29.815 -18.127 -2.136 11.488
M11
-5.613
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-29.819 -13.020 -1.651 11.413 43.370
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 346.178 9.371 36.940 < 2e-16 ***
y 28.991 5.375 5.394 2.20e-06 ***
M1 -34.207 12.900 -2.652 0.01088 *
M2 -38.922 12.900 -3.017 0.00411 **
M3 11.611 12.900 0.900 0.37266
M4 -6.089 12.900 -0.472 0.63909
M5 19.008 12.900 1.474 0.14728
M6 25.133 12.900 1.948 0.05736 .
M7 -29.815 12.900 -2.311 0.02525 *
M8 -18.127 12.900 -1.405 0.16652
M9 -2.136 12.900 -0.166 0.86922
M10 11.488 12.900 0.891 0.37770
M11 -5.613 12.900 -0.435 0.66546
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 20.4 on 47 degrees of freedom
Multiple R-squared: 0.6504, Adjusted R-squared: 0.5611
F-statistic: 7.286 on 12 and 47 DF, p-value: 2.650e-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.28192400 0.56384801 0.7180760
[2,] 0.22256854 0.44513709 0.7774315
[3,] 0.21855390 0.43710781 0.7814461
[4,] 0.21518978 0.43037956 0.7848102
[5,] 0.30288726 0.60577452 0.6971127
[6,] 0.23975602 0.47951203 0.7602440
[7,] 0.20845773 0.41691545 0.7915423
[8,] 0.15174701 0.30349402 0.8482530
[9,] 0.10403277 0.20806554 0.8959672
[10,] 0.10687998 0.21375996 0.8931200
[11,] 0.08405020 0.16810039 0.9159498
[12,] 0.05789570 0.11579141 0.9421043
[13,] 0.04290241 0.08580482 0.9570976
[14,] 0.27690172 0.55380344 0.7230983
[15,] 0.54544183 0.90911634 0.4545582
[16,] 0.78024135 0.43951730 0.2197586
[17,] 0.75757818 0.48484364 0.2424218
[18,] 0.72075490 0.55849021 0.2792451
[19,] 0.70470261 0.59059477 0.2952974
[20,] 0.71339975 0.57320050 0.2866002
[21,] 0.61903784 0.76192433 0.3809622
[22,] 0.51217216 0.97565569 0.4878278
[23,] 0.41464103 0.82928205 0.5853590
[24,] 0.38356558 0.76713117 0.6164344
[25,] 0.72312834 0.55374333 0.2768717
[26,] 0.60434034 0.79131933 0.3956597
[27,] 0.68295734 0.63408533 0.3170427
[28,] 0.61811597 0.76376805 0.3818840
[29,] 0.44875152 0.89750305 0.5512485
> postscript(file="/var/www/html/rcomp/tmp/1k2kh1227449466.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/2lcbd1227449466.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/3ocy91227449466.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/41xw81227449466.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/55tr81227449466.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 7
-12.341817 -1.311217 24.462783 8.756583 -29.819217 2.305783 -3.751617
8 9 10 11 12 13 14
-15.819217 -6.881817 -26.190417 9.537583 -1.051417 -14.715817 -11.277217
15 16 17 18 19 20 21
3.217783 -18.286417 -10.249217 -21.879217 -25.384617 21.524783 -16.417817
22 23 24 25 26 27 28
-8.289417 -3.788417 -7.044417 11.349183 11.603783 15.793783 -7.059417
29 30 31 32 33 34 35
43.369783 43.334783 -24.849617 20.805783 5.325183 18.098583 23.570583
36 37 38 39 40 41 42
3.351583 7.203725 -3.391675 -26.229675 -22.111875 -1.362675 -28.282675
43 44 45 46 47 48 49
25.671925 -14.370675 -5.691275 4.128125 -25.771875 -12.569875 8.504725
50 51 52 53 54 55 56
4.376325 -17.244675 38.701125 -1.938675 4.521325 28.313925 -12.140675
57 58 59 60
23.665725 12.253125 -3.547875 17.314125
> postscript(file="/var/www/html/rcomp/tmp/6fagt1227449466.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 -12.341817 NA
1 -1.311217 -12.341817
2 24.462783 -1.311217
3 8.756583 24.462783
4 -29.819217 8.756583
5 2.305783 -29.819217
6 -3.751617 2.305783
7 -15.819217 -3.751617
8 -6.881817 -15.819217
9 -26.190417 -6.881817
10 9.537583 -26.190417
11 -1.051417 9.537583
12 -14.715817 -1.051417
13 -11.277217 -14.715817
14 3.217783 -11.277217
15 -18.286417 3.217783
16 -10.249217 -18.286417
17 -21.879217 -10.249217
18 -25.384617 -21.879217
19 21.524783 -25.384617
20 -16.417817 21.524783
21 -8.289417 -16.417817
22 -3.788417 -8.289417
23 -7.044417 -3.788417
24 11.349183 -7.044417
25 11.603783 11.349183
26 15.793783 11.603783
27 -7.059417 15.793783
28 43.369783 -7.059417
29 43.334783 43.369783
30 -24.849617 43.334783
31 20.805783 -24.849617
32 5.325183 20.805783
33 18.098583 5.325183
34 23.570583 18.098583
35 3.351583 23.570583
36 7.203725 3.351583
37 -3.391675 7.203725
38 -26.229675 -3.391675
39 -22.111875 -26.229675
40 -1.362675 -22.111875
41 -28.282675 -1.362675
42 25.671925 -28.282675
43 -14.370675 25.671925
44 -5.691275 -14.370675
45 4.128125 -5.691275
46 -25.771875 4.128125
47 -12.569875 -25.771875
48 8.504725 -12.569875
49 4.376325 8.504725
50 -17.244675 4.376325
51 38.701125 -17.244675
52 -1.938675 38.701125
53 4.521325 -1.938675
54 28.313925 4.521325
55 -12.140675 28.313925
56 23.665725 -12.140675
57 12.253125 23.665725
58 -3.547875 12.253125
59 17.314125 -3.547875
60 NA 17.314125
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -1.311217 -12.341817
[2,] 24.462783 -1.311217
[3,] 8.756583 24.462783
[4,] -29.819217 8.756583
[5,] 2.305783 -29.819217
[6,] -3.751617 2.305783
[7,] -15.819217 -3.751617
[8,] -6.881817 -15.819217
[9,] -26.190417 -6.881817
[10,] 9.537583 -26.190417
[11,] -1.051417 9.537583
[12,] -14.715817 -1.051417
[13,] -11.277217 -14.715817
[14,] 3.217783 -11.277217
[15,] -18.286417 3.217783
[16,] -10.249217 -18.286417
[17,] -21.879217 -10.249217
[18,] -25.384617 -21.879217
[19,] 21.524783 -25.384617
[20,] -16.417817 21.524783
[21,] -8.289417 -16.417817
[22,] -3.788417 -8.289417
[23,] -7.044417 -3.788417
[24,] 11.349183 -7.044417
[25,] 11.603783 11.349183
[26,] 15.793783 11.603783
[27,] -7.059417 15.793783
[28,] 43.369783 -7.059417
[29,] 43.334783 43.369783
[30,] -24.849617 43.334783
[31,] 20.805783 -24.849617
[32,] 5.325183 20.805783
[33,] 18.098583 5.325183
[34,] 23.570583 18.098583
[35,] 3.351583 23.570583
[36,] 7.203725 3.351583
[37,] -3.391675 7.203725
[38,] -26.229675 -3.391675
[39,] -22.111875 -26.229675
[40,] -1.362675 -22.111875
[41,] -28.282675 -1.362675
[42,] 25.671925 -28.282675
[43,] -14.370675 25.671925
[44,] -5.691275 -14.370675
[45,] 4.128125 -5.691275
[46,] -25.771875 4.128125
[47,] -12.569875 -25.771875
[48,] 8.504725 -12.569875
[49,] 4.376325 8.504725
[50,] -17.244675 4.376325
[51,] 38.701125 -17.244675
[52,] -1.938675 38.701125
[53,] 4.521325 -1.938675
[54,] 28.313925 4.521325
[55,] -12.140675 28.313925
[56,] 23.665725 -12.140675
[57,] 12.253125 23.665725
[58,] -3.547875 12.253125
[59,] 17.314125 -3.547875
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -1.311217 -12.341817
2 24.462783 -1.311217
3 8.756583 24.462783
4 -29.819217 8.756583
5 2.305783 -29.819217
6 -3.751617 2.305783
7 -15.819217 -3.751617
8 -6.881817 -15.819217
9 -26.190417 -6.881817
10 9.537583 -26.190417
11 -1.051417 9.537583
12 -14.715817 -1.051417
13 -11.277217 -14.715817
14 3.217783 -11.277217
15 -18.286417 3.217783
16 -10.249217 -18.286417
17 -21.879217 -10.249217
18 -25.384617 -21.879217
19 21.524783 -25.384617
20 -16.417817 21.524783
21 -8.289417 -16.417817
22 -3.788417 -8.289417
23 -7.044417 -3.788417
24 11.349183 -7.044417
25 11.603783 11.349183
26 15.793783 11.603783
27 -7.059417 15.793783
28 43.369783 -7.059417
29 43.334783 43.369783
30 -24.849617 43.334783
31 20.805783 -24.849617
32 5.325183 20.805783
33 18.098583 5.325183
34 23.570583 18.098583
35 3.351583 23.570583
36 7.203725 3.351583
37 -3.391675 7.203725
38 -26.229675 -3.391675
39 -22.111875 -26.229675
40 -1.362675 -22.111875
41 -28.282675 -1.362675
42 25.671925 -28.282675
43 -14.370675 25.671925
44 -5.691275 -14.370675
45 4.128125 -5.691275
46 -25.771875 4.128125
47 -12.569875 -25.771875
48 8.504725 -12.569875
49 4.376325 8.504725
50 -17.244675 4.376325
51 38.701125 -17.244675
52 -1.938675 38.701125
53 4.521325 -1.938675
54 28.313925 4.521325
55 -12.140675 28.313925
56 23.665725 -12.140675
57 12.253125 23.665725
58 -3.547875 12.253125
59 17.314125 -3.547875
> 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/7aaka1227449466.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/8pzmw1227449466.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/93rpt1227449466.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/10pef71227449466.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/115bpz1227449466.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/12e8nj1227449466.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/1332d01227449466.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/140ffl1227449466.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/15sv2p1227449467.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/16ixo91227449467.tab")
+ }
>
> system("convert tmp/1k2kh1227449466.ps tmp/1k2kh1227449466.png")
> system("convert tmp/2lcbd1227449466.ps tmp/2lcbd1227449466.png")
> system("convert tmp/3ocy91227449466.ps tmp/3ocy91227449466.png")
> system("convert tmp/41xw81227449466.ps tmp/41xw81227449466.png")
> system("convert tmp/55tr81227449466.ps tmp/55tr81227449466.png")
> system("convert tmp/6fagt1227449466.ps tmp/6fagt1227449466.png")
> system("convert tmp/7aaka1227449466.ps tmp/7aaka1227449466.png")
> system("convert tmp/8pzmw1227449466.ps tmp/8pzmw1227449466.png")
> system("convert tmp/93rpt1227449466.ps tmp/93rpt1227449466.png")
> system("convert tmp/10pef71227449466.ps tmp/10pef71227449466.png")
>
>
> proc.time()
user system elapsed
2.405 1.557 2.800