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(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,0,510,0,514,0,517,0,508,0,493,0,490,0,469,0,478,1,528,1,534,1,518,1,506,1,502,1,516,1,528,1,533,1,536,1,537,1,524,1,536,1,587,1,597,1,581,1,564,1,558,0,575,0,580,0,575,0,563,0,552,0,537,0,545,0,601),dim=c(2,55),dimnames=list(c('X','Y'),1:55))
> y <- array(NA,dim=c(2,55),dimnames=list(c('X','Y'),1:55))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = 'Linear Trend'
> par2 = 'Include Monthly Dummies'
> par1 = '2'
> #'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 t
1 593 0 1 0 0 0 0 0 0 0 0 0 0 1
2 590 0 0 1 0 0 0 0 0 0 0 0 0 2
3 580 0 0 0 1 0 0 0 0 0 0 0 0 3
4 574 0 0 0 0 1 0 0 0 0 0 0 0 4
5 573 0 0 0 0 0 1 0 0 0 0 0 0 5
6 573 0 0 0 0 0 0 1 0 0 0 0 0 6
7 620 0 0 0 0 0 0 0 1 0 0 0 0 7
8 626 0 0 0 0 0 0 0 0 1 0 0 0 8
9 620 0 0 0 0 0 0 0 0 0 1 0 0 9
10 588 0 0 0 0 0 0 0 0 0 0 1 0 10
11 566 0 0 0 0 0 0 0 0 0 0 0 1 11
12 557 0 0 0 0 0 0 0 0 0 0 0 0 12
13 561 0 1 0 0 0 0 0 0 0 0 0 0 13
14 549 0 0 1 0 0 0 0 0 0 0 0 0 14
15 532 0 0 0 1 0 0 0 0 0 0 0 0 15
16 526 0 0 0 0 1 0 0 0 0 0 0 0 16
17 511 0 0 0 0 0 1 0 0 0 0 0 0 17
18 499 0 0 0 0 0 0 1 0 0 0 0 0 18
19 555 0 0 0 0 0 0 0 1 0 0 0 0 19
20 565 0 0 0 0 0 0 0 0 1 0 0 0 20
21 542 0 0 0 0 0 0 0 0 0 1 0 0 21
22 527 0 0 0 0 0 0 0 0 0 0 1 0 22
23 510 0 0 0 0 0 0 0 0 0 0 0 1 23
24 514 0 0 0 0 0 0 0 0 0 0 0 0 24
25 517 0 1 0 0 0 0 0 0 0 0 0 0 25
26 508 0 0 1 0 0 0 0 0 0 0 0 0 26
27 493 0 0 0 1 0 0 0 0 0 0 0 0 27
28 490 0 0 0 0 1 0 0 0 0 0 0 0 28
29 469 0 0 0 0 0 1 0 0 0 0 0 0 29
30 478 0 0 0 0 0 0 1 0 0 0 0 0 30
31 528 1 0 0 0 0 0 0 1 0 0 0 0 31
32 534 1 0 0 0 0 0 0 0 1 0 0 0 32
33 518 1 0 0 0 0 0 0 0 0 1 0 0 33
34 506 1 0 0 0 0 0 0 0 0 0 1 0 34
35 502 1 0 0 0 0 0 0 0 0 0 0 1 35
36 516 1 0 0 0 0 0 0 0 0 0 0 0 36
37 528 1 1 0 0 0 0 0 0 0 0 0 0 37
38 533 1 0 1 0 0 0 0 0 0 0 0 0 38
39 536 1 0 0 1 0 0 0 0 0 0 0 0 39
40 537 1 0 0 0 1 0 0 0 0 0 0 0 40
41 524 1 0 0 0 0 1 0 0 0 0 0 0 41
42 536 1 0 0 0 0 0 1 0 0 0 0 0 42
43 587 1 0 0 0 0 0 0 1 0 0 0 0 43
44 597 1 0 0 0 0 0 0 0 1 0 0 0 44
45 581 1 0 0 0 0 0 0 0 0 1 0 0 45
46 564 1 0 0 0 0 0 0 0 0 0 1 0 46
47 558 1 0 0 0 0 0 0 0 0 0 0 1 47
48 575 0 0 0 0 0 0 0 0 0 0 0 0 48
49 580 0 1 0 0 0 0 0 0 0 0 0 0 49
50 575 0 0 1 0 0 0 0 0 0 0 0 0 50
51 563 0 0 0 1 0 0 0 0 0 0 0 0 51
52 552 0 0 0 0 1 0 0 0 0 0 0 0 52
53 537 0 0 0 0 0 1 0 0 0 0 0 0 53
54 545 0 0 0 0 0 0 1 0 0 0 0 0 54
55 601 0 0 0 0 0 0 0 1 0 0 0 0 55
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X M1 M2 M3 M4
548.9808 -15.6216 13.7564 9.1089 -0.9386 -5.7861
M5 M6 M7 M8 M9 M10
-18.6336 -15.0811 40.1958 43.2954 28.1979 9.3504
M11 t
-2.7471 -0.1525
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-56.92 -31.30 9.85 26.68 44.19
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 548.9808 19.8543 27.651 <2e-16 ***
X -15.6216 12.2635 -1.274 0.2099
M1 13.7564 23.6366 0.582 0.5638
M2 9.1089 23.6184 0.386 0.7017
M3 -0.9386 23.6052 -0.040 0.9685
M4 -5.7861 23.5970 -0.245 0.8075
M5 -18.6336 23.5939 -0.790 0.4342
M6 -15.0811 23.5957 -0.639 0.5263
M7 40.1958 23.6489 1.700 0.0968 .
M8 43.2954 25.1710 1.720 0.0930 .
M9 28.1979 25.1343 1.122 0.2684
M10 9.3504 25.1021 0.372 0.7114
M11 -2.7471 25.0747 -0.110 0.9133
t -0.1525 0.3438 -0.444 0.6596
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 35.16 on 41 degrees of freedom
Multiple R-squared: 0.3, Adjusted R-squared: 0.07799
F-statistic: 1.351 on 13 and 41 DF, p-value: 0.2246
> 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.12773973 0.2554795 0.87226027
[2,] 0.21105199 0.4221040 0.78894801
[3,] 0.17098811 0.3419762 0.82901189
[4,] 0.12953820 0.2590764 0.87046180
[5,] 0.18847927 0.3769585 0.81152073
[6,] 0.18608923 0.3721785 0.81391077
[7,] 0.18276783 0.3655357 0.81723217
[8,] 0.20109229 0.4021846 0.79890771
[9,] 0.32085206 0.6417041 0.67914794
[10,] 0.35457196 0.7091439 0.64542804
[11,] 0.31684860 0.6336972 0.68315140
[12,] 0.32066624 0.6413325 0.67933376
[13,] 0.27757407 0.5551481 0.72242593
[14,] 0.41875187 0.8375037 0.58124813
[15,] 0.32723410 0.6544682 0.67276590
[16,] 0.23278971 0.4655794 0.76721029
[17,] 0.15669014 0.3133803 0.84330986
[18,] 0.10273661 0.2054732 0.89726339
[19,] 0.08894528 0.1778906 0.91105472
[20,] 0.31212482 0.6242496 0.68787518
[21,] 0.70242408 0.5951518 0.29757592
[22,] 0.94007457 0.1198509 0.05992543
> postscript(file="/var/www/html/rcomp/tmp/1ja0b1291031272.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/2t1iw1291031272.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/3t1iw1291031272.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/4t1iw1291031272.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/5t1iw1291031272.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 = 55
Frequency = 1
1 2 3 4 5 6
30.4153649 32.2153649 32.4153649 31.4153649 43.4153649 40.0153649
7 8 9 10 11 12
31.8910371 34.9439501 44.1939501 31.1939501 21.4439501 9.8493598
13 14 15 16 17 18
0.2455186 -6.9544814 -13.7544814 -14.7544814 -16.7544814 -32.1544814
19 20 21 22 23 24
-31.2788092 -24.2258963 -31.9758963 -27.9758963 -32.7258963 -31.3204866
25 26 27 28 29 30
-41.9243278 -46.1243278 -50.9243278 -48.9243278 -56.9243278 -51.3243278
31 32 33 34 35 36
-40.8270166 -37.7741037 -38.5241037 -31.5241037 -23.2741037 -11.8686940
37 38 39 40 41 42
-13.4725352 -3.6725352 9.5274648 15.5274648 15.5274648 24.1274648
43 44 45 46 47 48
20.0031370 27.0560499 26.3060499 28.3060499 34.5560499 33.3398207
49 50 51 52 53 54
24.7359795 24.5359795 22.7359795 16.7359795 14.7359795 19.3359795
55
20.2116517
> postscript(file="/var/www/html/rcomp/tmp/6mshz1291031272.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 = 55
Frequency = 1
lag(myerror, k = 1) myerror
0 30.4153649 NA
1 32.2153649 30.4153649
2 32.4153649 32.2153649
3 31.4153649 32.4153649
4 43.4153649 31.4153649
5 40.0153649 43.4153649
6 31.8910371 40.0153649
7 34.9439501 31.8910371
8 44.1939501 34.9439501
9 31.1939501 44.1939501
10 21.4439501 31.1939501
11 9.8493598 21.4439501
12 0.2455186 9.8493598
13 -6.9544814 0.2455186
14 -13.7544814 -6.9544814
15 -14.7544814 -13.7544814
16 -16.7544814 -14.7544814
17 -32.1544814 -16.7544814
18 -31.2788092 -32.1544814
19 -24.2258963 -31.2788092
20 -31.9758963 -24.2258963
21 -27.9758963 -31.9758963
22 -32.7258963 -27.9758963
23 -31.3204866 -32.7258963
24 -41.9243278 -31.3204866
25 -46.1243278 -41.9243278
26 -50.9243278 -46.1243278
27 -48.9243278 -50.9243278
28 -56.9243278 -48.9243278
29 -51.3243278 -56.9243278
30 -40.8270166 -51.3243278
31 -37.7741037 -40.8270166
32 -38.5241037 -37.7741037
33 -31.5241037 -38.5241037
34 -23.2741037 -31.5241037
35 -11.8686940 -23.2741037
36 -13.4725352 -11.8686940
37 -3.6725352 -13.4725352
38 9.5274648 -3.6725352
39 15.5274648 9.5274648
40 15.5274648 15.5274648
41 24.1274648 15.5274648
42 20.0031370 24.1274648
43 27.0560499 20.0031370
44 26.3060499 27.0560499
45 28.3060499 26.3060499
46 34.5560499 28.3060499
47 33.3398207 34.5560499
48 24.7359795 33.3398207
49 24.5359795 24.7359795
50 22.7359795 24.5359795
51 16.7359795 22.7359795
52 14.7359795 16.7359795
53 19.3359795 14.7359795
54 20.2116517 19.3359795
55 NA 20.2116517
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 32.2153649 30.4153649
[2,] 32.4153649 32.2153649
[3,] 31.4153649 32.4153649
[4,] 43.4153649 31.4153649
[5,] 40.0153649 43.4153649
[6,] 31.8910371 40.0153649
[7,] 34.9439501 31.8910371
[8,] 44.1939501 34.9439501
[9,] 31.1939501 44.1939501
[10,] 21.4439501 31.1939501
[11,] 9.8493598 21.4439501
[12,] 0.2455186 9.8493598
[13,] -6.9544814 0.2455186
[14,] -13.7544814 -6.9544814
[15,] -14.7544814 -13.7544814
[16,] -16.7544814 -14.7544814
[17,] -32.1544814 -16.7544814
[18,] -31.2788092 -32.1544814
[19,] -24.2258963 -31.2788092
[20,] -31.9758963 -24.2258963
[21,] -27.9758963 -31.9758963
[22,] -32.7258963 -27.9758963
[23,] -31.3204866 -32.7258963
[24,] -41.9243278 -31.3204866
[25,] -46.1243278 -41.9243278
[26,] -50.9243278 -46.1243278
[27,] -48.9243278 -50.9243278
[28,] -56.9243278 -48.9243278
[29,] -51.3243278 -56.9243278
[30,] -40.8270166 -51.3243278
[31,] -37.7741037 -40.8270166
[32,] -38.5241037 -37.7741037
[33,] -31.5241037 -38.5241037
[34,] -23.2741037 -31.5241037
[35,] -11.8686940 -23.2741037
[36,] -13.4725352 -11.8686940
[37,] -3.6725352 -13.4725352
[38,] 9.5274648 -3.6725352
[39,] 15.5274648 9.5274648
[40,] 15.5274648 15.5274648
[41,] 24.1274648 15.5274648
[42,] 20.0031370 24.1274648
[43,] 27.0560499 20.0031370
[44,] 26.3060499 27.0560499
[45,] 28.3060499 26.3060499
[46,] 34.5560499 28.3060499
[47,] 33.3398207 34.5560499
[48,] 24.7359795 33.3398207
[49,] 24.5359795 24.7359795
[50,] 22.7359795 24.5359795
[51,] 16.7359795 22.7359795
[52,] 14.7359795 16.7359795
[53,] 19.3359795 14.7359795
[54,] 20.2116517 19.3359795
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 32.2153649 30.4153649
2 32.4153649 32.2153649
3 31.4153649 32.4153649
4 43.4153649 31.4153649
5 40.0153649 43.4153649
6 31.8910371 40.0153649
7 34.9439501 31.8910371
8 44.1939501 34.9439501
9 31.1939501 44.1939501
10 21.4439501 31.1939501
11 9.8493598 21.4439501
12 0.2455186 9.8493598
13 -6.9544814 0.2455186
14 -13.7544814 -6.9544814
15 -14.7544814 -13.7544814
16 -16.7544814 -14.7544814
17 -32.1544814 -16.7544814
18 -31.2788092 -32.1544814
19 -24.2258963 -31.2788092
20 -31.9758963 -24.2258963
21 -27.9758963 -31.9758963
22 -32.7258963 -27.9758963
23 -31.3204866 -32.7258963
24 -41.9243278 -31.3204866
25 -46.1243278 -41.9243278
26 -50.9243278 -46.1243278
27 -48.9243278 -50.9243278
28 -56.9243278 -48.9243278
29 -51.3243278 -56.9243278
30 -40.8270166 -51.3243278
31 -37.7741037 -40.8270166
32 -38.5241037 -37.7741037
33 -31.5241037 -38.5241037
34 -23.2741037 -31.5241037
35 -11.8686940 -23.2741037
36 -13.4725352 -11.8686940
37 -3.6725352 -13.4725352
38 9.5274648 -3.6725352
39 15.5274648 9.5274648
40 15.5274648 15.5274648
41 24.1274648 15.5274648
42 20.0031370 24.1274648
43 27.0560499 20.0031370
44 26.3060499 27.0560499
45 28.3060499 26.3060499
46 34.5560499 28.3060499
47 33.3398207 34.5560499
48 24.7359795 33.3398207
49 24.5359795 24.7359795
50 22.7359795 24.5359795
51 16.7359795 22.7359795
52 14.7359795 16.7359795
53 19.3359795 14.7359795
54 20.2116517 19.3359795
> 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/7xjy21291031272.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/8xjy21291031272.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/9xjy21291031272.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/10qbg51291031272.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/11bbea1291031272.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/12euug1291031272.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/13ld9a1291031272.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/14j7gj1291031272.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/15hm7j1291031272.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/16vena1291031272.tab")
+ }
> try(system("convert tmp/1ja0b1291031272.ps tmp/1ja0b1291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/2t1iw1291031272.ps tmp/2t1iw1291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/3t1iw1291031272.ps tmp/3t1iw1291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/4t1iw1291031272.ps tmp/4t1iw1291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/5t1iw1291031272.ps tmp/5t1iw1291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/6mshz1291031272.ps tmp/6mshz1291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/7xjy21291031272.ps tmp/7xjy21291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/8xjy21291031272.ps tmp/8xjy21291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/9xjy21291031272.ps tmp/9xjy21291031272.png",intern=TRUE))
character(0)
> try(system("convert tmp/10qbg51291031272.ps tmp/10qbg51291031272.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.296 1.573 5.758