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.
Natural language support but running in an English locale
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(13768040.14
+ ,14731798.37
+ ,17487530.67
+ ,16471559.62
+ ,16198106.13
+ ,15213975.95
+ ,17535166.38
+ ,17637387.4
+ ,16571771.60
+ ,17972385.83
+ ,16198892.67
+ ,16896235.55
+ ,16554237.93
+ ,16697955.94
+ ,19554176.37
+ ,19691579.52
+ ,15903762.33
+ ,15930700.75
+ ,18003781.65
+ ,17444615.98
+ ,18329610.38
+ ,17699369.88
+ ,16260733.42
+ ,15189796.81
+ ,14851949.20
+ ,15672722.75
+ ,18174068.44
+ ,17180794.3
+ ,18406552.23
+ ,17664893.45
+ ,18466459.42
+ ,17862884.98
+ ,16016524.60
+ ,16162288.88
+ ,17428458.32
+ ,17463628.82
+ ,17167191.42
+ ,16772112.17
+ ,19629987.60
+ ,19106861.48
+ ,17183629.01
+ ,16721314.25
+ ,18344657.85
+ ,18161267.85
+ ,19301440.71
+ ,18509941.2
+ ,18147463.68
+ ,17802737.97
+ ,16192909.22
+ ,16409869.75
+ ,18374420.60
+ ,17967742.04
+ ,20515191.95
+ ,20286602.27
+ ,18957217.20
+ ,19537280.81
+ ,16471529.53
+ ,18021889.62
+ ,18746813.27
+ ,20194317.23
+ ,19009453.59
+ ,19049596.62
+ ,19211178.55
+ ,20244720.94
+ ,20547653.75
+ ,21473302.24
+ ,19325754.03
+ ,19673603.19
+ ,20605542.58
+ ,21053177.29
+ ,20056915.06
+ ,20159479.84
+ ,16141449.72
+ ,18203628.31
+ ,20359793.22
+ ,21289464.94
+ ,19711553.27
+ ,20432335.71
+ ,15638580.70
+ ,17180395.07
+ ,14384486.00
+ ,15816786.32
+ ,13855616.12
+ ,15071819.75
+ ,14308336.46
+ ,14521120.61
+ ,15290621.44
+ ,15668789.39
+ ,14423755.53
+ ,14346884.11
+ ,13779681.49
+ ,13881008.13
+ ,15686348.94
+ ,15465943.69
+ ,14733828.17
+ ,14238232.92
+ ,12522497.94
+ ,13557713.21
+ ,16189383.57
+ ,16127590.29
+ ,16059123.25
+ ,16793894.2
+ ,16007123.26
+ ,16014007.43
+ ,15806842.33
+ ,16867867.15
+ ,15159951.13
+ ,16014583.21
+ ,15692144.17
+ ,15878594.85
+ ,18908869.11
+ ,18664899.14
+ ,16969881.42
+ ,17962530.06
+ ,16997477.78
+ ,17332692.2
+ ,19858875.65
+ ,19542066.35
+ ,17681170.13
+ ,17203555.19)
+ ,dim=c(2
+ ,60)
+ ,dimnames=list(c('Uitvoer'
+ ,'Invoer')
+ ,1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('Uitvoer','Invoer'),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 = 'Linear Trend'
> par2 = 'Do not include Seasonal 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
Invoer Uitvoer t
1 14731798 13768040 1
2 16471560 17487531 2
3 15213976 16198106 3
4 17637387 17535166 4
5 17972386 16571772 5
6 16896236 16198893 6
7 16697956 16554238 7
8 19691580 19554176 8
9 15930701 15903762 9
10 17444616 18003782 10
11 17699370 18329610 11
12 15189797 16260733 12
13 15672723 14851949 13
14 17180794 18174068 14
15 17664893 18406552 15
16 17862885 18466459 16
17 16162289 16016525 17
18 17463629 17428458 18
19 16772112 17167191 19
20 19106861 19629988 20
21 16721314 17183629 21
22 18161268 18344658 22
23 18509941 19301441 23
24 17802738 18147464 24
25 16409870 16192909 25
26 17967742 18374421 26
27 20286602 20515192 27
28 19537281 18957217 28
29 18021890 16471530 29
30 20194317 18746813 30
31 19049597 19009454 31
32 20244721 19211179 32
33 21473302 20547654 33
34 19673603 19325754 34
35 21053177 20605543 35
36 20159480 20056915 36
37 18203628 16141450 37
38 21289465 20359793 38
39 20432336 19711553 39
40 17180395 15638581 40
41 15816786 14384486 41
42 15071820 13855616 42
43 14521121 14308336 43
44 15668789 15290621 44
45 14346884 14423756 45
46 13881008 13779681 46
47 15465944 15686349 47
48 14238233 14733828 48
49 13557713 12522498 49
50 16127590 16189384 50
51 16793894 16059123 51
52 16014007 16007123 52
53 16867867 15806842 53
54 16014583 15159951 54
55 15878595 15692144 55
56 18664899 18908869 56
57 17962530 16969881 57
58 17332692 16997478 58
59 19542066 19858876 59
60 17203555 17681170 60
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Uitvoer t
1.286e+06 9.195e-01 1.032e+04
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1171512 -555526 -53322 508458 1694004
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.286e+06 8.547e+05 1.504 0.1380
Uitvoer 9.195e-01 4.714e-02 19.508 <2e-16 ***
t 1.032e+04 5.369e+03 1.922 0.0596 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 711600 on 57 degrees of freedom
Multiple R-squared: 0.8701, Adjusted R-squared: 0.8655
F-statistic: 190.9 on 2 and 57 DF, p-value: < 2.2e-16
> 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.7956903 0.4086193 0.20430967
[2,] 0.7973007 0.4053986 0.20269928
[3,] 0.7129423 0.5741154 0.28705770
[4,] 0.7748532 0.4502935 0.22514676
[5,] 0.7519717 0.4960565 0.24802826
[6,] 0.6994654 0.6010693 0.30053463
[7,] 0.7829915 0.4340170 0.21700852
[8,] 0.7643670 0.4712660 0.23563300
[9,] 0.7446300 0.5107401 0.25537003
[10,] 0.6898162 0.6203676 0.31018378
[11,] 0.6298239 0.7403521 0.37017607
[12,] 0.5517224 0.8965553 0.44827763
[13,] 0.4886002 0.9772005 0.51139976
[14,] 0.4309188 0.8618376 0.56908119
[15,] 0.3935989 0.7871978 0.60640111
[16,] 0.3667285 0.7334571 0.63327145
[17,] 0.3319263 0.6638526 0.66807372
[18,] 0.3574721 0.7149442 0.64252791
[19,] 0.3609924 0.7219848 0.63900762
[20,] 0.3347101 0.6694201 0.66528994
[21,] 0.3952697 0.7905393 0.60473033
[22,] 0.4823342 0.9646684 0.51766582
[23,] 0.5482514 0.9034972 0.45174862
[24,] 0.6432839 0.7134322 0.35671612
[25,] 0.7658949 0.4682102 0.23410511
[26,] 0.7651295 0.4697410 0.23487052
[27,] 0.7606164 0.4787673 0.23938364
[28,] 0.7593092 0.4813815 0.24069076
[29,] 0.7113801 0.5772397 0.28861986
[30,] 0.6555917 0.6888166 0.34440829
[31,] 0.6519051 0.6961898 0.34809488
[32,] 0.7474035 0.5051930 0.25259650
[33,] 0.7048253 0.5903494 0.29517470
[34,] 0.6404207 0.7191586 0.35957932
[35,] 0.7324875 0.5350251 0.26751254
[36,] 0.8404594 0.3190812 0.15954062
[37,] 0.9120251 0.1759497 0.08797486
[38,] 0.9125654 0.1748692 0.08743461
[39,] 0.9210324 0.1579353 0.07896764
[40,] 0.9093645 0.1812710 0.09063550
[41,] 0.8868269 0.2263461 0.11317305
[42,] 0.8544127 0.2911745 0.14558725
[43,] 0.9383864 0.1232272 0.06161362
[44,] 0.8978288 0.2043423 0.10217117
[45,] 0.8909693 0.2180613 0.10903067
[46,] 0.8175839 0.3648323 0.18241615
[47,] 0.8533470 0.2933060 0.14665299
[48,] 0.7605950 0.4788101 0.23940505
[49,] 0.6009712 0.7980575 0.39902877
> postscript(file="/var/www/html/freestat/rcomp/tmp/1a4au1291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index')
> points(x[,1]-mysum$resid)
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/2leax1291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/3leax1291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/4leax1291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals')
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/5wn901291146188.ps",horizontal=F,onefile=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
776032.38 -914587.41 -996868.66 186799.76 1397317.16 653708.16
7 8 9 10 11 12
118369.72 343238.41 -71414.67 -498781.21 -553945.77 -1171512.09
13 14 15 16 17 18
596467.16 -960460.53 -700449.43 -567862.25 -26069.96 -33318.99
19 20 21 22 23 24
-494921.31 -435026.08 -581473.25 -219402.29 -760807.97 -417252.39
25 26 27 28 29 30
-23233.08 -481574.28 -141467.20 531444.94 1291316.78 1361307.48
31 32 33 34 35 36
-35230.03 964088.90 953465.12 266979.61 459471.85 59916.02
37 38 39 40 41 42
1694004.12 890765.79 619371.54 1102197.98 881405.96 612413.92
43 44 45 46 47 48
-364880.16 -130739.54 -665883.88 -549855.43 -728415.13 -1090605.56
49 50 51 52 53 54
251866.87 -560267.03 215491.03 -526901.75 500795.88 232006.74
55 56 57 58 59 60
-403651.49 -585436.68 484768.17 -180764.32 -612757.39 -959194.25
> postscript(file="/var/www/html/freestat/rcomp/tmp/6wn901291146188.ps",horizontal=F,onefile=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 776032.38 NA
1 -914587.41 776032.38
2 -996868.66 -914587.41
3 186799.76 -996868.66
4 1397317.16 186799.76
5 653708.16 1397317.16
6 118369.72 653708.16
7 343238.41 118369.72
8 -71414.67 343238.41
9 -498781.21 -71414.67
10 -553945.77 -498781.21
11 -1171512.09 -553945.77
12 596467.16 -1171512.09
13 -960460.53 596467.16
14 -700449.43 -960460.53
15 -567862.25 -700449.43
16 -26069.96 -567862.25
17 -33318.99 -26069.96
18 -494921.31 -33318.99
19 -435026.08 -494921.31
20 -581473.25 -435026.08
21 -219402.29 -581473.25
22 -760807.97 -219402.29
23 -417252.39 -760807.97
24 -23233.08 -417252.39
25 -481574.28 -23233.08
26 -141467.20 -481574.28
27 531444.94 -141467.20
28 1291316.78 531444.94
29 1361307.48 1291316.78
30 -35230.03 1361307.48
31 964088.90 -35230.03
32 953465.12 964088.90
33 266979.61 953465.12
34 459471.85 266979.61
35 59916.02 459471.85
36 1694004.12 59916.02
37 890765.79 1694004.12
38 619371.54 890765.79
39 1102197.98 619371.54
40 881405.96 1102197.98
41 612413.92 881405.96
42 -364880.16 612413.92
43 -130739.54 -364880.16
44 -665883.88 -130739.54
45 -549855.43 -665883.88
46 -728415.13 -549855.43
47 -1090605.56 -728415.13
48 251866.87 -1090605.56
49 -560267.03 251866.87
50 215491.03 -560267.03
51 -526901.75 215491.03
52 500795.88 -526901.75
53 232006.74 500795.88
54 -403651.49 232006.74
55 -585436.68 -403651.49
56 484768.17 -585436.68
57 -180764.32 484768.17
58 -612757.39 -180764.32
59 -959194.25 -612757.39
60 NA -959194.25
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -914587.41 776032.38
[2,] -996868.66 -914587.41
[3,] 186799.76 -996868.66
[4,] 1397317.16 186799.76
[5,] 653708.16 1397317.16
[6,] 118369.72 653708.16
[7,] 343238.41 118369.72
[8,] -71414.67 343238.41
[9,] -498781.21 -71414.67
[10,] -553945.77 -498781.21
[11,] -1171512.09 -553945.77
[12,] 596467.16 -1171512.09
[13,] -960460.53 596467.16
[14,] -700449.43 -960460.53
[15,] -567862.25 -700449.43
[16,] -26069.96 -567862.25
[17,] -33318.99 -26069.96
[18,] -494921.31 -33318.99
[19,] -435026.08 -494921.31
[20,] -581473.25 -435026.08
[21,] -219402.29 -581473.25
[22,] -760807.97 -219402.29
[23,] -417252.39 -760807.97
[24,] -23233.08 -417252.39
[25,] -481574.28 -23233.08
[26,] -141467.20 -481574.28
[27,] 531444.94 -141467.20
[28,] 1291316.78 531444.94
[29,] 1361307.48 1291316.78
[30,] -35230.03 1361307.48
[31,] 964088.90 -35230.03
[32,] 953465.12 964088.90
[33,] 266979.61 953465.12
[34,] 459471.85 266979.61
[35,] 59916.02 459471.85
[36,] 1694004.12 59916.02
[37,] 890765.79 1694004.12
[38,] 619371.54 890765.79
[39,] 1102197.98 619371.54
[40,] 881405.96 1102197.98
[41,] 612413.92 881405.96
[42,] -364880.16 612413.92
[43,] -130739.54 -364880.16
[44,] -665883.88 -130739.54
[45,] -549855.43 -665883.88
[46,] -728415.13 -549855.43
[47,] -1090605.56 -728415.13
[48,] 251866.87 -1090605.56
[49,] -560267.03 251866.87
[50,] 215491.03 -560267.03
[51,] -526901.75 215491.03
[52,] 500795.88 -526901.75
[53,] 232006.74 500795.88
[54,] -403651.49 232006.74
[55,] -585436.68 -403651.49
[56,] 484768.17 -585436.68
[57,] -180764.32 484768.17
[58,] -612757.39 -180764.32
[59,] -959194.25 -612757.39
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -914587.41 776032.38
2 -996868.66 -914587.41
3 186799.76 -996868.66
4 1397317.16 186799.76
5 653708.16 1397317.16
6 118369.72 653708.16
7 343238.41 118369.72
8 -71414.67 343238.41
9 -498781.21 -71414.67
10 -553945.77 -498781.21
11 -1171512.09 -553945.77
12 596467.16 -1171512.09
13 -960460.53 596467.16
14 -700449.43 -960460.53
15 -567862.25 -700449.43
16 -26069.96 -567862.25
17 -33318.99 -26069.96
18 -494921.31 -33318.99
19 -435026.08 -494921.31
20 -581473.25 -435026.08
21 -219402.29 -581473.25
22 -760807.97 -219402.29
23 -417252.39 -760807.97
24 -23233.08 -417252.39
25 -481574.28 -23233.08
26 -141467.20 -481574.28
27 531444.94 -141467.20
28 1291316.78 531444.94
29 1361307.48 1291316.78
30 -35230.03 1361307.48
31 964088.90 -35230.03
32 953465.12 964088.90
33 266979.61 953465.12
34 459471.85 266979.61
35 59916.02 459471.85
36 1694004.12 59916.02
37 890765.79 1694004.12
38 619371.54 890765.79
39 1102197.98 619371.54
40 881405.96 1102197.98
41 612413.92 881405.96
42 -364880.16 612413.92
43 -130739.54 -364880.16
44 -665883.88 -130739.54
45 -549855.43 -665883.88
46 -728415.13 -549855.43
47 -1090605.56 -728415.13
48 251866.87 -1090605.56
49 -560267.03 251866.87
50 215491.03 -560267.03
51 -526901.75 215491.03
52 500795.88 -526901.75
53 232006.74 500795.88
54 -403651.49 232006.74
55 -585436.68 -403651.49
56 484768.17 -585436.68
57 -180764.32 484768.17
58 -612757.39 -180764.32
59 -959194.25 -612757.39
> plot(z,main=paste('Residual Lag plot, lowess, and regression line'), ylab='values of Residuals', xlab='lagged values of Residuals')
> lines(lowess(z))
> abline(lm(z))
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/7pw831291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/8pw831291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/9h5p61291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
> plot(mylm, las = 1, sub='Residual Diagnostics')
> par(opar)
> dev.off()
null device
1
> if (n > n25) {
+ postscript(file="/var/www/html/freestat/rcomp/tmp/10h5p61291146188.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
+ plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint')
+ grid()
+ dev.off()
+ }
null device
1
>
> #Note: the /var/www/html/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/html/freestat/rcomp/createtable")
>
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Estimated Regression Equation', 1, TRUE)
> a<-table.row.end(a)
> myeq <- colnames(x)[1]
> myeq <- paste(myeq, '[t] = ', sep='')
> for (i in 1:k){
+ if (mysum$coefficients[i,1] > 0) myeq <- paste(myeq, '+', '')
+ myeq <- paste(myeq, mysum$coefficients[i,1], sep=' ')
+ if (rownames(mysum$coefficients)[i] != '(Intercept)') {
+ myeq <- paste(myeq, rownames(mysum$coefficients)[i], sep='')
+ if (rownames(mysum$coefficients)[i] != 't') myeq <- paste(myeq, '[t]', sep='')
+ }
+ }
> myeq <- paste(myeq, ' + e[t]')
> a<-table.row.start(a)
> a<-table.element(a, myeq)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/1136ou1291146188.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,hyperlink('http://www.xycoon.com/ols1.htm','Multiple Linear Regression - Ordinary Least Squares',''), 6, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Variable',header=TRUE)
> a<-table.element(a,'Parameter',header=TRUE)
> a<-table.element(a,'S.D.',header=TRUE)
> a<-table.element(a,'T-STAT
H0: parameter = 0',header=TRUE)
> a<-table.element(a,'2-tail p-value',header=TRUE)
> a<-table.element(a,'1-tail p-value',header=TRUE)
> a<-table.row.end(a)
> for (i in 1:k){
+ a<-table.row.start(a)
+ a<-table.element(a,rownames(mysum$coefficients)[i],header=TRUE)
+ a<-table.element(a,mysum$coefficients[i,1])
+ a<-table.element(a, round(mysum$coefficients[i,2],6))
+ a<-table.element(a, round(mysum$coefficients[i,3],4))
+ a<-table.element(a, round(mysum$coefficients[i,4],6))
+ a<-table.element(a, round(mysum$coefficients[i,4]/2,6))
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/126o4h1291146188.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Regression Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple R',1,TRUE)
> a<-table.element(a, sqrt(mysum$r.squared))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'R-squared',1,TRUE)
> a<-table.element(a, mysum$r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Adjusted R-squared',1,TRUE)
> a<-table.element(a, mysum$adj.r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (value)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[3])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'p-value',1,TRUE)
> a<-table.element(a, 1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3]))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Residual Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Residual Standard Deviation',1,TRUE)
> a<-table.element(a, mysum$sigma)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Sum Squared Residuals',1,TRUE)
> a<-table.element(a, sum(myerror*myerror))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/135i4o1291146189.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Actuals, Interpolation, and Residuals', 4, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Time or Index', 1, TRUE)
> a<-table.element(a, 'Actuals', 1, TRUE)
> a<-table.element(a, 'Interpolation
Forecast', 1, TRUE)
> a<-table.element(a, 'Residuals
Prediction Error', 1, TRUE)
> a<-table.row.end(a)
> for (i in 1:n) {
+ a<-table.row.start(a)
+ a<-table.element(a,i, 1, TRUE)
+ a<-table.element(a,x[i])
+ a<-table.element(a,x[i]-mysum$resid[i])
+ a<-table.element(a,mysum$resid[i])
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/14grlr1291146189.tab")
> if (n > n25) {
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'p-values',header=TRUE)
+ a<-table.element(a,'Alternative Hypothesis',3,header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'breakpoint index',header=TRUE)
+ a<-table.element(a,'greater',header=TRUE)
+ a<-table.element(a,'2-sided',header=TRUE)
+ a<-table.element(a,'less',header=TRUE)
+ a<-table.row.end(a)
+ for (mypoint in kp3:nmkm3) {
+ a<-table.row.start(a)
+ a<-table.element(a,mypoint,header=TRUE)
+ a<-table.element(a,gqarr[mypoint-kp3+1,1])
+ a<-table.element(a,gqarr[mypoint-kp3+1,2])
+ a<-table.element(a,gqarr[mypoint-kp3+1,3])
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/www/html/freestat/rcomp/tmp/151rkf1291146189.tab")
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Meta Analysis of Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'Description',header=TRUE)
+ a<-table.element(a,'# significant tests',header=TRUE)
+ a<-table.element(a,'% significant tests',header=TRUE)
+ a<-table.element(a,'OK/NOK',header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'1% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant1)
+ a<-table.element(a,numsignificant1/numgqtests)
+ if (numsignificant1/numgqtests < 0.01) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'5% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant5)
+ a<-table.element(a,numsignificant5/numgqtests)
+ if (numsignificant5/numgqtests < 0.05) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'10% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant10)
+ a<-table.element(a,numsignificant10/numgqtests)
+ if (numsignificant10/numgqtests < 0.1) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.end(a)
+ table.save(a,file="/var/www/html/freestat/rcomp/tmp/16x1i61291146189.tab")
+ }
>
> try(system("convert tmp/1a4au1291146188.ps tmp/1a4au1291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/2leax1291146188.ps tmp/2leax1291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/3leax1291146188.ps tmp/3leax1291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/4leax1291146188.ps tmp/4leax1291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/5wn901291146188.ps tmp/5wn901291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/6wn901291146188.ps tmp/6wn901291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/7pw831291146188.ps tmp/7pw831291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/8pw831291146188.ps tmp/8pw831291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/9h5p61291146188.ps tmp/9h5p61291146188.png",intern=TRUE))
character(0)
> try(system("convert tmp/10h5p61291146188.ps tmp/10h5p61291146188.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.982 2.544 4.543