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(31514
+ ,-9
+ ,0
+ ,8.3
+ ,1.2
+ ,27071
+ ,-13
+ ,4
+ ,8.2
+ ,1.7
+ ,29462
+ ,-18
+ ,5
+ ,8
+ ,1.8
+ ,26105
+ ,-11
+ ,-7
+ ,7.9
+ ,1.5
+ ,22397
+ ,-9
+ ,-2
+ ,7.6
+ ,1
+ ,23843
+ ,-10
+ ,1
+ ,7.6
+ ,1.6
+ ,21705
+ ,-13
+ ,3
+ ,8.3
+ ,1.5
+ ,18089
+ ,-11
+ ,-2
+ ,8.4
+ ,1.8
+ ,20764
+ ,-5
+ ,-6
+ ,8.4
+ ,1.8
+ ,25316
+ ,-15
+ ,10
+ ,8.4
+ ,1.6
+ ,17704
+ ,-6
+ ,-9
+ ,8.4
+ ,1.9
+ ,15548
+ ,-6
+ ,0
+ ,8.6
+ ,1.7
+ ,28029
+ ,-3
+ ,-3
+ ,8.9
+ ,1.6
+ ,29383
+ ,-1
+ ,-2
+ ,8.8
+ ,1.3
+ ,36438
+ ,-3
+ ,2
+ ,8.3
+ ,1.1
+ ,32034
+ ,-4
+ ,1
+ ,7.5
+ ,1.9
+ ,22679
+ ,-6
+ ,2
+ ,7.2
+ ,2.6
+ ,24319
+ ,0
+ ,-6
+ ,7.4
+ ,2.3
+ ,18004
+ ,-4
+ ,4
+ ,8.8
+ ,2.4
+ ,17537
+ ,-2
+ ,-2
+ ,9.3
+ ,2.2
+ ,20366
+ ,-2
+ ,0
+ ,9.3
+ ,2
+ ,22782
+ ,-6
+ ,4
+ ,8.7
+ ,2.9
+ ,19169
+ ,-7
+ ,1
+ ,8.2
+ ,2.6
+ ,13807
+ ,-6
+ ,-1
+ ,8.3
+ ,2.3
+ ,29743
+ ,-6
+ ,0
+ ,8.5
+ ,2.3
+ ,25591
+ ,-3
+ ,-3
+ ,8.6
+ ,2.6
+ ,29096
+ ,-2
+ ,-1
+ ,8.5
+ ,3.1
+ ,26482
+ ,-5
+ ,3
+ ,8.2
+ ,2.8
+ ,22405
+ ,-11
+ ,6
+ ,8.1
+ ,2.5
+ ,27044
+ ,-11
+ ,0
+ ,7.9
+ ,2.9
+ ,17970
+ ,-11
+ ,0
+ ,8.6
+ ,3.1
+ ,18730
+ ,-10
+ ,-1
+ ,8.7
+ ,3.1
+ ,19684
+ ,-14
+ ,4
+ ,8.7
+ ,3.2
+ ,19785
+ ,-8
+ ,-6
+ ,8.5
+ ,2.5
+ ,18479
+ ,-9
+ ,1
+ ,8.4
+ ,2.6
+ ,10698
+ ,-5
+ ,-4
+ ,8.5
+ ,2.9
+ ,31956
+ ,-1
+ ,-4
+ ,8.7
+ ,2.6
+ ,29506
+ ,-2
+ ,1
+ ,8.7
+ ,2.4
+ ,34506
+ ,-5
+ ,3
+ ,8.6
+ ,1.7
+ ,27165
+ ,-4
+ ,-1
+ ,8.5
+ ,2
+ ,26736
+ ,-6
+ ,2
+ ,8.3
+ ,2.2
+ ,23691
+ ,-2
+ ,-4
+ ,8
+ ,1.9
+ ,18157
+ ,-2
+ ,0
+ ,8.2
+ ,1.6
+ ,17328
+ ,-2
+ ,0
+ ,8.1
+ ,1.6
+ ,18205
+ ,-2
+ ,0
+ ,8.1
+ ,1.2
+ ,20995
+ ,2
+ ,-4
+ ,8
+ ,1.2
+ ,17382
+ ,1
+ ,1
+ ,7.9
+ ,1.5
+ ,9367
+ ,-8
+ ,9
+ ,7.9
+ ,1.6
+ ,31124
+ ,-1
+ ,-7
+ ,8
+ ,1.7
+ ,26551
+ ,1
+ ,-2
+ ,8
+ ,1.8
+ ,30651
+ ,-1
+ ,2
+ ,7.9
+ ,1.8
+ ,25859
+ ,2
+ ,-3
+ ,8
+ ,1.8
+ ,25100
+ ,2
+ ,0
+ ,7.7
+ ,1.3
+ ,25778
+ ,1
+ ,1
+ ,7.2
+ ,1.3
+ ,20418
+ ,-1
+ ,2
+ ,7.5
+ ,1.4
+ ,18688
+ ,-2
+ ,1
+ ,7.3
+ ,1.1
+ ,20424
+ ,-2
+ ,0
+ ,7
+ ,1.5
+ ,24776
+ ,-1
+ ,-1
+ ,7
+ ,2.2
+ ,19814
+ ,-8
+ ,7
+ ,7
+ ,2.9
+ ,12738
+ ,-4
+ ,-4
+ ,7.2
+ ,3.1
+ ,31566
+ ,-6
+ ,2
+ ,7.3
+ ,3.5
+ ,30111
+ ,-3
+ ,-3
+ ,7.1
+ ,3.6
+ ,30019
+ ,-3
+ ,0
+ ,6.8
+ ,4.4
+ ,31934
+ ,-7
+ ,4
+ ,6.4
+ ,4.2
+ ,25826
+ ,-9
+ ,2
+ ,6.1
+ ,5.2
+ ,26835
+ ,-11
+ ,2
+ ,6.5
+ ,5.8
+ ,20205
+ ,-13
+ ,2
+ ,7.7
+ ,5.9
+ ,17789
+ ,-11
+ ,-2
+ ,7.9
+ ,5.4
+ ,20520
+ ,-9
+ ,-2
+ ,7.5
+ ,5.5
+ ,22518
+ ,-17
+ ,8
+ ,6.9
+ ,4.7
+ ,15572
+ ,-22
+ ,5
+ ,6.6
+ ,3.1
+ ,11509
+ ,-25
+ ,3
+ ,6.9
+ ,2.6
+ ,25447
+ ,-20
+ ,-5
+ ,7.7
+ ,2.3
+ ,24090
+ ,-24
+ ,4
+ ,8
+ ,1.9
+ ,27786
+ ,-24
+ ,0
+ ,8
+ ,0.6
+ ,26195
+ ,-22
+ ,-2
+ ,7.7
+ ,0.6
+ ,20516
+ ,-19
+ ,-3
+ ,7.3
+ ,-0.4
+ ,22759
+ ,-18
+ ,-1
+ ,7.4
+ ,-1.1
+ ,19028
+ ,-17
+ ,-1
+ ,8.1
+ ,-1.7
+ ,16971
+ ,-11
+ ,-6
+ ,8.3
+ ,-0.8
+ ,20036
+ ,-11
+ ,0
+ ,8.1
+ ,-1.2
+ ,22485
+ ,-12
+ ,1
+ ,7.9
+ ,-1
+ ,18730
+ ,-10
+ ,-2
+ ,7.9
+ ,-0.1
+ ,14538
+ ,-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 31514 -9 0
2 27071 -13 4
3 29462 -18 5
4 26105 -11 -7
5 22397 -9 -2
6 23843 -10 1
7 21705 -13 3
8 18089 -11 -2
9 20764 -5 -6
10 25316 -15 10
11 17704 -6 -9
12 15548 -6 0
13 28029 -3 -3
14 29383 -1 -2
15 36438 -3 2
16 32034 -4 1
17 22679 -6 2
18 24319 0 -6
19 18004 -4 4
20 17537 -2 -2
21 20366 -2 0
22 22782 -6 4
23 19169 -7 1
24 13807 -6 -1
25 29743 -6 0
26 25591 -3 -3
27 29096 -2 -1
28 26482 -5 3
29 22405 -11 6
30 27044 -11 0
31 17970 -11 0
32 18730 -10 -1
33 19684 -14 4
34 19785 -8 -6
35 18479 -9 1
36 10698 -5 -4
37 31956 -1 -4
38 29506 -2 1
39 34506 -5 3
40 27165 -4 -1
41 26736 -6 2
42 23691 -2 -4
43 18157 -2 0
44 17328 -2 0
45 18205 -2 0
46 20995 2 -4
47 17382 1 1
48 9367 -8 9
49 31124 -1 -7
50 26551 1 -2
51 30651 -1 2
52 25859 2 -3
53 25100 2 0
54 25778 1 1
55 20418 -1 2
56 18688 -2 1
57 20424 -2 0
58 24776 -1 -1
59 19814 -8 7
60 12738 -4 -4
61 31566 -6 2
62 30111 -3 -3
63 30019 -3 0
64 31934 -7 4
65 25826 -9 2
66 26835 -11 2
67 20205 -13 2
68 17789 -11 -2
69 20520 -9 -2
70 22518 -17 8
71 15572 -22 5
72 11509 -25 3
73 25447 -20 -5
74 24090 -24 4
75 27786 -24 0
76 26195 -22 -2
77 20516 -19 -3
78 22759 -18 -1
79 19028 -17 -1
80 16971 -11 -6
81 20036 -11 0
82 22485 -12 1
83 18730 -10 -2
84 14538 -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
29192.30 174.39
Evolutie_consumentenvertrouwen Totaal_Werkloosheid
42.88 -649.24
Algemene_index
135.26
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-13903.5 -4079.4 -184.3 4541.3 12923.1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 29192.30 8666.02 3.369 0.00117 **
Consumentenvertrouwen 174.39 101.42 1.720 0.08944 .
Evolutie_consumentenvertrouwen 42.88 184.36 0.233 0.81668
Totaal_Werkloosheid -649.24 1024.48 -0.634 0.52809
Algemene_index 135.26 458.34 0.295 0.76869
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5721 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/1zetz1292690628.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/2zetz1292690628.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/3r5tl1292690628.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/4r5tl1292690628.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/5kwso1292690628.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
9117.62325 5068.11365 8144.80781 4057.27605 -341.03941 1069.55965
7 8 9 10 11 12
-163.03174 -3889.07144 -2088.89446 3547.99661 -4859.39544 -7244.40098
13 14 15 16 17 18
5050.36521 5988.36084 12923.05530 8108.72082 -1229.83295 -122.71568
19 20 21 22 23 24
-5273.52571 -5480.35990 -2710.06514 -279.30041 -3873.32036 -9218.45094
25 26 27 28 29 30
6804.51953 2282.33343 5394.63338 2978.09209 -205.55339 4506.76509
31 32 33 34 35 36
-4139.81533 -3446.40188 -2022.76211 -2574.48243 -4084.69248 -12324.51113
37 38 39 40 41 42
8406.35741 5943.40625 11410.57416 3961.19655 3595.43952 -44.04351
43 44 45 46 47 48
-5579.13091 -6473.05537 -5541.95202 -3342.92055 -7101.42535 -13903.47359
49 50 51 52 53 54
7370.25407 2220.55698 6432.89763 1397.04596 382.26638 867.15507
55 56 57 58 59 60
-4005.69688 -5607.70036 -4077.69865 48.10947 -4130.87274 -11329.97033
61 62 63 64 65 66
7600.35897 5693.20808 5169.59260 7377.99048 1374.49458 2910.81636
67 68 69 70 71 72
-2604.83694 -5000.62393 -2891.62658 -208.63547 -6132.41269 -9324.08475
73 74 75 76 77 78
4644.96878 3848.49729 7891.84702 5843.05160 -440.67785 1701.78108
79 80 81 82 83 84
-1667.98211 -4548.81028 -1816.82662 606.78379 -3490.09232 -6904.69966
> postscript(file="/var/www/html/rcomp/tmp/6kwso1292690628.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 9117.62325 NA
1 5068.11365 9117.62325
2 8144.80781 5068.11365
3 4057.27605 8144.80781
4 -341.03941 4057.27605
5 1069.55965 -341.03941
6 -163.03174 1069.55965
7 -3889.07144 -163.03174
8 -2088.89446 -3889.07144
9 3547.99661 -2088.89446
10 -4859.39544 3547.99661
11 -7244.40098 -4859.39544
12 5050.36521 -7244.40098
13 5988.36084 5050.36521
14 12923.05530 5988.36084
15 8108.72082 12923.05530
16 -1229.83295 8108.72082
17 -122.71568 -1229.83295
18 -5273.52571 -122.71568
19 -5480.35990 -5273.52571
20 -2710.06514 -5480.35990
21 -279.30041 -2710.06514
22 -3873.32036 -279.30041
23 -9218.45094 -3873.32036
24 6804.51953 -9218.45094
25 2282.33343 6804.51953
26 5394.63338 2282.33343
27 2978.09209 5394.63338
28 -205.55339 2978.09209
29 4506.76509 -205.55339
30 -4139.81533 4506.76509
31 -3446.40188 -4139.81533
32 -2022.76211 -3446.40188
33 -2574.48243 -2022.76211
34 -4084.69248 -2574.48243
35 -12324.51113 -4084.69248
36 8406.35741 -12324.51113
37 5943.40625 8406.35741
38 11410.57416 5943.40625
39 3961.19655 11410.57416
40 3595.43952 3961.19655
41 -44.04351 3595.43952
42 -5579.13091 -44.04351
43 -6473.05537 -5579.13091
44 -5541.95202 -6473.05537
45 -3342.92055 -5541.95202
46 -7101.42535 -3342.92055
47 -13903.47359 -7101.42535
48 7370.25407 -13903.47359
49 2220.55698 7370.25407
50 6432.89763 2220.55698
51 1397.04596 6432.89763
52 382.26638 1397.04596
53 867.15507 382.26638
54 -4005.69688 867.15507
55 -5607.70036 -4005.69688
56 -4077.69865 -5607.70036
57 48.10947 -4077.69865
58 -4130.87274 48.10947
59 -11329.97033 -4130.87274
60 7600.35897 -11329.97033
61 5693.20808 7600.35897
62 5169.59260 5693.20808
63 7377.99048 5169.59260
64 1374.49458 7377.99048
65 2910.81636 1374.49458
66 -2604.83694 2910.81636
67 -5000.62393 -2604.83694
68 -2891.62658 -5000.62393
69 -208.63547 -2891.62658
70 -6132.41269 -208.63547
71 -9324.08475 -6132.41269
72 4644.96878 -9324.08475
73 3848.49729 4644.96878
74 7891.84702 3848.49729
75 5843.05160 7891.84702
76 -440.67785 5843.05160
77 1701.78108 -440.67785
78 -1667.98211 1701.78108
79 -4548.81028 -1667.98211
80 -1816.82662 -4548.81028
81 606.78379 -1816.82662
82 -3490.09232 606.78379
83 -6904.69966 -3490.09232
84 NA -6904.69966
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 5068.11365 9117.62325
[2,] 8144.80781 5068.11365
[3,] 4057.27605 8144.80781
[4,] -341.03941 4057.27605
[5,] 1069.55965 -341.03941
[6,] -163.03174 1069.55965
[7,] -3889.07144 -163.03174
[8,] -2088.89446 -3889.07144
[9,] 3547.99661 -2088.89446
[10,] -4859.39544 3547.99661
[11,] -7244.40098 -4859.39544
[12,] 5050.36521 -7244.40098
[13,] 5988.36084 5050.36521
[14,] 12923.05530 5988.36084
[15,] 8108.72082 12923.05530
[16,] -1229.83295 8108.72082
[17,] -122.71568 -1229.83295
[18,] -5273.52571 -122.71568
[19,] -5480.35990 -5273.52571
[20,] -2710.06514 -5480.35990
[21,] -279.30041 -2710.06514
[22,] -3873.32036 -279.30041
[23,] -9218.45094 -3873.32036
[24,] 6804.51953 -9218.45094
[25,] 2282.33343 6804.51953
[26,] 5394.63338 2282.33343
[27,] 2978.09209 5394.63338
[28,] -205.55339 2978.09209
[29,] 4506.76509 -205.55339
[30,] -4139.81533 4506.76509
[31,] -3446.40188 -4139.81533
[32,] -2022.76211 -3446.40188
[33,] -2574.48243 -2022.76211
[34,] -4084.69248 -2574.48243
[35,] -12324.51113 -4084.69248
[36,] 8406.35741 -12324.51113
[37,] 5943.40625 8406.35741
[38,] 11410.57416 5943.40625
[39,] 3961.19655 11410.57416
[40,] 3595.43952 3961.19655
[41,] -44.04351 3595.43952
[42,] -5579.13091 -44.04351
[43,] -6473.05537 -5579.13091
[44,] -5541.95202 -6473.05537
[45,] -3342.92055 -5541.95202
[46,] -7101.42535 -3342.92055
[47,] -13903.47359 -7101.42535
[48,] 7370.25407 -13903.47359
[49,] 2220.55698 7370.25407
[50,] 6432.89763 2220.55698
[51,] 1397.04596 6432.89763
[52,] 382.26638 1397.04596
[53,] 867.15507 382.26638
[54,] -4005.69688 867.15507
[55,] -5607.70036 -4005.69688
[56,] -4077.69865 -5607.70036
[57,] 48.10947 -4077.69865
[58,] -4130.87274 48.10947
[59,] -11329.97033 -4130.87274
[60,] 7600.35897 -11329.97033
[61,] 5693.20808 7600.35897
[62,] 5169.59260 5693.20808
[63,] 7377.99048 5169.59260
[64,] 1374.49458 7377.99048
[65,] 2910.81636 1374.49458
[66,] -2604.83694 2910.81636
[67,] -5000.62393 -2604.83694
[68,] -2891.62658 -5000.62393
[69,] -208.63547 -2891.62658
[70,] -6132.41269 -208.63547
[71,] -9324.08475 -6132.41269
[72,] 4644.96878 -9324.08475
[73,] 3848.49729 4644.96878
[74,] 7891.84702 3848.49729
[75,] 5843.05160 7891.84702
[76,] -440.67785 5843.05160
[77,] 1701.78108 -440.67785
[78,] -1667.98211 1701.78108
[79,] -4548.81028 -1667.98211
[80,] -1816.82662 -4548.81028
[81,] 606.78379 -1816.82662
[82,] -3490.09232 606.78379
[83,] -6904.69966 -3490.09232
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 5068.11365 9117.62325
2 8144.80781 5068.11365
3 4057.27605 8144.80781
4 -341.03941 4057.27605
5 1069.55965 -341.03941
6 -163.03174 1069.55965
7 -3889.07144 -163.03174
8 -2088.89446 -3889.07144
9 3547.99661 -2088.89446
10 -4859.39544 3547.99661
11 -7244.40098 -4859.39544
12 5050.36521 -7244.40098
13 5988.36084 5050.36521
14 12923.05530 5988.36084
15 8108.72082 12923.05530
16 -1229.83295 8108.72082
17 -122.71568 -1229.83295
18 -5273.52571 -122.71568
19 -5480.35990 -5273.52571
20 -2710.06514 -5480.35990
21 -279.30041 -2710.06514
22 -3873.32036 -279.30041
23 -9218.45094 -3873.32036
24 6804.51953 -9218.45094
25 2282.33343 6804.51953
26 5394.63338 2282.33343
27 2978.09209 5394.63338
28 -205.55339 2978.09209
29 4506.76509 -205.55339
30 -4139.81533 4506.76509
31 -3446.40188 -4139.81533
32 -2022.76211 -3446.40188
33 -2574.48243 -2022.76211
34 -4084.69248 -2574.48243
35 -12324.51113 -4084.69248
36 8406.35741 -12324.51113
37 5943.40625 8406.35741
38 11410.57416 5943.40625
39 3961.19655 11410.57416
40 3595.43952 3961.19655
41 -44.04351 3595.43952
42 -5579.13091 -44.04351
43 -6473.05537 -5579.13091
44 -5541.95202 -6473.05537
45 -3342.92055 -5541.95202
46 -7101.42535 -3342.92055
47 -13903.47359 -7101.42535
48 7370.25407 -13903.47359
49 2220.55698 7370.25407
50 6432.89763 2220.55698
51 1397.04596 6432.89763
52 382.26638 1397.04596
53 867.15507 382.26638
54 -4005.69688 867.15507
55 -5607.70036 -4005.69688
56 -4077.69865 -5607.70036
57 48.10947 -4077.69865
58 -4130.87274 48.10947
59 -11329.97033 -4130.87274
60 7600.35897 -11329.97033
61 5693.20808 7600.35897
62 5169.59260 5693.20808
63 7377.99048 5169.59260
64 1374.49458 7377.99048
65 2910.81636 1374.49458
66 -2604.83694 2910.81636
67 -5000.62393 -2604.83694
68 -2891.62658 -5000.62393
69 -208.63547 -2891.62658
70 -6132.41269 -208.63547
71 -9324.08475 -6132.41269
72 4644.96878 -9324.08475
73 3848.49729 4644.96878
74 7891.84702 3848.49729
75 5843.05160 7891.84702
76 -440.67785 5843.05160
77 1701.78108 -440.67785
78 -1667.98211 1701.78108
79 -4548.81028 -1667.98211
80 -1816.82662 -4548.81028
81 606.78379 -1816.82662
82 -3490.09232 606.78379
83 -6904.69966 -3490.09232
> 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/7d6r81292690628.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/8d6r81292690628.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/9d6r81292690628.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/105x8t1292690628.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/119f7z1292690628.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/12ugnn1292690628.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/13q8le1292690628.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/14c8jj1292690628.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/15x9i81292690628.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/161rzd1292690628.tab")
+ }
>
> try(system("convert tmp/1zetz1292690628.ps tmp/1zetz1292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/2zetz1292690628.ps tmp/2zetz1292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/3r5tl1292690628.ps tmp/3r5tl1292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/4r5tl1292690628.ps tmp/4r5tl1292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/5kwso1292690628.ps tmp/5kwso1292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/6kwso1292690628.ps tmp/6kwso1292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/7d6r81292690628.ps tmp/7d6r81292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/8d6r81292690628.ps tmp/8d6r81292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/9d6r81292690628.ps tmp/9d6r81292690628.png",intern=TRUE))
character(0)
> try(system("convert tmp/105x8t1292690628.ps tmp/105x8t1292690628.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.858 1.705 6.717