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(16198.9,16896.2,16554.2,16698,19554.2,19691.6,15903.8,15930.7,18003.8,17444.6,18329.6,17699.4,16260.7,15189.8,14851.9,15672.7,18174.1,17180.8,18406.6,17664.9,18466.5,17862.9,16016.5,16162.3,17428.5,17463.6,17167.2,16772.1,19630,19106.9,17183.6,16721.3,18344.7,18161.3,19301.4,18509.9,18147.5,17802.7,16192.9,16409.9,18374.4,17967.7,20515.2,20286.6,18957.2,19537.3,16471.5,18021.9,18746.8,20194.3,19009.5,19049.6,19211.2,20244.7,20547.7,21473.3,19325.8,19673.6,20605.5,21053.2,20056.9,20159.5,16141.4,18203.6,20359.8,21289.5,19711.6,20432.3,15638.6,17180.4,14384.5,15816.8,13855.6,15071.8,14308.3,14521.1,15290.6,15668.8,14423.8,14346.9,13779.7,13881,15686.3,15465.9,14733.8,14238.2,12522.5,13557.7,16189.4,16127.6,16059.1,16793.9,16007.1,16014,15806.8,16867.9,15160,16014.6,15692.1,15878.6,18908.9,18664.9,16969.9,17962.5,16997.5,17332.7,19858.9,19542.1,17681.2,17203.6),dim=c(2,55),dimnames=list(c('uitvoer','invoer'),1:55))
> y <- array(NA,dim=c(2,55),dimnames=list(c('uitvoer','invoer'),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 = '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
uitvoer invoer M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 16198.9 16896.2 1 0 0 0 0 0 0 0 0 0 0
2 16554.2 16698.0 0 1 0 0 0 0 0 0 0 0 0
3 19554.2 19691.6 0 0 1 0 0 0 0 0 0 0 0
4 15903.8 15930.7 0 0 0 1 0 0 0 0 0 0 0
5 18003.8 17444.6 0 0 0 0 1 0 0 0 0 0 0
6 18329.6 17699.4 0 0 0 0 0 1 0 0 0 0 0
7 16260.7 15189.8 0 0 0 0 0 0 1 0 0 0 0
8 14851.9 15672.7 0 0 0 0 0 0 0 1 0 0 0
9 18174.1 17180.8 0 0 0 0 0 0 0 0 1 0 0
10 18406.6 17664.9 0 0 0 0 0 0 0 0 0 1 0
11 18466.5 17862.9 0 0 0 0 0 0 0 0 0 0 1
12 16016.5 16162.3 0 0 0 0 0 0 0 0 0 0 0
13 17428.5 17463.6 1 0 0 0 0 0 0 0 0 0 0
14 17167.2 16772.1 0 1 0 0 0 0 0 0 0 0 0
15 19630.0 19106.9 0 0 1 0 0 0 0 0 0 0 0
16 17183.6 16721.3 0 0 0 1 0 0 0 0 0 0 0
17 18344.7 18161.3 0 0 0 0 1 0 0 0 0 0 0
18 19301.4 18509.9 0 0 0 0 0 1 0 0 0 0 0
19 18147.5 17802.7 0 0 0 0 0 0 1 0 0 0 0
20 16192.9 16409.9 0 0 0 0 0 0 0 1 0 0 0
21 18374.4 17967.7 0 0 0 0 0 0 0 0 1 0 0
22 20515.2 20286.6 0 0 0 0 0 0 0 0 0 1 0
23 18957.2 19537.3 0 0 0 0 0 0 0 0 0 0 1
24 16471.5 18021.9 0 0 0 0 0 0 0 0 0 0 0
25 18746.8 20194.3 1 0 0 0 0 0 0 0 0 0 0
26 19009.5 19049.6 0 1 0 0 0 0 0 0 0 0 0
27 19211.2 20244.7 0 0 1 0 0 0 0 0 0 0 0
28 20547.7 21473.3 0 0 0 1 0 0 0 0 0 0 0
29 19325.8 19673.6 0 0 0 0 1 0 0 0 0 0 0
30 20605.5 21053.2 0 0 0 0 0 1 0 0 0 0 0
31 20056.9 20159.5 0 0 0 0 0 0 1 0 0 0 0
32 16141.4 18203.6 0 0 0 0 0 0 0 1 0 0 0
33 20359.8 21289.5 0 0 0 0 0 0 0 0 1 0 0
34 19711.6 20432.3 0 0 0 0 0 0 0 0 0 1 0
35 15638.6 17180.4 0 0 0 0 0 0 0 0 0 0 1
36 14384.5 15816.8 0 0 0 0 0 0 0 0 0 0 0
37 13855.6 15071.8 1 0 0 0 0 0 0 0 0 0 0
38 14308.3 14521.1 0 1 0 0 0 0 0 0 0 0 0
39 15290.6 15668.8 0 0 1 0 0 0 0 0 0 0 0
40 14423.8 14346.9 0 0 0 1 0 0 0 0 0 0 0
41 13779.7 13881.0 0 0 0 0 1 0 0 0 0 0 0
42 15686.3 15465.9 0 0 0 0 0 1 0 0 0 0 0
43 14733.8 14238.2 0 0 0 0 0 0 1 0 0 0 0
44 12522.5 13557.7 0 0 0 0 0 0 0 1 0 0 0
45 16189.4 16127.6 0 0 0 0 0 0 0 0 1 0 0
46 16059.1 16793.9 0 0 0 0 0 0 0 0 0 1 0
47 16007.1 16014.0 0 0 0 0 0 0 0 0 0 0 1
48 15806.8 16867.9 0 0 0 0 0 0 0 0 0 0 0
49 15160.0 16014.6 1 0 0 0 0 0 0 0 0 0 0
50 15692.1 15878.6 0 1 0 0 0 0 0 0 0 0 0
51 18908.9 18664.9 0 0 1 0 0 0 0 0 0 0 0
52 16969.9 17962.5 0 0 0 1 0 0 0 0 0 0 0
53 16997.5 17332.7 0 0 0 0 1 0 0 0 0 0 0
54 19858.9 19542.1 0 0 0 0 0 1 0 0 0 0 0
55 17681.2 17203.6 0 0 0 0 0 0 1 0 0 0 0
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) invoer M1 M2 M3 M4
748.0722 0.8926 241.3890 995.4584 1101.3108 827.4088
M5 M6 M7 M8 M9 M10
1101.5054 1536.1848 1526.3054 -67.6232 1333.3850 1149.1965
M11
766.1374
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1210.791 -379.801 2.483 341.833 1007.911
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 748.07215 789.35027 0.948 0.348700
invoer 0.89260 0.04399 20.291 < 2e-16 ***
M1 241.38901 385.22372 0.627 0.534299
M2 995.45841 384.84417 2.587 0.013247 *
M3 1101.31079 394.32318 2.793 0.007833 **
M4 827.40883 385.61475 2.146 0.037724 *
M5 1101.50544 385.64854 2.856 0.006637 **
M6 1536.18476 392.31175 3.916 0.000325 ***
M7 1526.30537 384.90158 3.965 0.000280 ***
M8 -67.62316 406.97625 -0.166 0.868828
M9 1333.38499 410.42413 3.249 0.002284 **
M10 1149.19652 415.77962 2.764 0.008443 **
M11 766.13740 407.67852 1.879 0.067158 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 573.6 on 42 degrees of freedom
Multiple R-squared: 0.9372, Adjusted R-squared: 0.9193
F-statistic: 52.25 on 12 and 42 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.5268496 0.9463008 0.4731504
[2,] 0.4500496 0.9000992 0.5499504
[3,] 0.3312893 0.6625787 0.6687107
[4,] 0.2780630 0.5561260 0.7219370
[5,] 0.4172149 0.8344298 0.5827851
[6,] 0.4087961 0.8175921 0.5912039
[7,] 0.4394497 0.8788995 0.5605503
[8,] 0.5244526 0.9510948 0.4755474
[9,] 0.5662495 0.8675010 0.4337505
[10,] 0.4590491 0.9180981 0.5409509
[11,] 0.4783914 0.9567828 0.5216086
[12,] 0.6253327 0.7493346 0.3746673
[13,] 0.5435883 0.9128234 0.4564117
[14,] 0.4644300 0.9288600 0.5355700
[15,] 0.3921140 0.7842279 0.6078860
[16,] 0.3019192 0.6038384 0.6980808
[17,] 0.3590038 0.7180077 0.6409962
[18,] 0.3455416 0.6910832 0.6544584
[19,] 0.3018622 0.6037244 0.6981378
[20,] 0.8821297 0.2357406 0.1178703
[21,] 0.8556455 0.2887091 0.1443545
[22,] 0.8388314 0.3223372 0.1611686
[23,] 0.7449621 0.5100758 0.2550379
[24,] 0.7889035 0.4221930 0.2110965
> postscript(file="/var/www/html/freestat/rcomp/tmp/12nuh1290170821.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/freestat/rcomp/tmp/22nuh1290170821.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/freestat/rcomp/tmp/3vfu21290170821.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/freestat/rcomp/tmp/4vfu21290170821.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/freestat/rcomp/tmp/5vfu21290170821.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
127.933355 -93.923224 128.144593 108.616401 583.216472 246.903315
7 8 9 10 11 12
427.945333 182.038538 757.104136 741.686165 1007.910991 841.999655
13 14 15 16 17 18
851.073547 452.935303 725.846337 682.728837 284.391861 495.253061
19 20 21 22 23 24
-17.522613 865.015678 255.019182 510.163361 4.045776 -362.874613
25 26 27 28 29 30
-268.042382 262.344550 -708.551072 -194.794371 -84.383303 -470.790101
31 32 33 34 35 36
-211.796346 -787.536416 -724.611115 -423.488092 -1210.791232 -481.607917
37 38 39 40 41 42
-586.911809 -396.727778 -544.714279 42.312285 -460.023161 -402.780221
43 44 45 46 47 48
-249.558909 -259.517800 -287.512202 -828.361434 198.834465 2.482875
49 50 51 52 53 54
-124.052710 -224.628852 399.274422 -638.863151 -323.201870 131.413946
55
50.932535
> postscript(file="/var/www/html/freestat/rcomp/tmp/6n6bn1290170821.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 127.933355 NA
1 -93.923224 127.933355
2 128.144593 -93.923224
3 108.616401 128.144593
4 583.216472 108.616401
5 246.903315 583.216472
6 427.945333 246.903315
7 182.038538 427.945333
8 757.104136 182.038538
9 741.686165 757.104136
10 1007.910991 741.686165
11 841.999655 1007.910991
12 851.073547 841.999655
13 452.935303 851.073547
14 725.846337 452.935303
15 682.728837 725.846337
16 284.391861 682.728837
17 495.253061 284.391861
18 -17.522613 495.253061
19 865.015678 -17.522613
20 255.019182 865.015678
21 510.163361 255.019182
22 4.045776 510.163361
23 -362.874613 4.045776
24 -268.042382 -362.874613
25 262.344550 -268.042382
26 -708.551072 262.344550
27 -194.794371 -708.551072
28 -84.383303 -194.794371
29 -470.790101 -84.383303
30 -211.796346 -470.790101
31 -787.536416 -211.796346
32 -724.611115 -787.536416
33 -423.488092 -724.611115
34 -1210.791232 -423.488092
35 -481.607917 -1210.791232
36 -586.911809 -481.607917
37 -396.727778 -586.911809
38 -544.714279 -396.727778
39 42.312285 -544.714279
40 -460.023161 42.312285
41 -402.780221 -460.023161
42 -249.558909 -402.780221
43 -259.517800 -249.558909
44 -287.512202 -259.517800
45 -828.361434 -287.512202
46 198.834465 -828.361434
47 2.482875 198.834465
48 -124.052710 2.482875
49 -224.628852 -124.052710
50 399.274422 -224.628852
51 -638.863151 399.274422
52 -323.201870 -638.863151
53 131.413946 -323.201870
54 50.932535 131.413946
55 NA 50.932535
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -93.923224 127.933355
[2,] 128.144593 -93.923224
[3,] 108.616401 128.144593
[4,] 583.216472 108.616401
[5,] 246.903315 583.216472
[6,] 427.945333 246.903315
[7,] 182.038538 427.945333
[8,] 757.104136 182.038538
[9,] 741.686165 757.104136
[10,] 1007.910991 741.686165
[11,] 841.999655 1007.910991
[12,] 851.073547 841.999655
[13,] 452.935303 851.073547
[14,] 725.846337 452.935303
[15,] 682.728837 725.846337
[16,] 284.391861 682.728837
[17,] 495.253061 284.391861
[18,] -17.522613 495.253061
[19,] 865.015678 -17.522613
[20,] 255.019182 865.015678
[21,] 510.163361 255.019182
[22,] 4.045776 510.163361
[23,] -362.874613 4.045776
[24,] -268.042382 -362.874613
[25,] 262.344550 -268.042382
[26,] -708.551072 262.344550
[27,] -194.794371 -708.551072
[28,] -84.383303 -194.794371
[29,] -470.790101 -84.383303
[30,] -211.796346 -470.790101
[31,] -787.536416 -211.796346
[32,] -724.611115 -787.536416
[33,] -423.488092 -724.611115
[34,] -1210.791232 -423.488092
[35,] -481.607917 -1210.791232
[36,] -586.911809 -481.607917
[37,] -396.727778 -586.911809
[38,] -544.714279 -396.727778
[39,] 42.312285 -544.714279
[40,] -460.023161 42.312285
[41,] -402.780221 -460.023161
[42,] -249.558909 -402.780221
[43,] -259.517800 -249.558909
[44,] -287.512202 -259.517800
[45,] -828.361434 -287.512202
[46,] 198.834465 -828.361434
[47,] 2.482875 198.834465
[48,] -124.052710 2.482875
[49,] -224.628852 -124.052710
[50,] 399.274422 -224.628852
[51,] -638.863151 399.274422
[52,] -323.201870 -638.863151
[53,] 131.413946 -323.201870
[54,] 50.932535 131.413946
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -93.923224 127.933355
2 128.144593 -93.923224
3 108.616401 128.144593
4 583.216472 108.616401
5 246.903315 583.216472
6 427.945333 246.903315
7 182.038538 427.945333
8 757.104136 182.038538
9 741.686165 757.104136
10 1007.910991 741.686165
11 841.999655 1007.910991
12 851.073547 841.999655
13 452.935303 851.073547
14 725.846337 452.935303
15 682.728837 725.846337
16 284.391861 682.728837
17 495.253061 284.391861
18 -17.522613 495.253061
19 865.015678 -17.522613
20 255.019182 865.015678
21 510.163361 255.019182
22 4.045776 510.163361
23 -362.874613 4.045776
24 -268.042382 -362.874613
25 262.344550 -268.042382
26 -708.551072 262.344550
27 -194.794371 -708.551072
28 -84.383303 -194.794371
29 -470.790101 -84.383303
30 -211.796346 -470.790101
31 -787.536416 -211.796346
32 -724.611115 -787.536416
33 -423.488092 -724.611115
34 -1210.791232 -423.488092
35 -481.607917 -1210.791232
36 -586.911809 -481.607917
37 -396.727778 -586.911809
38 -544.714279 -396.727778
39 42.312285 -544.714279
40 -460.023161 42.312285
41 -402.780221 -460.023161
42 -249.558909 -402.780221
43 -259.517800 -249.558909
44 -287.512202 -259.517800
45 -828.361434 -287.512202
46 198.834465 -828.361434
47 2.482875 198.834465
48 -124.052710 2.482875
49 -224.628852 -124.052710
50 399.274422 -224.628852
51 -638.863151 399.274422
52 -323.201870 -638.863151
53 131.413946 -323.201870
54 50.932535 131.413946
> 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/7yfs81290170821.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/freestat/rcomp/tmp/8yfs81290170821.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/freestat/rcomp/tmp/9yfs81290170821.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/freestat/rcomp/tmp/10rprt1290170821.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/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/11c7qy1290170821.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/12x7om1290170821.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/13uzmv1290170821.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/14fi2j1290170821.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/15iijp1290170821.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/16eagx1290170821.tab")
+ }
>
> try(system("convert tmp/12nuh1290170821.ps tmp/12nuh1290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/22nuh1290170821.ps tmp/22nuh1290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/3vfu21290170821.ps tmp/3vfu21290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/4vfu21290170821.ps tmp/4vfu21290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/5vfu21290170821.ps tmp/5vfu21290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/6n6bn1290170821.ps tmp/6n6bn1290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/7yfs81290170821.ps tmp/7yfs81290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/8yfs81290170821.ps tmp/8yfs81290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/9yfs81290170821.ps tmp/9yfs81290170821.png",intern=TRUE))
character(0)
> try(system("convert tmp/10rprt1290170821.ps tmp/10rprt1290170821.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.817 2.487 20.607