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(562,0,561,0,555,0,544,0,537,0,543,0,594,0,611,0,613,0,611,0,594,0,595,0,591,0,589,0,584,0,573,0,567,0,569,0,621,0,629,0,628,0,612,0,595,0,597,0,593,0,590,0,580,0,574,0,573,0,573,0,620,0,626,0,620,0,588,0,566,0,557,0,561,0,549,0,532,0,526,0,511,0,499,0,555,0,565,0,542,0,527,1,510,1,514,1,517,1,508,1,493,1,490,1,469,1,478,1,528,1,534,1,518,1,506,1,502,1,516,1,528,1),dim=c(2,61),dimnames=list(c('Y','X'),1:61))
> y <- array(NA,dim=c(2,61),dimnames=list(c('Y','X'),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
Y X M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 562 0 1 0 0 0 0 0 0 0 0 0 0
2 561 0 0 1 0 0 0 0 0 0 0 0 0
3 555 0 0 0 1 0 0 0 0 0 0 0 0
4 544 0 0 0 0 1 0 0 0 0 0 0 0
5 537 0 0 0 0 0 1 0 0 0 0 0 0
6 543 0 0 0 0 0 0 1 0 0 0 0 0
7 594 0 0 0 0 0 0 0 1 0 0 0 0
8 611 0 0 0 0 0 0 0 0 1 0 0 0
9 613 0 0 0 0 0 0 0 0 0 1 0 0
10 611 0 0 0 0 0 0 0 0 0 0 1 0
11 594 0 0 0 0 0 0 0 0 0 0 0 1
12 595 0 0 0 0 0 0 0 0 0 0 0 0
13 591 0 1 0 0 0 0 0 0 0 0 0 0
14 589 0 0 1 0 0 0 0 0 0 0 0 0
15 584 0 0 0 1 0 0 0 0 0 0 0 0
16 573 0 0 0 0 1 0 0 0 0 0 0 0
17 567 0 0 0 0 0 1 0 0 0 0 0 0
18 569 0 0 0 0 0 0 1 0 0 0 0 0
19 621 0 0 0 0 0 0 0 1 0 0 0 0
20 629 0 0 0 0 0 0 0 0 1 0 0 0
21 628 0 0 0 0 0 0 0 0 0 1 0 0
22 612 0 0 0 0 0 0 0 0 0 0 1 0
23 595 0 0 0 0 0 0 0 0 0 0 0 1
24 597 0 0 0 0 0 0 0 0 0 0 0 0
25 593 0 1 0 0 0 0 0 0 0 0 0 0
26 590 0 0 1 0 0 0 0 0 0 0 0 0
27 580 0 0 0 1 0 0 0 0 0 0 0 0
28 574 0 0 0 0 1 0 0 0 0 0 0 0
29 573 0 0 0 0 0 1 0 0 0 0 0 0
30 573 0 0 0 0 0 0 1 0 0 0 0 0
31 620 0 0 0 0 0 0 0 1 0 0 0 0
32 626 0 0 0 0 0 0 0 0 1 0 0 0
33 620 0 0 0 0 0 0 0 0 0 1 0 0
34 588 0 0 0 0 0 0 0 0 0 0 1 0
35 566 0 0 0 0 0 0 0 0 0 0 0 1
36 557 0 0 0 0 0 0 0 0 0 0 0 0
37 561 0 1 0 0 0 0 0 0 0 0 0 0
38 549 0 0 1 0 0 0 0 0 0 0 0 0
39 532 0 0 0 1 0 0 0 0 0 0 0 0
40 526 0 0 0 0 1 0 0 0 0 0 0 0
41 511 0 0 0 0 0 1 0 0 0 0 0 0
42 499 0 0 0 0 0 0 1 0 0 0 0 0
43 555 0 0 0 0 0 0 0 1 0 0 0 0
44 565 0 0 0 0 0 0 0 0 1 0 0 0
45 542 0 0 0 0 0 0 0 0 0 1 0 0
46 527 1 0 0 0 0 0 0 0 0 0 1 0
47 510 1 0 0 0 0 0 0 0 0 0 0 1
48 514 1 0 0 0 0 0 0 0 0 0 0 0
49 517 1 1 0 0 0 0 0 0 0 0 0 0
50 508 1 0 1 0 0 0 0 0 0 0 0 0
51 493 1 0 0 1 0 0 0 0 0 0 0 0
52 490 1 0 0 0 1 0 0 0 0 0 0 0
53 469 1 0 0 0 0 1 0 0 0 0 0 0
54 478 1 0 0 0 0 0 1 0 0 0 0 0
55 528 1 0 0 0 0 0 0 1 0 0 0 0
56 534 1 0 0 0 0 0 0 0 1 0 0 0
57 518 1 0 0 0 0 0 0 0 0 1 0 0
58 506 1 0 0 0 0 0 0 0 0 0 1 0
59 502 1 0 0 0 0 0 0 0 0 0 0 1
60 516 1 0 0 0 0 0 0 0 0 0 0 0
61 528 1 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) X M1 M2 M3 M4
584.372 -71.429 -1.895 -10.686 -21.286 -28.686
M5 M6 M7 M8 M9 M10
-38.686 -37.686 13.514 22.914 14.114 13.000
M11
-2.400
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-56.486 -9.372 3.058 16.314 29.514
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 584.372 10.563 55.324 < 2e-16 ***
X -71.429 6.781 -10.534 4.48e-14 ***
M1 -1.895 13.830 -0.137 0.8916
M2 -10.686 14.501 -0.737 0.4648
M3 -21.286 14.501 -1.468 0.1486
M4 -28.686 14.501 -1.978 0.0537 .
M5 -38.686 14.501 -2.668 0.0104 *
M6 -37.686 14.501 -2.599 0.0124 *
M7 13.514 14.501 0.932 0.3560
M8 22.914 14.501 1.580 0.1206
M9 14.114 14.501 0.973 0.3353
M10 13.000 14.437 0.900 0.3724
M11 -2.400 14.437 -0.166 0.8687
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 22.83 on 48 degrees of freedom
Multiple R-squared: 0.7616, Adjusted R-squared: 0.702
F-statistic: 12.78 on 12 and 48 DF, p-value: 3.679e-11
> 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.49636004 0.992720084 0.503639958
[2,] 0.43834596 0.876691914 0.561654043
[3,] 0.37286210 0.745724192 0.627137904
[4,] 0.32678770 0.653575404 0.673212298
[5,] 0.25794727 0.515894549 0.742052725
[6,] 0.21698373 0.433967465 0.783016268
[7,] 0.15114018 0.302280369 0.848859816
[8,] 0.10140623 0.202812470 0.898593765
[9,] 0.06641265 0.132825306 0.933587347
[10,] 0.04843506 0.096870121 0.951564939
[11,] 0.03778809 0.075576177 0.962211911
[12,] 0.02990658 0.059813164 0.970093418
[13,] 0.02678089 0.053561780 0.973219110
[14,] 0.04367927 0.087358532 0.956320734
[15,] 0.08029617 0.160592338 0.919703831
[16,] 0.13946640 0.278932794 0.860533603
[17,] 0.26476601 0.529532027 0.735233987
[18,] 0.89461545 0.210769108 0.105384554
[19,] 0.98367191 0.032656186 0.016328093
[20,] 0.99792014 0.004159722 0.002079861
[21,] 0.99846360 0.003072798 0.001536399
[22,] 0.99735368 0.005292637 0.002646319
[23,] 0.99682766 0.006344680 0.003172340
[24,] 0.99654403 0.006911943 0.003455971
[25,] 0.99467483 0.010650349 0.005325174
[26,] 0.99616575 0.007668509 0.003834254
[27,] 0.99441286 0.011174279 0.005587140
[28,] 0.98741450 0.025170995 0.012585498
[29,] 0.97115785 0.057684305 0.028842153
[30,] 0.93412343 0.131753133 0.065876566
> postscript(file="/var/www/html/rcomp/tmp/1uj071258570299.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/2g9h51258570299.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/3ot8i1258570299.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/4ot6g1258570299.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/5w3dt1258570299.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
-20.4764706 -12.6858824 -8.0858824 -11.6858824 -8.6858824 -3.6858824
7 8 9 10 11 12
-3.8858824 3.7141176 14.5141176 13.6282353 12.0282353 10.6282353
13 14 15 16 17 18
8.5235294 15.3141176 20.9141176 17.3141176 21.3141176 22.3141176
19 20 21 22 23 24
23.1141176 21.7141176 29.5141176 14.6282353 13.0282353 12.6282353
25 26 27 28 29 30
10.5235294 16.3141176 16.9141176 18.3141176 27.3141176 26.3141176
31 32 33 34 35 36
22.1141176 18.7141176 21.5141176 -9.3717647 -15.9717647 -27.3717647
37 38 39 40 41 42
-21.4764706 -24.6858824 -31.0858824 -29.6858824 -34.6858824 -47.6858824
43 44 45 46 47 48
-42.8858824 -42.2858824 -56.4858824 1.0576471 -0.5423529 1.0576471
49 50 51 52 53 54
5.9529412 5.7435294 1.3435294 5.7435294 -5.2564706 2.7435294
55 56 57 58 59 60
1.5435294 -1.8564706 -9.0564706 -19.9423529 -8.5423529 3.0576471
61
16.9529412
> postscript(file="/var/www/html/rcomp/tmp/6t6qo1258570299.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 -20.4764706 NA
1 -12.6858824 -20.4764706
2 -8.0858824 -12.6858824
3 -11.6858824 -8.0858824
4 -8.6858824 -11.6858824
5 -3.6858824 -8.6858824
6 -3.8858824 -3.6858824
7 3.7141176 -3.8858824
8 14.5141176 3.7141176
9 13.6282353 14.5141176
10 12.0282353 13.6282353
11 10.6282353 12.0282353
12 8.5235294 10.6282353
13 15.3141176 8.5235294
14 20.9141176 15.3141176
15 17.3141176 20.9141176
16 21.3141176 17.3141176
17 22.3141176 21.3141176
18 23.1141176 22.3141176
19 21.7141176 23.1141176
20 29.5141176 21.7141176
21 14.6282353 29.5141176
22 13.0282353 14.6282353
23 12.6282353 13.0282353
24 10.5235294 12.6282353
25 16.3141176 10.5235294
26 16.9141176 16.3141176
27 18.3141176 16.9141176
28 27.3141176 18.3141176
29 26.3141176 27.3141176
30 22.1141176 26.3141176
31 18.7141176 22.1141176
32 21.5141176 18.7141176
33 -9.3717647 21.5141176
34 -15.9717647 -9.3717647
35 -27.3717647 -15.9717647
36 -21.4764706 -27.3717647
37 -24.6858824 -21.4764706
38 -31.0858824 -24.6858824
39 -29.6858824 -31.0858824
40 -34.6858824 -29.6858824
41 -47.6858824 -34.6858824
42 -42.8858824 -47.6858824
43 -42.2858824 -42.8858824
44 -56.4858824 -42.2858824
45 1.0576471 -56.4858824
46 -0.5423529 1.0576471
47 1.0576471 -0.5423529
48 5.9529412 1.0576471
49 5.7435294 5.9529412
50 1.3435294 5.7435294
51 5.7435294 1.3435294
52 -5.2564706 5.7435294
53 2.7435294 -5.2564706
54 1.5435294 2.7435294
55 -1.8564706 1.5435294
56 -9.0564706 -1.8564706
57 -19.9423529 -9.0564706
58 -8.5423529 -19.9423529
59 3.0576471 -8.5423529
60 16.9529412 3.0576471
61 NA 16.9529412
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -12.6858824 -20.4764706
[2,] -8.0858824 -12.6858824
[3,] -11.6858824 -8.0858824
[4,] -8.6858824 -11.6858824
[5,] -3.6858824 -8.6858824
[6,] -3.8858824 -3.6858824
[7,] 3.7141176 -3.8858824
[8,] 14.5141176 3.7141176
[9,] 13.6282353 14.5141176
[10,] 12.0282353 13.6282353
[11,] 10.6282353 12.0282353
[12,] 8.5235294 10.6282353
[13,] 15.3141176 8.5235294
[14,] 20.9141176 15.3141176
[15,] 17.3141176 20.9141176
[16,] 21.3141176 17.3141176
[17,] 22.3141176 21.3141176
[18,] 23.1141176 22.3141176
[19,] 21.7141176 23.1141176
[20,] 29.5141176 21.7141176
[21,] 14.6282353 29.5141176
[22,] 13.0282353 14.6282353
[23,] 12.6282353 13.0282353
[24,] 10.5235294 12.6282353
[25,] 16.3141176 10.5235294
[26,] 16.9141176 16.3141176
[27,] 18.3141176 16.9141176
[28,] 27.3141176 18.3141176
[29,] 26.3141176 27.3141176
[30,] 22.1141176 26.3141176
[31,] 18.7141176 22.1141176
[32,] 21.5141176 18.7141176
[33,] -9.3717647 21.5141176
[34,] -15.9717647 -9.3717647
[35,] -27.3717647 -15.9717647
[36,] -21.4764706 -27.3717647
[37,] -24.6858824 -21.4764706
[38,] -31.0858824 -24.6858824
[39,] -29.6858824 -31.0858824
[40,] -34.6858824 -29.6858824
[41,] -47.6858824 -34.6858824
[42,] -42.8858824 -47.6858824
[43,] -42.2858824 -42.8858824
[44,] -56.4858824 -42.2858824
[45,] 1.0576471 -56.4858824
[46,] -0.5423529 1.0576471
[47,] 1.0576471 -0.5423529
[48,] 5.9529412 1.0576471
[49,] 5.7435294 5.9529412
[50,] 1.3435294 5.7435294
[51,] 5.7435294 1.3435294
[52,] -5.2564706 5.7435294
[53,] 2.7435294 -5.2564706
[54,] 1.5435294 2.7435294
[55,] -1.8564706 1.5435294
[56,] -9.0564706 -1.8564706
[57,] -19.9423529 -9.0564706
[58,] -8.5423529 -19.9423529
[59,] 3.0576471 -8.5423529
[60,] 16.9529412 3.0576471
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -12.6858824 -20.4764706
2 -8.0858824 -12.6858824
3 -11.6858824 -8.0858824
4 -8.6858824 -11.6858824
5 -3.6858824 -8.6858824
6 -3.8858824 -3.6858824
7 3.7141176 -3.8858824
8 14.5141176 3.7141176
9 13.6282353 14.5141176
10 12.0282353 13.6282353
11 10.6282353 12.0282353
12 8.5235294 10.6282353
13 15.3141176 8.5235294
14 20.9141176 15.3141176
15 17.3141176 20.9141176
16 21.3141176 17.3141176
17 22.3141176 21.3141176
18 23.1141176 22.3141176
19 21.7141176 23.1141176
20 29.5141176 21.7141176
21 14.6282353 29.5141176
22 13.0282353 14.6282353
23 12.6282353 13.0282353
24 10.5235294 12.6282353
25 16.3141176 10.5235294
26 16.9141176 16.3141176
27 18.3141176 16.9141176
28 27.3141176 18.3141176
29 26.3141176 27.3141176
30 22.1141176 26.3141176
31 18.7141176 22.1141176
32 21.5141176 18.7141176
33 -9.3717647 21.5141176
34 -15.9717647 -9.3717647
35 -27.3717647 -15.9717647
36 -21.4764706 -27.3717647
37 -24.6858824 -21.4764706
38 -31.0858824 -24.6858824
39 -29.6858824 -31.0858824
40 -34.6858824 -29.6858824
41 -47.6858824 -34.6858824
42 -42.8858824 -47.6858824
43 -42.2858824 -42.8858824
44 -56.4858824 -42.2858824
45 1.0576471 -56.4858824
46 -0.5423529 1.0576471
47 1.0576471 -0.5423529
48 5.9529412 1.0576471
49 5.7435294 5.9529412
50 1.3435294 5.7435294
51 5.7435294 1.3435294
52 -5.2564706 5.7435294
53 2.7435294 -5.2564706
54 1.5435294 2.7435294
55 -1.8564706 1.5435294
56 -9.0564706 -1.8564706
57 -19.9423529 -9.0564706
58 -8.5423529 -19.9423529
59 3.0576471 -8.5423529
60 16.9529412 3.0576471
> 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/7msq01258570299.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/8frqj1258570299.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/9anhb1258570299.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/10khzt1258570299.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/11wbpf1258570299.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/12eqqf1258570299.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/13f5ck1258570300.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/14hk6k1258570300.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/15dmpr1258570300.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/16v3zz1258570300.tab")
+ }
>
> system("convert tmp/1uj071258570299.ps tmp/1uj071258570299.png")
> system("convert tmp/2g9h51258570299.ps tmp/2g9h51258570299.png")
> system("convert tmp/3ot8i1258570299.ps tmp/3ot8i1258570299.png")
> system("convert tmp/4ot6g1258570299.ps tmp/4ot6g1258570299.png")
> system("convert tmp/5w3dt1258570299.ps tmp/5w3dt1258570299.png")
> system("convert tmp/6t6qo1258570299.ps tmp/6t6qo1258570299.png")
> system("convert tmp/7msq01258570299.ps tmp/7msq01258570299.png")
> system("convert tmp/8frqj1258570299.ps tmp/8frqj1258570299.png")
> system("convert tmp/9anhb1258570299.ps tmp/9anhb1258570299.png")
> system("convert tmp/10khzt1258570299.ps tmp/10khzt1258570299.png")
>
>
> proc.time()
user system elapsed
2.394 1.539 2.928