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(31.514
+ ,-9
+ ,0
+ ,8.3
+ ,1.2
+ ,27.071
+ ,-13
+ ,4
+ ,8.2
+ ,1.7
+ ,29.462
+ ,-18
+ ,5
+ ,8
+ ,1.8
+ ,26.105
+ ,-11
+ ,-7
+ ,7.9
+ ,1.5
+ ,22.397
+ ,-9
+ ,-2
+ ,7.6
+ ,1
+ ,23.843
+ ,-10
+ ,1
+ ,7.6
+ ,1.6
+ ,21.705
+ ,-13
+ ,3
+ ,8.3
+ ,1.5
+ ,18.089
+ ,-11
+ ,-2
+ ,8.4
+ ,1.8
+ ,20.764
+ ,-5
+ ,-6
+ ,8.4
+ ,1.8
+ ,25.316
+ ,-15
+ ,10
+ ,8.4
+ ,1.6
+ ,17.704
+ ,-6
+ ,-9
+ ,8.4
+ ,1.9
+ ,15.548
+ ,-6
+ ,0
+ ,8.6
+ ,1.7
+ ,28.029
+ ,-3
+ ,-3
+ ,8.9
+ ,1.6
+ ,29.383
+ ,-1
+ ,-2
+ ,8.8
+ ,1.3
+ ,36.438
+ ,-3
+ ,2
+ ,8.3
+ ,1.1
+ ,32.034
+ ,-4
+ ,1
+ ,7.5
+ ,1.9
+ ,22.679
+ ,-6
+ ,2
+ ,7.2
+ ,2.6
+ ,24.319
+ ,0
+ ,-6
+ ,7.4
+ ,2.3
+ ,18.004
+ ,-4
+ ,4
+ ,8.8
+ ,2.4
+ ,17.537
+ ,-2
+ ,-2
+ ,9.3
+ ,2.2
+ ,20.366
+ ,-2
+ ,0
+ ,9.3
+ ,2
+ ,22.782
+ ,-6
+ ,4
+ ,8.7
+ ,2.9
+ ,19.169
+ ,-7
+ ,1
+ ,8.2
+ ,2.6
+ ,13.807
+ ,-6
+ ,-1
+ ,8.3
+ ,2.3
+ ,29.743
+ ,-6
+ ,0
+ ,8.5
+ ,2.3
+ ,25.591
+ ,-3
+ ,-3
+ ,8.6
+ ,2.6
+ ,29.096
+ ,-2
+ ,-1
+ ,8.5
+ ,3.1
+ ,26.482
+ ,-5
+ ,3
+ ,8.2
+ ,2.8
+ ,22.405
+ ,-11
+ ,6
+ ,8.1
+ ,2.5
+ ,27.044
+ ,-11
+ ,0
+ ,7.9
+ ,2.9
+ ,17.970
+ ,-11
+ ,0
+ ,8.6
+ ,3.1
+ ,18.730
+ ,-10
+ ,-1
+ ,8.7
+ ,3.1
+ ,19.684
+ ,-14
+ ,4
+ ,8.7
+ ,3.2
+ ,19.785
+ ,-8
+ ,-6
+ ,8.5
+ ,2.5
+ ,18.479
+ ,-9
+ ,1
+ ,8.4
+ ,2.6
+ ,10.698
+ ,-5
+ ,-4
+ ,8.5
+ ,2.9
+ ,31.956
+ ,-1
+ ,-4
+ ,8.7
+ ,2.6
+ ,29.506
+ ,-2
+ ,1
+ ,8.7
+ ,2.4
+ ,34.506
+ ,-5
+ ,3
+ ,8.6
+ ,1.7
+ ,27.165
+ ,-4
+ ,-1
+ ,8.5
+ ,2
+ ,26.736
+ ,-6
+ ,2
+ ,8.3
+ ,2.2
+ ,23.691
+ ,-2
+ ,-4
+ ,8
+ ,1.9
+ ,18.157
+ ,-2
+ ,0
+ ,8.2
+ ,1.6
+ ,17.328
+ ,-2
+ ,0
+ ,8.1
+ ,1.6
+ ,18.205
+ ,-2
+ ,0
+ ,8.1
+ ,1.2
+ ,20.995
+ ,2
+ ,-4
+ ,8
+ ,1.2
+ ,17.382
+ ,1
+ ,1
+ ,7.9
+ ,1.5
+ ,9.367
+ ,-8
+ ,9
+ ,7.9
+ ,1.6
+ ,31.124
+ ,-1
+ ,-7
+ ,8
+ ,1.7
+ ,26.551
+ ,1
+ ,-2
+ ,8
+ ,1.8
+ ,30.651
+ ,-1
+ ,2
+ ,7.9
+ ,1.8
+ ,25.859
+ ,2
+ ,-3
+ ,8
+ ,1.8
+ ,25.100
+ ,2
+ ,0
+ ,7.7
+ ,1.3
+ ,25.778
+ ,1
+ ,1
+ ,7.2
+ ,1.3
+ ,20.418
+ ,-1
+ ,2
+ ,7.5
+ ,1.4
+ ,18.688
+ ,-2
+ ,1
+ ,7.3
+ ,1.1
+ ,20.424
+ ,-2
+ ,0
+ ,7
+ ,1.5
+ ,24.776
+ ,-1
+ ,-1
+ ,7
+ ,2.2
+ ,19.814
+ ,-8
+ ,7
+ ,7
+ ,2.9
+ ,12.738
+ ,-4
+ ,-4
+ ,7.2
+ ,3.1
+ ,31.566
+ ,-6
+ ,2
+ ,7.3
+ ,3.5
+ ,30.111
+ ,-3
+ ,-3
+ ,7.1
+ ,3.6
+ ,30.019
+ ,-3
+ ,0
+ ,6.8
+ ,4.4
+ ,31.934
+ ,-7
+ ,4
+ ,6.4
+ ,4.2
+ ,25.826
+ ,-9
+ ,2
+ ,6.1
+ ,5.2
+ ,26.835
+ ,-11
+ ,2
+ ,6.5
+ ,5.8
+ ,20.205
+ ,-13
+ ,2
+ ,7.7
+ ,5.9
+ ,17.789
+ ,-11
+ ,-2
+ ,7.9
+ ,5.4
+ ,20.520
+ ,-9
+ ,-2
+ ,7.5
+ ,5.5
+ ,22.518
+ ,-17
+ ,8
+ ,6.9
+ ,4.7
+ ,15.572
+ ,-22
+ ,5
+ ,6.6
+ ,3.1
+ ,11.509
+ ,-25
+ ,3
+ ,6.9
+ ,2.6
+ ,25.447
+ ,-20
+ ,-5
+ ,7.7
+ ,2.3
+ ,24.090
+ ,-24
+ ,4
+ ,8
+ ,1.9
+ ,27.786
+ ,-24
+ ,0
+ ,8
+ ,0.6
+ ,26.195
+ ,-22
+ ,-2
+ ,7.7
+ ,0.6
+ ,20.516
+ ,-19
+ ,-3
+ ,7.3
+ ,-0.4
+ ,22.759
+ ,-18
+ ,-1
+ ,7.4
+ ,-1.1
+ ,19.028
+ ,-17
+ ,-1
+ ,8.1
+ ,-1.7
+ ,16.971
+ ,-11
+ ,-6
+ ,8.3
+ ,-0.8
+ ,20.036
+ ,-11
+ ,0
+ ,8.1
+ ,-1.2
+ ,22.485
+ ,-12
+ ,1
+ ,7.9
+ ,-1
+ ,18.730
+ ,-10
+ ,-2
+ ,7.9
+ ,-0.1
+ ,14.538
+ ,-15
+ ,5
+ ,8.3
+ ,0.3)
+ ,dim=c(5
+ ,84)
+ ,dimnames=list(c('Inschrijvingen'
+ ,'Consumentenvertrouwen'
+ ,'Evolutie_consumentenvertrouwen'
+ ,'Totaal_Werkloosheid'
+ ,'Algemene_index')
+ ,1:84))
> y <- array(NA,dim=c(5,84),dimnames=list(c('Inschrijvingen','Consumentenvertrouwen','Evolutie_consumentenvertrouwen','Totaal_Werkloosheid','Algemene_index'),1:84))
> 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 = 'Do not include Seasonal 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
Inschrijvingen Consumentenvertrouwen Evolutie_consumentenvertrouwen
1 31.514 -9 0
2 27.071 -13 4
3 29.462 -18 5
4 26.105 -11 -7
5 22.397 -9 -2
6 23.843 -10 1
7 21.705 -13 3
8 18.089 -11 -2
9 20.764 -5 -6
10 25.316 -15 10
11 17.704 -6 -9
12 15.548 -6 0
13 28.029 -3 -3
14 29.383 -1 -2
15 36.438 -3 2
16 32.034 -4 1
17 22.679 -6 2
18 24.319 0 -6
19 18.004 -4 4
20 17.537 -2 -2
21 20.366 -2 0
22 22.782 -6 4
23 19.169 -7 1
24 13.807 -6 -1
25 29.743 -6 0
26 25.591 -3 -3
27 29.096 -2 -1
28 26.482 -5 3
29 22.405 -11 6
30 27.044 -11 0
31 17.970 -11 0
32 18.730 -10 -1
33 19.684 -14 4
34 19.785 -8 -6
35 18.479 -9 1
36 10.698 -5 -4
37 31.956 -1 -4
38 29.506 -2 1
39 34.506 -5 3
40 27.165 -4 -1
41 26.736 -6 2
42 23.691 -2 -4
43 18.157 -2 0
44 17.328 -2 0
45 18.205 -2 0
46 20.995 2 -4
47 17.382 1 1
48 9.367 -8 9
49 31.124 -1 -7
50 26.551 1 -2
51 30.651 -1 2
52 25.859 2 -3
53 25.100 2 0
54 25.778 1 1
55 20.418 -1 2
56 18.688 -2 1
57 20.424 -2 0
58 24.776 -1 -1
59 19.814 -8 7
60 12.738 -4 -4
61 31.566 -6 2
62 30.111 -3 -3
63 30.019 -3 0
64 31.934 -7 4
65 25.826 -9 2
66 26.835 -11 2
67 20.205 -13 2
68 17.789 -11 -2
69 20.520 -9 -2
70 22.518 -17 8
71 15.572 -22 5
72 11.509 -25 3
73 25.447 -20 -5
74 24.090 -24 4
75 27.786 -24 0
76 26.195 -22 -2
77 20.516 -19 -3
78 22.759 -18 -1
79 19.028 -17 -1
80 16.971 -11 -6
81 20.036 -11 0
82 22.485 -12 1
83 18.730 -10 -2
84 14.538 -15 5
Totaal_Werkloosheid Algemene_index
1 8.3 1.2
2 8.2 1.7
3 8.0 1.8
4 7.9 1.5
5 7.6 1.0
6 7.6 1.6
7 8.3 1.5
8 8.4 1.8
9 8.4 1.8
10 8.4 1.6
11 8.4 1.9
12 8.6 1.7
13 8.9 1.6
14 8.8 1.3
15 8.3 1.1
16 7.5 1.9
17 7.2 2.6
18 7.4 2.3
19 8.8 2.4
20 9.3 2.2
21 9.3 2.0
22 8.7 2.9
23 8.2 2.6
24 8.3 2.3
25 8.5 2.3
26 8.6 2.6
27 8.5 3.1
28 8.2 2.8
29 8.1 2.5
30 7.9 2.9
31 8.6 3.1
32 8.7 3.1
33 8.7 3.2
34 8.5 2.5
35 8.4 2.6
36 8.5 2.9
37 8.7 2.6
38 8.7 2.4
39 8.6 1.7
40 8.5 2.0
41 8.3 2.2
42 8.0 1.9
43 8.2 1.6
44 8.1 1.6
45 8.1 1.2
46 8.0 1.2
47 7.9 1.5
48 7.9 1.6
49 8.0 1.7
50 8.0 1.8
51 7.9 1.8
52 8.0 1.8
53 7.7 1.3
54 7.2 1.3
55 7.5 1.4
56 7.3 1.1
57 7.0 1.5
58 7.0 2.2
59 7.0 2.9
60 7.2 3.1
61 7.3 3.5
62 7.1 3.6
63 6.8 4.4
64 6.4 4.2
65 6.1 5.2
66 6.5 5.8
67 7.7 5.9
68 7.9 5.4
69 7.5 5.5
70 6.9 4.7
71 6.6 3.1
72 6.9 2.6
73 7.7 2.3
74 8.0 1.9
75 8.0 0.6
76 7.7 0.6
77 7.3 -0.4
78 7.4 -1.1
79 8.1 -1.7
80 8.3 -0.8
81 8.1 -1.2
82 7.9 -1.0
83 7.9 -0.1
84 8.3 0.3
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Consumentenvertrouwen
29.19230 0.17439
Evolutie_consumentenvertrouwen Totaal_Werkloosheid
0.04288 -0.64924
Algemene_index
0.13526
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-13.9035 -4.0794 -0.1843 4.5413 12.9231
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 29.19230 8.66602 3.369 0.00117 **
Consumentenvertrouwen 0.17439 0.10142 1.720 0.08944 .
Evolutie_consumentenvertrouwen 0.04288 0.18436 0.233 0.81668
Totaal_Werkloosheid -0.64924 1.02448 -0.634 0.52809
Algemene_index 0.13526 0.45834 0.295 0.76869
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.721 on 79 degrees of freedom
Multiple R-squared: 0.04352, Adjusted R-squared: -0.004913
F-statistic: 0.8986 on 4 and 79 DF, p-value: 0.469
> 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.6269545 0.746091052 0.373045526
[2,] 0.4767619 0.953523853 0.523238073
[3,] 0.3478297 0.695659311 0.652170345
[4,] 0.2433457 0.486691435 0.756654283
[5,] 0.2076913 0.415382649 0.792308676
[6,] 0.3363623 0.672724685 0.663637658
[7,] 0.3119559 0.623911848 0.688044076
[8,] 0.4327049 0.865409864 0.567295068
[9,] 0.4698906 0.939781243 0.530109378
[10,] 0.3851097 0.770219326 0.614890337
[11,] 0.3027959 0.605591779 0.697204110
[12,] 0.2647291 0.529458283 0.735270858
[13,] 0.2078290 0.415657993 0.792171004
[14,] 0.1533837 0.306767495 0.846616252
[15,] 0.1516243 0.303248670 0.848375665
[16,] 0.1117517 0.223503484 0.888248258
[17,] 0.1525313 0.305062500 0.847468750
[18,] 0.2545993 0.509198529 0.745400736
[19,] 0.2733832 0.546766301 0.726616850
[20,] 0.3759066 0.751813196 0.624093402
[21,] 0.3315632 0.663126447 0.668436777
[22,] 0.2735452 0.547090394 0.726454803
[23,] 0.2863890 0.572778018 0.713610991
[24,] 0.2359541 0.471908168 0.764045916
[25,] 0.1910760 0.382151923 0.808924039
[26,] 0.1484516 0.296903113 0.851548443
[27,] 0.1175188 0.235037510 0.882481245
[28,] 0.0991478 0.198295608 0.900852196
[29,] 0.2265287 0.453057354 0.773471323
[30,] 0.3310775 0.662154943 0.668922529
[31,] 0.3281952 0.656390455 0.671804772
[32,] 0.5120217 0.975956607 0.487978303
[33,] 0.4855291 0.971058282 0.514470859
[34,] 0.4747756 0.949551300 0.525224350
[35,] 0.4157915 0.831583050 0.584208475
[36,] 0.4792359 0.958471750 0.520764125
[37,] 0.5455025 0.908994979 0.454497489
[38,] 0.5796577 0.840684674 0.420342337
[39,] 0.5525757 0.894848629 0.447424314
[40,] 0.5952318 0.809536420 0.404768210
[41,] 0.8278437 0.344312612 0.172156306
[42,] 0.8446667 0.310666538 0.155333269
[43,] 0.8087042 0.382591503 0.191295752
[44,] 0.8366036 0.326792892 0.163396446
[45,] 0.7994865 0.401027078 0.200513539
[46,] 0.7550399 0.489920233 0.244960116
[47,] 0.7068094 0.586381148 0.293190574
[48,] 0.6601572 0.679685531 0.339842765
[49,] 0.6410032 0.717993570 0.358996785
[50,] 0.6011897 0.797620689 0.398810345
[51,] 0.5302975 0.939404965 0.469702483
[52,] 0.4825526 0.965105178 0.517447411
[53,] 0.7225920 0.554815983 0.277407992
[54,] 0.7888863 0.422227489 0.211113745
[55,] 0.7763643 0.447271359 0.223635680
[56,] 0.7658221 0.468355705 0.234177853
[57,] 0.8541577 0.291684651 0.145842325
[58,] 0.8337613 0.332477384 0.166238692
[59,] 0.8968196 0.206360880 0.103180440
[60,] 0.8507571 0.298485704 0.149242852
[61,] 0.8213211 0.357357750 0.178678875
[62,] 0.7487734 0.502453286 0.251226643
[63,] 0.8439505 0.312099025 0.156049512
[64,] 0.7920024 0.415995270 0.207997635
[65,] 0.9588484 0.082303191 0.041151596
[66,] 0.9210717 0.157856650 0.078928325
[67,] 0.8522438 0.295512320 0.147756160
[68,] 0.8427263 0.314547324 0.157273662
[69,] 0.9969254 0.006149255 0.003074627
> postscript(file="/var/www/html/rcomp/tmp/17d7p1292676443.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/rcomp/tmp/27d7p1292676443.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/rcomp/tmp/3hnoa1292676443.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/rcomp/tmp/4hnoa1292676443.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/rcomp/tmp/5hnoa1292676443.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 = 84
Frequency = 1
1 2 3 4 5 6
9.11762325 5.06811365 8.14480781 4.05727605 -0.34103941 1.06955965
7 8 9 10 11 12
-0.16303174 -3.88907144 -2.08889446 3.54799661 -4.85939544 -7.24440098
13 14 15 16 17 18
5.05036521 5.98836084 12.92305530 8.10872082 -1.22983295 -0.12271568
19 20 21 22 23 24
-5.27352571 -5.48035990 -2.71006514 -0.27930041 -3.87332036 -9.21845094
25 26 27 28 29 30
6.80451953 2.28233343 5.39463338 2.97809209 -0.20555339 4.50676509
31 32 33 34 35 36
-4.13981533 -3.44640188 -2.02276211 -2.57448243 -4.08469248 -12.32451113
37 38 39 40 41 42
8.40635741 5.94340625 11.41057416 3.96119655 3.59543952 -0.04404351
43 44 45 46 47 48
-5.57913091 -6.47305537 -5.54195202 -3.34292055 -7.10142535 -13.90347359
49 50 51 52 53 54
7.37025407 2.22055698 6.43289763 1.39704596 0.38226638 0.86715507
55 56 57 58 59 60
-4.00569688 -5.60770036 -4.07769865 0.04810947 -4.13087274 -11.32997033
61 62 63 64 65 66
7.60035897 5.69320808 5.16959260 7.37799048 1.37449458 2.91081636
67 68 69 70 71 72
-2.60483694 -5.00062393 -2.89162658 -0.20863547 -6.13241269 -9.32408475
73 74 75 76 77 78
4.64496878 3.84849729 7.89184702 5.84305160 -0.44067785 1.70178108
79 80 81 82 83 84
-1.66798211 -4.54881028 -1.81682662 0.60678379 -3.49009232 -6.90469966
> postscript(file="/var/www/html/rcomp/tmp/6aw6d1292676443.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 = 84
Frequency = 1
lag(myerror, k = 1) myerror
0 9.11762325 NA
1 5.06811365 9.11762325
2 8.14480781 5.06811365
3 4.05727605 8.14480781
4 -0.34103941 4.05727605
5 1.06955965 -0.34103941
6 -0.16303174 1.06955965
7 -3.88907144 -0.16303174
8 -2.08889446 -3.88907144
9 3.54799661 -2.08889446
10 -4.85939544 3.54799661
11 -7.24440098 -4.85939544
12 5.05036521 -7.24440098
13 5.98836084 5.05036521
14 12.92305530 5.98836084
15 8.10872082 12.92305530
16 -1.22983295 8.10872082
17 -0.12271568 -1.22983295
18 -5.27352571 -0.12271568
19 -5.48035990 -5.27352571
20 -2.71006514 -5.48035990
21 -0.27930041 -2.71006514
22 -3.87332036 -0.27930041
23 -9.21845094 -3.87332036
24 6.80451953 -9.21845094
25 2.28233343 6.80451953
26 5.39463338 2.28233343
27 2.97809209 5.39463338
28 -0.20555339 2.97809209
29 4.50676509 -0.20555339
30 -4.13981533 4.50676509
31 -3.44640188 -4.13981533
32 -2.02276211 -3.44640188
33 -2.57448243 -2.02276211
34 -4.08469248 -2.57448243
35 -12.32451113 -4.08469248
36 8.40635741 -12.32451113
37 5.94340625 8.40635741
38 11.41057416 5.94340625
39 3.96119655 11.41057416
40 3.59543952 3.96119655
41 -0.04404351 3.59543952
42 -5.57913091 -0.04404351
43 -6.47305537 -5.57913091
44 -5.54195202 -6.47305537
45 -3.34292055 -5.54195202
46 -7.10142535 -3.34292055
47 -13.90347359 -7.10142535
48 7.37025407 -13.90347359
49 2.22055698 7.37025407
50 6.43289763 2.22055698
51 1.39704596 6.43289763
52 0.38226638 1.39704596
53 0.86715507 0.38226638
54 -4.00569688 0.86715507
55 -5.60770036 -4.00569688
56 -4.07769865 -5.60770036
57 0.04810947 -4.07769865
58 -4.13087274 0.04810947
59 -11.32997033 -4.13087274
60 7.60035897 -11.32997033
61 5.69320808 7.60035897
62 5.16959260 5.69320808
63 7.37799048 5.16959260
64 1.37449458 7.37799048
65 2.91081636 1.37449458
66 -2.60483694 2.91081636
67 -5.00062393 -2.60483694
68 -2.89162658 -5.00062393
69 -0.20863547 -2.89162658
70 -6.13241269 -0.20863547
71 -9.32408475 -6.13241269
72 4.64496878 -9.32408475
73 3.84849729 4.64496878
74 7.89184702 3.84849729
75 5.84305160 7.89184702
76 -0.44067785 5.84305160
77 1.70178108 -0.44067785
78 -1.66798211 1.70178108
79 -4.54881028 -1.66798211
80 -1.81682662 -4.54881028
81 0.60678379 -1.81682662
82 -3.49009232 0.60678379
83 -6.90469966 -3.49009232
84 NA -6.90469966
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 5.06811365 9.11762325
[2,] 8.14480781 5.06811365
[3,] 4.05727605 8.14480781
[4,] -0.34103941 4.05727605
[5,] 1.06955965 -0.34103941
[6,] -0.16303174 1.06955965
[7,] -3.88907144 -0.16303174
[8,] -2.08889446 -3.88907144
[9,] 3.54799661 -2.08889446
[10,] -4.85939544 3.54799661
[11,] -7.24440098 -4.85939544
[12,] 5.05036521 -7.24440098
[13,] 5.98836084 5.05036521
[14,] 12.92305530 5.98836084
[15,] 8.10872082 12.92305530
[16,] -1.22983295 8.10872082
[17,] -0.12271568 -1.22983295
[18,] -5.27352571 -0.12271568
[19,] -5.48035990 -5.27352571
[20,] -2.71006514 -5.48035990
[21,] -0.27930041 -2.71006514
[22,] -3.87332036 -0.27930041
[23,] -9.21845094 -3.87332036
[24,] 6.80451953 -9.21845094
[25,] 2.28233343 6.80451953
[26,] 5.39463338 2.28233343
[27,] 2.97809209 5.39463338
[28,] -0.20555339 2.97809209
[29,] 4.50676509 -0.20555339
[30,] -4.13981533 4.50676509
[31,] -3.44640188 -4.13981533
[32,] -2.02276211 -3.44640188
[33,] -2.57448243 -2.02276211
[34,] -4.08469248 -2.57448243
[35,] -12.32451113 -4.08469248
[36,] 8.40635741 -12.32451113
[37,] 5.94340625 8.40635741
[38,] 11.41057416 5.94340625
[39,] 3.96119655 11.41057416
[40,] 3.59543952 3.96119655
[41,] -0.04404351 3.59543952
[42,] -5.57913091 -0.04404351
[43,] -6.47305537 -5.57913091
[44,] -5.54195202 -6.47305537
[45,] -3.34292055 -5.54195202
[46,] -7.10142535 -3.34292055
[47,] -13.90347359 -7.10142535
[48,] 7.37025407 -13.90347359
[49,] 2.22055698 7.37025407
[50,] 6.43289763 2.22055698
[51,] 1.39704596 6.43289763
[52,] 0.38226638 1.39704596
[53,] 0.86715507 0.38226638
[54,] -4.00569688 0.86715507
[55,] -5.60770036 -4.00569688
[56,] -4.07769865 -5.60770036
[57,] 0.04810947 -4.07769865
[58,] -4.13087274 0.04810947
[59,] -11.32997033 -4.13087274
[60,] 7.60035897 -11.32997033
[61,] 5.69320808 7.60035897
[62,] 5.16959260 5.69320808
[63,] 7.37799048 5.16959260
[64,] 1.37449458 7.37799048
[65,] 2.91081636 1.37449458
[66,] -2.60483694 2.91081636
[67,] -5.00062393 -2.60483694
[68,] -2.89162658 -5.00062393
[69,] -0.20863547 -2.89162658
[70,] -6.13241269 -0.20863547
[71,] -9.32408475 -6.13241269
[72,] 4.64496878 -9.32408475
[73,] 3.84849729 4.64496878
[74,] 7.89184702 3.84849729
[75,] 5.84305160 7.89184702
[76,] -0.44067785 5.84305160
[77,] 1.70178108 -0.44067785
[78,] -1.66798211 1.70178108
[79,] -4.54881028 -1.66798211
[80,] -1.81682662 -4.54881028
[81,] 0.60678379 -1.81682662
[82,] -3.49009232 0.60678379
[83,] -6.90469966 -3.49009232
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 5.06811365 9.11762325
2 8.14480781 5.06811365
3 4.05727605 8.14480781
4 -0.34103941 4.05727605
5 1.06955965 -0.34103941
6 -0.16303174 1.06955965
7 -3.88907144 -0.16303174
8 -2.08889446 -3.88907144
9 3.54799661 -2.08889446
10 -4.85939544 3.54799661
11 -7.24440098 -4.85939544
12 5.05036521 -7.24440098
13 5.98836084 5.05036521
14 12.92305530 5.98836084
15 8.10872082 12.92305530
16 -1.22983295 8.10872082
17 -0.12271568 -1.22983295
18 -5.27352571 -0.12271568
19 -5.48035990 -5.27352571
20 -2.71006514 -5.48035990
21 -0.27930041 -2.71006514
22 -3.87332036 -0.27930041
23 -9.21845094 -3.87332036
24 6.80451953 -9.21845094
25 2.28233343 6.80451953
26 5.39463338 2.28233343
27 2.97809209 5.39463338
28 -0.20555339 2.97809209
29 4.50676509 -0.20555339
30 -4.13981533 4.50676509
31 -3.44640188 -4.13981533
32 -2.02276211 -3.44640188
33 -2.57448243 -2.02276211
34 -4.08469248 -2.57448243
35 -12.32451113 -4.08469248
36 8.40635741 -12.32451113
37 5.94340625 8.40635741
38 11.41057416 5.94340625
39 3.96119655 11.41057416
40 3.59543952 3.96119655
41 -0.04404351 3.59543952
42 -5.57913091 -0.04404351
43 -6.47305537 -5.57913091
44 -5.54195202 -6.47305537
45 -3.34292055 -5.54195202
46 -7.10142535 -3.34292055
47 -13.90347359 -7.10142535
48 7.37025407 -13.90347359
49 2.22055698 7.37025407
50 6.43289763 2.22055698
51 1.39704596 6.43289763
52 0.38226638 1.39704596
53 0.86715507 0.38226638
54 -4.00569688 0.86715507
55 -5.60770036 -4.00569688
56 -4.07769865 -5.60770036
57 0.04810947 -4.07769865
58 -4.13087274 0.04810947
59 -11.32997033 -4.13087274
60 7.60035897 -11.32997033
61 5.69320808 7.60035897
62 5.16959260 5.69320808
63 7.37799048 5.16959260
64 1.37449458 7.37799048
65 2.91081636 1.37449458
66 -2.60483694 2.91081636
67 -5.00062393 -2.60483694
68 -2.89162658 -5.00062393
69 -0.20863547 -2.89162658
70 -6.13241269 -0.20863547
71 -9.32408475 -6.13241269
72 4.64496878 -9.32408475
73 3.84849729 4.64496878
74 7.89184702 3.84849729
75 5.84305160 7.89184702
76 -0.44067785 5.84305160
77 1.70178108 -0.44067785
78 -1.66798211 1.70178108
79 -4.54881028 -1.66798211
80 -1.81682662 -4.54881028
81 0.60678379 -1.81682662
82 -3.49009232 0.60678379
83 -6.90469966 -3.49009232
> 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/73nng1292676443.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/rcomp/tmp/83nng1292676443.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/rcomp/tmp/93nng1292676443.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/rcomp/tmp/10vemj1292676443.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/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/11zfk61292676443.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/12aok91292676443.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/13g7zl1292676443.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/14ryy61292676443.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/15uzwu1292676443.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/168rc31292676443.tab")
+ }
>
> try(system("convert tmp/17d7p1292676443.ps tmp/17d7p1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/27d7p1292676443.ps tmp/27d7p1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/3hnoa1292676443.ps tmp/3hnoa1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/4hnoa1292676443.ps tmp/4hnoa1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/5hnoa1292676443.ps tmp/5hnoa1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/6aw6d1292676443.ps tmp/6aw6d1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/73nng1292676443.ps tmp/73nng1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/83nng1292676443.ps tmp/83nng1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/93nng1292676443.ps tmp/93nng1292676443.png",intern=TRUE))
character(0)
> try(system("convert tmp/10vemj1292676443.ps tmp/10vemj1292676443.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.808 1.647 14.138