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(4.8,19.2,5.5,26.6,5.4,26.6,5.9,31.4,5.8,31.2,5.1,26.4,4.1,20.7,4.4,20.7,3.6,15,3.5,13.3,3.1,8.7,2.9,10.2,2.2,4.3,1.4,-0.1,1.2,-4.6,1.3,-3.9,1.3,-3.5,1.3,-3.4,1.8,-2.5,1.8,-1.1,1.8,0.3,1.7,-0.9,2.1,3.6,2,2.7,1.7,-0.2,1.9,-1,2.3,5.8,2.4,6.4,2.5,9.6,2.8,13.2,2.6,10.6,2.2,10.9,2.8,12.9,2.8,15.9,2.8,12.2,2.3,9.1,2.2,9,3,17.4,2.9,14.7,2.7,17,2.7,13.7,2.3,9.5,2.4,14.8,2.8,13.6,2.3,12.6,2,8.9,1.9,10.2,2.3,12.7,2.7,16,1.8,10.4,2,9.9,2.1,9.5,2,8.6,2.4,10,1.7,3.5,1,-4.2,1.2,-4.4,1.4,-1.5,1.7,-0.1,1.8,0.8),dim=c(2,60),dimnames=list(c('Inflatie_België','Inflatie_energiedragers'),1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('Inflatie_België','Inflatie_energiedragers'),1:60))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = '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
Inflatie_Belgi\353 Inflatie_energiedragers M1 M2 M3 M4 M5 M6 M7 M8 M9 M10
1 4.8 19.2 1 0 0 0 0 0 0 0 0 0
2 5.5 26.6 0 1 0 0 0 0 0 0 0 0
3 5.4 26.6 0 0 1 0 0 0 0 0 0 0
4 5.9 31.4 0 0 0 1 0 0 0 0 0 0
5 5.8 31.2 0 0 0 0 1 0 0 0 0 0
6 5.1 26.4 0 0 0 0 0 1 0 0 0 0
7 4.1 20.7 0 0 0 0 0 0 1 0 0 0
8 4.4 20.7 0 0 0 0 0 0 0 1 0 0
9 3.6 15.0 0 0 0 0 0 0 0 0 1 0
10 3.5 13.3 0 0 0 0 0 0 0 0 0 1
11 3.1 8.7 0 0 0 0 0 0 0 0 0 0
12 2.9 10.2 0 0 0 0 0 0 0 0 0 0
13 2.2 4.3 1 0 0 0 0 0 0 0 0 0
14 1.4 -0.1 0 1 0 0 0 0 0 0 0 0
15 1.2 -4.6 0 0 1 0 0 0 0 0 0 0
16 1.3 -3.9 0 0 0 1 0 0 0 0 0 0
17 1.3 -3.5 0 0 0 0 1 0 0 0 0 0
18 1.3 -3.4 0 0 0 0 0 1 0 0 0 0
19 1.8 -2.5 0 0 0 0 0 0 1 0 0 0
20 1.8 -1.1 0 0 0 0 0 0 0 1 0 0
21 1.8 0.3 0 0 0 0 0 0 0 0 1 0
22 1.7 -0.9 0 0 0 0 0 0 0 0 0 1
23 2.1 3.6 0 0 0 0 0 0 0 0 0 0
24 2.0 2.7 0 0 0 0 0 0 0 0 0 0
25 1.7 -0.2 1 0 0 0 0 0 0 0 0 0
26 1.9 -1.0 0 1 0 0 0 0 0 0 0 0
27 2.3 5.8 0 0 1 0 0 0 0 0 0 0
28 2.4 6.4 0 0 0 1 0 0 0 0 0 0
29 2.5 9.6 0 0 0 0 1 0 0 0 0 0
30 2.8 13.2 0 0 0 0 0 1 0 0 0 0
31 2.6 10.6 0 0 0 0 0 0 1 0 0 0
32 2.2 10.9 0 0 0 0 0 0 0 1 0 0
33 2.8 12.9 0 0 0 0 0 0 0 0 1 0
34 2.8 15.9 0 0 0 0 0 0 0 0 0 1
35 2.8 12.2 0 0 0 0 0 0 0 0 0 0
36 2.3 9.1 0 0 0 0 0 0 0 0 0 0
37 2.2 9.0 1 0 0 0 0 0 0 0 0 0
38 3.0 17.4 0 1 0 0 0 0 0 0 0 0
39 2.9 14.7 0 0 1 0 0 0 0 0 0 0
40 2.7 17.0 0 0 0 1 0 0 0 0 0 0
41 2.7 13.7 0 0 0 0 1 0 0 0 0 0
42 2.3 9.5 0 0 0 0 0 1 0 0 0 0
43 2.4 14.8 0 0 0 0 0 0 1 0 0 0
44 2.8 13.6 0 0 0 0 0 0 0 1 0 0
45 2.3 12.6 0 0 0 0 0 0 0 0 1 0
46 2.0 8.9 0 0 0 0 0 0 0 0 0 1
47 1.9 10.2 0 0 0 0 0 0 0 0 0 0
48 2.3 12.7 0 0 0 0 0 0 0 0 0 0
49 2.7 16.0 1 0 0 0 0 0 0 0 0 0
50 1.8 10.4 0 1 0 0 0 0 0 0 0 0
51 2.0 9.9 0 0 1 0 0 0 0 0 0 0
52 2.1 9.5 0 0 0 1 0 0 0 0 0 0
53 2.0 8.6 0 0 0 0 1 0 0 0 0 0
54 2.4 10.0 0 0 0 0 0 1 0 0 0 0
55 1.7 3.5 0 0 0 0 0 0 1 0 0 0
56 1.0 -4.2 0 0 0 0 0 0 0 1 0 0
57 1.2 -4.4 0 0 0 0 0 0 0 0 1 0
58 1.4 -1.5 0 0 0 0 0 0 0 0 0 1
59 1.7 -0.1 0 0 0 0 0 0 0 0 0 0
60 1.8 0.8 0 0 0 0 0 0 0 0 0 0
M11
1 0
2 0
3 0
4 0
5 0
6 0
7 0
8 0
9 0
10 0
11 1
12 0
13 0
14 0
15 0
16 0
17 0
18 0
19 0
20 0
21 0
22 0
23 1
24 0
25 0
26 0
27 0
28 0
29 0
30 0
31 0
32 0
33 0
34 0
35 1
36 0
37 0
38 0
39 0
40 0
41 0
42 0
43 0
44 0
45 0
46 0
47 1
48 0
49 0
50 0
51 0
52 0
53 0
54 0
55 0
56 0
57 0
58 0
59 1
60 0
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Inflatie_energiedragers M1
1.42012 0.11829 0.15717
M2 M3 M4
0.03887 0.10017 0.03090
M5 M6 M7
0.02983 0.04209 -0.01444
M8 M9 M10
0.07590 0.05871 0.01527
M11
0.08129
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.8892 -0.3884 0.1038 0.2918 0.9515
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.420116 0.244267 5.814 5.14e-07 ***
Inflatie_energiedragers 0.118294 0.007805 15.155 < 2e-16 ***
M1 0.157169 0.337030 0.466 0.643
M2 0.038875 0.337583 0.115 0.909
M3 0.100168 0.337470 0.297 0.768
M4 0.030898 0.338675 0.091 0.928
M5 0.029825 0.338534 0.088 0.930
M6 0.042094 0.337912 0.125 0.901
M7 -0.014441 0.336924 -0.043 0.966
M8 0.075902 0.336508 0.226 0.823
M9 0.058707 0.336440 0.174 0.862
M10 0.015268 0.336438 0.045 0.964
M11 0.081293 0.336440 0.242 0.810
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.532 on 47 degrees of freedom
Multiple R-squared: 0.8366, Adjusted R-squared: 0.7949
F-statistic: 20.05 on 12 and 47 DF, p-value: 1.497e-14
> 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.6986167 0.60276661 0.30138330
[2,] 0.6213910 0.75721796 0.37860898
[3,] 0.5420208 0.91595838 0.45797919
[4,] 0.8020345 0.39593105 0.19796553
[5,] 0.7570306 0.48593883 0.24296942
[6,] 0.6684735 0.66305301 0.33152651
[7,] 0.5749275 0.85014509 0.42507255
[8,] 0.5340147 0.93197053 0.46598526
[9,] 0.4419728 0.88394554 0.55802723
[10,] 0.3779623 0.75592467 0.62203767
[11,] 0.4468324 0.89366487 0.55316757
[12,] 0.4486810 0.89736193 0.55131904
[13,] 0.4980286 0.99605722 0.50197139
[14,] 0.5537579 0.89248410 0.44624205
[15,] 0.6296413 0.74071736 0.37035868
[16,] 0.7158238 0.56835238 0.28417619
[17,] 0.8663654 0.26726926 0.13363463
[18,] 0.8797872 0.24042551 0.12021276
[19,] 0.9199939 0.16001221 0.08000611
[20,] 0.9354009 0.12919829 0.06459914
[21,] 0.9069243 0.18615147 0.09307574
[22,] 0.8975543 0.20489132 0.10244566
[23,] 0.9612045 0.07759090 0.03879545
[24,] 0.9762079 0.04758411 0.02379206
[25,] 0.9691804 0.06163929 0.03081964
[26,] 0.9588149 0.08237019 0.04118509
[27,] 0.9136305 0.17273894 0.08636947
[28,] 0.8525457 0.29490864 0.14745432
[29,] 0.9519019 0.09619627 0.04809814
> postscript(file="/var/www/html/rcomp/tmp/1qqxy1227717843.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/226171227717843.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/3u3l91227717843.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/46pbw1227717843.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/50rfk1227717843.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 = 60
Frequency = 1
1 2 3 4 5
0.9514795091 0.8944007730 0.7331079335 0.7345685655 0.6593003076
6 7 8 9 10
0.5148403887 0.2456487277 0.4553060121 0.3467737746 0.4913117165
11 12 13 14 15
0.5694374765 0.2732899872 0.1140534415 -0.0471613750 0.2238667718
16 17 18 19 20
0.3103309691 0.2640865797 0.2399882535 0.6900591459 0.4341054568
21 22 23 24 25
0.2856889965 0.3710801622 0.1727345943 0.2604916310 0.1463744277
26 27 28 29 30
0.5593028223 0.0936138257 0.1919073783 -0.0855589582 -0.2236847182
31 32 33 34 35
-0.0595863920 -0.5854171733 -0.2048097651 -0.5162515200 -0.1445899573
36 37 38 39 40
-0.1965871050 -0.4419262553 -0.5172985439 -0.3591987916 -0.7620042783
41 42 43 44 45
-0.3705625235 -0.2859985739 -0.7564193125 -0.3048097651 -0.6693216994
46 47 48 49 50
-0.4881966524 -0.8080028522 -0.6224438941 -0.7699811229 -0.8892436763
51 52 53 54 55
-0.6913897395 -0.4748026345 -0.4672654057 -0.2451453501 -0.1197021691
56 57 58 59 60
0.0008154696 0.2416686933 0.1420562937 0.2104207386 0.2852493808
> postscript(file="/var/www/html/rcomp/tmp/6ue201227717843.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 = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 0.9514795091 NA
1 0.8944007730 0.9514795091
2 0.7331079335 0.8944007730
3 0.7345685655 0.7331079335
4 0.6593003076 0.7345685655
5 0.5148403887 0.6593003076
6 0.2456487277 0.5148403887
7 0.4553060121 0.2456487277
8 0.3467737746 0.4553060121
9 0.4913117165 0.3467737746
10 0.5694374765 0.4913117165
11 0.2732899872 0.5694374765
12 0.1140534415 0.2732899872
13 -0.0471613750 0.1140534415
14 0.2238667718 -0.0471613750
15 0.3103309691 0.2238667718
16 0.2640865797 0.3103309691
17 0.2399882535 0.2640865797
18 0.6900591459 0.2399882535
19 0.4341054568 0.6900591459
20 0.2856889965 0.4341054568
21 0.3710801622 0.2856889965
22 0.1727345943 0.3710801622
23 0.2604916310 0.1727345943
24 0.1463744277 0.2604916310
25 0.5593028223 0.1463744277
26 0.0936138257 0.5593028223
27 0.1919073783 0.0936138257
28 -0.0855589582 0.1919073783
29 -0.2236847182 -0.0855589582
30 -0.0595863920 -0.2236847182
31 -0.5854171733 -0.0595863920
32 -0.2048097651 -0.5854171733
33 -0.5162515200 -0.2048097651
34 -0.1445899573 -0.5162515200
35 -0.1965871050 -0.1445899573
36 -0.4419262553 -0.1965871050
37 -0.5172985439 -0.4419262553
38 -0.3591987916 -0.5172985439
39 -0.7620042783 -0.3591987916
40 -0.3705625235 -0.7620042783
41 -0.2859985739 -0.3705625235
42 -0.7564193125 -0.2859985739
43 -0.3048097651 -0.7564193125
44 -0.6693216994 -0.3048097651
45 -0.4881966524 -0.6693216994
46 -0.8080028522 -0.4881966524
47 -0.6224438941 -0.8080028522
48 -0.7699811229 -0.6224438941
49 -0.8892436763 -0.7699811229
50 -0.6913897395 -0.8892436763
51 -0.4748026345 -0.6913897395
52 -0.4672654057 -0.4748026345
53 -0.2451453501 -0.4672654057
54 -0.1197021691 -0.2451453501
55 0.0008154696 -0.1197021691
56 0.2416686933 0.0008154696
57 0.1420562937 0.2416686933
58 0.2104207386 0.1420562937
59 0.2852493808 0.2104207386
60 NA 0.2852493808
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 0.8944007730 0.9514795091
[2,] 0.7331079335 0.8944007730
[3,] 0.7345685655 0.7331079335
[4,] 0.6593003076 0.7345685655
[5,] 0.5148403887 0.6593003076
[6,] 0.2456487277 0.5148403887
[7,] 0.4553060121 0.2456487277
[8,] 0.3467737746 0.4553060121
[9,] 0.4913117165 0.3467737746
[10,] 0.5694374765 0.4913117165
[11,] 0.2732899872 0.5694374765
[12,] 0.1140534415 0.2732899872
[13,] -0.0471613750 0.1140534415
[14,] 0.2238667718 -0.0471613750
[15,] 0.3103309691 0.2238667718
[16,] 0.2640865797 0.3103309691
[17,] 0.2399882535 0.2640865797
[18,] 0.6900591459 0.2399882535
[19,] 0.4341054568 0.6900591459
[20,] 0.2856889965 0.4341054568
[21,] 0.3710801622 0.2856889965
[22,] 0.1727345943 0.3710801622
[23,] 0.2604916310 0.1727345943
[24,] 0.1463744277 0.2604916310
[25,] 0.5593028223 0.1463744277
[26,] 0.0936138257 0.5593028223
[27,] 0.1919073783 0.0936138257
[28,] -0.0855589582 0.1919073783
[29,] -0.2236847182 -0.0855589582
[30,] -0.0595863920 -0.2236847182
[31,] -0.5854171733 -0.0595863920
[32,] -0.2048097651 -0.5854171733
[33,] -0.5162515200 -0.2048097651
[34,] -0.1445899573 -0.5162515200
[35,] -0.1965871050 -0.1445899573
[36,] -0.4419262553 -0.1965871050
[37,] -0.5172985439 -0.4419262553
[38,] -0.3591987916 -0.5172985439
[39,] -0.7620042783 -0.3591987916
[40,] -0.3705625235 -0.7620042783
[41,] -0.2859985739 -0.3705625235
[42,] -0.7564193125 -0.2859985739
[43,] -0.3048097651 -0.7564193125
[44,] -0.6693216994 -0.3048097651
[45,] -0.4881966524 -0.6693216994
[46,] -0.8080028522 -0.4881966524
[47,] -0.6224438941 -0.8080028522
[48,] -0.7699811229 -0.6224438941
[49,] -0.8892436763 -0.7699811229
[50,] -0.6913897395 -0.8892436763
[51,] -0.4748026345 -0.6913897395
[52,] -0.4672654057 -0.4748026345
[53,] -0.2451453501 -0.4672654057
[54,] -0.1197021691 -0.2451453501
[55,] 0.0008154696 -0.1197021691
[56,] 0.2416686933 0.0008154696
[57,] 0.1420562937 0.2416686933
[58,] 0.2104207386 0.1420562937
[59,] 0.2852493808 0.2104207386
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 0.8944007730 0.9514795091
2 0.7331079335 0.8944007730
3 0.7345685655 0.7331079335
4 0.6593003076 0.7345685655
5 0.5148403887 0.6593003076
6 0.2456487277 0.5148403887
7 0.4553060121 0.2456487277
8 0.3467737746 0.4553060121
9 0.4913117165 0.3467737746
10 0.5694374765 0.4913117165
11 0.2732899872 0.5694374765
12 0.1140534415 0.2732899872
13 -0.0471613750 0.1140534415
14 0.2238667718 -0.0471613750
15 0.3103309691 0.2238667718
16 0.2640865797 0.3103309691
17 0.2399882535 0.2640865797
18 0.6900591459 0.2399882535
19 0.4341054568 0.6900591459
20 0.2856889965 0.4341054568
21 0.3710801622 0.2856889965
22 0.1727345943 0.3710801622
23 0.2604916310 0.1727345943
24 0.1463744277 0.2604916310
25 0.5593028223 0.1463744277
26 0.0936138257 0.5593028223
27 0.1919073783 0.0936138257
28 -0.0855589582 0.1919073783
29 -0.2236847182 -0.0855589582
30 -0.0595863920 -0.2236847182
31 -0.5854171733 -0.0595863920
32 -0.2048097651 -0.5854171733
33 -0.5162515200 -0.2048097651
34 -0.1445899573 -0.5162515200
35 -0.1965871050 -0.1445899573
36 -0.4419262553 -0.1965871050
37 -0.5172985439 -0.4419262553
38 -0.3591987916 -0.5172985439
39 -0.7620042783 -0.3591987916
40 -0.3705625235 -0.7620042783
41 -0.2859985739 -0.3705625235
42 -0.7564193125 -0.2859985739
43 -0.3048097651 -0.7564193125
44 -0.6693216994 -0.3048097651
45 -0.4881966524 -0.6693216994
46 -0.8080028522 -0.4881966524
47 -0.6224438941 -0.8080028522
48 -0.7699811229 -0.6224438941
49 -0.8892436763 -0.7699811229
50 -0.6913897395 -0.8892436763
51 -0.4748026345 -0.6913897395
52 -0.4672654057 -0.4748026345
53 -0.2451453501 -0.4672654057
54 -0.1197021691 -0.2451453501
55 0.0008154696 -0.1197021691
56 0.2416686933 0.0008154696
57 0.1420562937 0.2416686933
58 0.2104207386 0.1420562937
59 0.2852493808 0.2104207386
> 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/75zoj1227717844.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/8pdct1227717844.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/9pkdi1227717844.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/10l7sq1227717844.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/1174pb1227717844.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/12x01c1227717844.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/13ynr61227717844.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/14me6w1227717844.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/150pp11227717844.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/164tfw1227717844.tab")
+ }
>
> system("convert tmp/1qqxy1227717843.ps tmp/1qqxy1227717843.png")
> system("convert tmp/226171227717843.ps tmp/226171227717843.png")
> system("convert tmp/3u3l91227717843.ps tmp/3u3l91227717843.png")
> system("convert tmp/46pbw1227717843.ps tmp/46pbw1227717843.png")
> system("convert tmp/50rfk1227717843.ps tmp/50rfk1227717843.png")
> system("convert tmp/6ue201227717843.ps tmp/6ue201227717843.png")
> system("convert tmp/75zoj1227717844.ps tmp/75zoj1227717844.png")
> system("convert tmp/8pdct1227717844.ps tmp/8pdct1227717844.png")
> system("convert tmp/9pkdi1227717844.ps tmp/9pkdi1227717844.png")
> system("convert tmp/10l7sq1227717844.ps tmp/10l7sq1227717844.png")
>
>
> proc.time()
user system elapsed
2.407 1.527 3.557