R version 2.8.0 (2008-10-20)
Copyright (C) 2008 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> x <- array(list(577992,0,565464,0,547344,0,554788,0,562325,0,560854,0,555332,1,543599,1,536662,1,542722,1,593530,1,610763,1,612613,1,611324,1,594167,1,595454,1,590865,1,589379,1,584428,1,573100,1,567456,1,569028,1,620735,1,628884,1,628232,1,612117,1,595404,1,597141,1,593408,1,590072,1,579799,1,574205,1,572775,1,572942,1,619567,1,625809,1,619916,1,587625,1,565742,1,557274,1,560576,1,548854,1,531673,1,525919,1,511038,1,498662,1,555362,1,564591,1,541657,1,527070,1,509846,1,514258,1,516922,1,507561,1,492622,1,490243,1,469357,1,477580,1,528379,1,533590,1,517945,1),dim=c(2,61),dimnames=list(c('Werkloosheid','Aanslag'),1:61))
> y <- array(NA,dim=c(2,61),dimnames=list(c('Werkloosheid','Aanslag'),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
Werkloosheid Aanslag M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 577992 0 1 0 0 0 0 0 0 0 0 0 0
2 565464 0 0 1 0 0 0 0 0 0 0 0 0
3 547344 0 0 0 1 0 0 0 0 0 0 0 0
4 554788 0 0 0 0 1 0 0 0 0 0 0 0
5 562325 0 0 0 0 0 1 0 0 0 0 0 0
6 560854 0 0 0 0 0 0 1 0 0 0 0 0
7 555332 1 0 0 0 0 0 0 1 0 0 0 0
8 543599 1 0 0 0 0 0 0 0 1 0 0 0
9 536662 1 0 0 0 0 0 0 0 0 1 0 0
10 542722 1 0 0 0 0 0 0 0 0 0 1 0
11 593530 1 0 0 0 0 0 0 0 0 0 0 1
12 610763 1 0 0 0 0 0 0 0 0 0 0 0
13 612613 1 1 0 0 0 0 0 0 0 0 0 0
14 611324 1 0 1 0 0 0 0 0 0 0 0 0
15 594167 1 0 0 1 0 0 0 0 0 0 0 0
16 595454 1 0 0 0 1 0 0 0 0 0 0 0
17 590865 1 0 0 0 0 1 0 0 0 0 0 0
18 589379 1 0 0 0 0 0 1 0 0 0 0 0
19 584428 1 0 0 0 0 0 0 1 0 0 0 0
20 573100 1 0 0 0 0 0 0 0 1 0 0 0
21 567456 1 0 0 0 0 0 0 0 0 1 0 0
22 569028 1 0 0 0 0 0 0 0 0 0 1 0
23 620735 1 0 0 0 0 0 0 0 0 0 0 1
24 628884 1 0 0 0 0 0 0 0 0 0 0 0
25 628232 1 1 0 0 0 0 0 0 0 0 0 0
26 612117 1 0 1 0 0 0 0 0 0 0 0 0
27 595404 1 0 0 1 0 0 0 0 0 0 0 0
28 597141 1 0 0 0 1 0 0 0 0 0 0 0
29 593408 1 0 0 0 0 1 0 0 0 0 0 0
30 590072 1 0 0 0 0 0 1 0 0 0 0 0
31 579799 1 0 0 0 0 0 0 1 0 0 0 0
32 574205 1 0 0 0 0 0 0 0 1 0 0 0
33 572775 1 0 0 0 0 0 0 0 0 1 0 0
34 572942 1 0 0 0 0 0 0 0 0 0 1 0
35 619567 1 0 0 0 0 0 0 0 0 0 0 1
36 625809 1 0 0 0 0 0 0 0 0 0 0 0
37 619916 1 1 0 0 0 0 0 0 0 0 0 0
38 587625 1 0 1 0 0 0 0 0 0 0 0 0
39 565742 1 0 0 1 0 0 0 0 0 0 0 0
40 557274 1 0 0 0 1 0 0 0 0 0 0 0
41 560576 1 0 0 0 0 1 0 0 0 0 0 0
42 548854 1 0 0 0 0 0 1 0 0 0 0 0
43 531673 1 0 0 0 0 0 0 1 0 0 0 0
44 525919 1 0 0 0 0 0 0 0 1 0 0 0
45 511038 1 0 0 0 0 0 0 0 0 1 0 0
46 498662 1 0 0 0 0 0 0 0 0 0 1 0
47 555362 1 0 0 0 0 0 0 0 0 0 0 1
48 564591 1 0 0 0 0 0 0 0 0 0 0 0
49 541657 1 1 0 0 0 0 0 0 0 0 0 0
50 527070 1 0 1 0 0 0 0 0 0 0 0 0
51 509846 1 0 0 1 0 0 0 0 0 0 0 0
52 514258 1 0 0 0 1 0 0 0 0 0 0 0
53 516922 1 0 0 0 0 1 0 0 0 0 0 0
54 507561 1 0 0 0 0 0 1 0 0 0 0 0
55 492622 1 0 0 0 0 0 0 1 0 0 0 0
56 490243 1 0 0 0 0 0 0 0 1 0 0 0
57 469357 1 0 0 0 0 0 0 0 0 1 0 0
58 477580 1 0 0 0 0 0 0 0 0 0 1 0
59 528379 1 0 0 0 0 0 0 0 0 0 0 1
60 533590 1 0 0 0 0 0 0 0 0 0 0 0
61 517945 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) Aanslag M1 M2 M3 M4
583322 9405 -8101 -10126 -28346 -27063
M5 M6 M7 M8 M9 M10
-26027 -31502 -43957 -51314 -61270 -60541
M11
-9213
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-66682 -28136 5204 31022 43605
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 583322 24703 23.613 <2e-16 ***
Aanslag 9405 17615 0.534 0.5958
M1 -8101 23633 -0.343 0.7333
M2 -10126 24745 -0.409 0.6842
M3 -28346 24745 -1.146 0.2577
M4 -27063 24745 -1.094 0.2795
M5 -26027 24745 -1.052 0.2981
M6 -31502 24745 -1.273 0.2091
M7 -43957 24493 -1.795 0.0790 .
M8 -51314 24493 -2.095 0.0415 *
M9 -61270 24493 -2.502 0.0158 *
M10 -60541 24493 -2.472 0.0170 *
M11 -9213 24493 -0.376 0.7085
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 38730 on 48 degrees of freedom
Multiple R-squared: 0.2505, Adjusted R-squared: 0.06315
F-statistic: 1.337 on 12 and 48 DF, p-value: 0.2298
> 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,] 9.190185e-04 1.838037e-03 0.9990810
[2,] 4.122365e-04 8.244730e-04 0.9995878
[3,] 9.839585e-05 1.967917e-04 0.9999016
[4,] 4.560046e-04 9.120092e-04 0.9995440
[5,] 5.695529e-04 1.139106e-03 0.9994304
[6,] 5.937964e-04 1.187593e-03 0.9994062
[7,] 4.151583e-04 8.303165e-04 0.9995848
[8,] 3.029171e-04 6.058343e-04 0.9996971
[9,] 1.486407e-04 2.972813e-04 0.9998514
[10,] 7.927937e-05 1.585587e-04 0.9999207
[11,] 2.896854e-05 5.793709e-05 0.9999710
[12,] 1.096184e-05 2.192368e-05 0.9999890
[13,] 4.238479e-06 8.476959e-06 0.9999958
[14,] 1.546980e-06 3.093960e-06 0.9999985
[15,] 6.521916e-07 1.304383e-06 0.9999993
[16,] 3.738408e-07 7.476816e-07 0.9999996
[17,] 3.095361e-07 6.190722e-07 0.9999997
[18,] 6.415536e-07 1.283107e-06 0.9999994
[19,] 1.483812e-06 2.967624e-06 0.9999985
[20,] 2.931254e-06 5.862508e-06 0.9999971
[21,] 6.951193e-06 1.390239e-05 0.9999930
[22,] 9.042739e-05 1.808548e-04 0.9999096
[23,] 3.549593e-04 7.099186e-04 0.9996450
[24,] 1.608824e-03 3.217648e-03 0.9983912
[25,] 7.066234e-03 1.413247e-02 0.9929338
[26,] 1.936663e-02 3.873325e-02 0.9806334
[27,] 5.535694e-02 1.107139e-01 0.9446431
[28,] 1.115675e-01 2.231350e-01 0.8884325
[29,] 1.614752e-01 3.229503e-01 0.8385248
[30,] 3.066131e-01 6.132261e-01 0.6933869
> postscript(file="/var/www/html/rcomp/tmp/1i1m51227561606.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/2uoae1227561606.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/35okx1227561606.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/45e861227561606.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/5lizm1227561606.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 7
2770.586 -7731.757 -7632.357 -1470.757 5030.043 9034.243 6561.200
8 9 10 11 12 13 14
2185.800 5204.400 10535.200 10015.400 18035.600 27986.283 28722.939
15 16 17 18 19 20 21
29785.339 29789.939 24164.739 28153.939 35657.200 31686.800 35998.400
22 23 24 25 26 27 28
36841.200 37220.400 36156.600 43605.283 29515.939 31022.339 31476.939
29 30 31 32 33 34 35
26707.739 28846.939 31028.200 32791.800 41317.400 40755.200 36052.400
36 37 38 39 40 41 42
33081.600 35289.283 5023.939 1360.339 -8390.061 -6124.261 -12371.061
43 44 45 46 47 48 49
-17097.800 -15494.200 -20419.600 -33524.800 -28152.600 -28136.400 -42969.717
50 51 52 53 54 55 56
-55531.061 -54535.661 -51406.061 -49778.261 -53664.061 -56148.800 -51170.200
57 58 59 60 61
-62100.600 -54606.800 -55135.600 -59137.400 -66681.717
> postscript(file="/var/www/html/rcomp/tmp/67asp1227561606.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 2770.586 NA
1 -7731.757 2770.586
2 -7632.357 -7731.757
3 -1470.757 -7632.357
4 5030.043 -1470.757
5 9034.243 5030.043
6 6561.200 9034.243
7 2185.800 6561.200
8 5204.400 2185.800
9 10535.200 5204.400
10 10015.400 10535.200
11 18035.600 10015.400
12 27986.283 18035.600
13 28722.939 27986.283
14 29785.339 28722.939
15 29789.939 29785.339
16 24164.739 29789.939
17 28153.939 24164.739
18 35657.200 28153.939
19 31686.800 35657.200
20 35998.400 31686.800
21 36841.200 35998.400
22 37220.400 36841.200
23 36156.600 37220.400
24 43605.283 36156.600
25 29515.939 43605.283
26 31022.339 29515.939
27 31476.939 31022.339
28 26707.739 31476.939
29 28846.939 26707.739
30 31028.200 28846.939
31 32791.800 31028.200
32 41317.400 32791.800
33 40755.200 41317.400
34 36052.400 40755.200
35 33081.600 36052.400
36 35289.283 33081.600
37 5023.939 35289.283
38 1360.339 5023.939
39 -8390.061 1360.339
40 -6124.261 -8390.061
41 -12371.061 -6124.261
42 -17097.800 -12371.061
43 -15494.200 -17097.800
44 -20419.600 -15494.200
45 -33524.800 -20419.600
46 -28152.600 -33524.800
47 -28136.400 -28152.600
48 -42969.717 -28136.400
49 -55531.061 -42969.717
50 -54535.661 -55531.061
51 -51406.061 -54535.661
52 -49778.261 -51406.061
53 -53664.061 -49778.261
54 -56148.800 -53664.061
55 -51170.200 -56148.800
56 -62100.600 -51170.200
57 -54606.800 -62100.600
58 -55135.600 -54606.800
59 -59137.400 -55135.600
60 -66681.717 -59137.400
61 NA -66681.717
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -7731.757 2770.586
[2,] -7632.357 -7731.757
[3,] -1470.757 -7632.357
[4,] 5030.043 -1470.757
[5,] 9034.243 5030.043
[6,] 6561.200 9034.243
[7,] 2185.800 6561.200
[8,] 5204.400 2185.800
[9,] 10535.200 5204.400
[10,] 10015.400 10535.200
[11,] 18035.600 10015.400
[12,] 27986.283 18035.600
[13,] 28722.939 27986.283
[14,] 29785.339 28722.939
[15,] 29789.939 29785.339
[16,] 24164.739 29789.939
[17,] 28153.939 24164.739
[18,] 35657.200 28153.939
[19,] 31686.800 35657.200
[20,] 35998.400 31686.800
[21,] 36841.200 35998.400
[22,] 37220.400 36841.200
[23,] 36156.600 37220.400
[24,] 43605.283 36156.600
[25,] 29515.939 43605.283
[26,] 31022.339 29515.939
[27,] 31476.939 31022.339
[28,] 26707.739 31476.939
[29,] 28846.939 26707.739
[30,] 31028.200 28846.939
[31,] 32791.800 31028.200
[32,] 41317.400 32791.800
[33,] 40755.200 41317.400
[34,] 36052.400 40755.200
[35,] 33081.600 36052.400
[36,] 35289.283 33081.600
[37,] 5023.939 35289.283
[38,] 1360.339 5023.939
[39,] -8390.061 1360.339
[40,] -6124.261 -8390.061
[41,] -12371.061 -6124.261
[42,] -17097.800 -12371.061
[43,] -15494.200 -17097.800
[44,] -20419.600 -15494.200
[45,] -33524.800 -20419.600
[46,] -28152.600 -33524.800
[47,] -28136.400 -28152.600
[48,] -42969.717 -28136.400
[49,] -55531.061 -42969.717
[50,] -54535.661 -55531.061
[51,] -51406.061 -54535.661
[52,] -49778.261 -51406.061
[53,] -53664.061 -49778.261
[54,] -56148.800 -53664.061
[55,] -51170.200 -56148.800
[56,] -62100.600 -51170.200
[57,] -54606.800 -62100.600
[58,] -55135.600 -54606.800
[59,] -59137.400 -55135.600
[60,] -66681.717 -59137.400
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -7731.757 2770.586
2 -7632.357 -7731.757
3 -1470.757 -7632.357
4 5030.043 -1470.757
5 9034.243 5030.043
6 6561.200 9034.243
7 2185.800 6561.200
8 5204.400 2185.800
9 10535.200 5204.400
10 10015.400 10535.200
11 18035.600 10015.400
12 27986.283 18035.600
13 28722.939 27986.283
14 29785.339 28722.939
15 29789.939 29785.339
16 24164.739 29789.939
17 28153.939 24164.739
18 35657.200 28153.939
19 31686.800 35657.200
20 35998.400 31686.800
21 36841.200 35998.400
22 37220.400 36841.200
23 36156.600 37220.400
24 43605.283 36156.600
25 29515.939 43605.283
26 31022.339 29515.939
27 31476.939 31022.339
28 26707.739 31476.939
29 28846.939 26707.739
30 31028.200 28846.939
31 32791.800 31028.200
32 41317.400 32791.800
33 40755.200 41317.400
34 36052.400 40755.200
35 33081.600 36052.400
36 35289.283 33081.600
37 5023.939 35289.283
38 1360.339 5023.939
39 -8390.061 1360.339
40 -6124.261 -8390.061
41 -12371.061 -6124.261
42 -17097.800 -12371.061
43 -15494.200 -17097.800
44 -20419.600 -15494.200
45 -33524.800 -20419.600
46 -28152.600 -33524.800
47 -28136.400 -28152.600
48 -42969.717 -28136.400
49 -55531.061 -42969.717
50 -54535.661 -55531.061
51 -51406.061 -54535.661
52 -49778.261 -51406.061
53 -53664.061 -49778.261
54 -56148.800 -53664.061
55 -51170.200 -56148.800
56 -62100.600 -51170.200
57 -54606.800 -62100.600
58 -55135.600 -54606.800
59 -59137.400 -55135.600
60 -66681.717 -59137.400
> 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/7ddw21227561606.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/8i1kl1227561606.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/9l5ts1227561606.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/10nu3v1227561606.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/11mo2w1227561606.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/129msj1227561606.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/13ozpu1227561606.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/14fmyo1227561606.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/15k3il1227561606.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/16gsoc1227561606.tab")
+ }
>
> system("convert tmp/1i1m51227561606.ps tmp/1i1m51227561606.png")
> system("convert tmp/2uoae1227561606.ps tmp/2uoae1227561606.png")
> system("convert tmp/35okx1227561606.ps tmp/35okx1227561606.png")
> system("convert tmp/45e861227561606.ps tmp/45e861227561606.png")
> system("convert tmp/5lizm1227561606.ps tmp/5lizm1227561606.png")
> system("convert tmp/67asp1227561606.ps tmp/67asp1227561606.png")
> system("convert tmp/7ddw21227561606.ps tmp/7ddw21227561606.png")
> system("convert tmp/8i1kl1227561606.ps tmp/8i1kl1227561606.png")
> system("convert tmp/9l5ts1227561606.ps tmp/9l5ts1227561606.png")
> system("convert tmp/10nu3v1227561606.ps tmp/10nu3v1227561606.png")
>
>
> proc.time()
user system elapsed
2.463 1.604 4.460