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(17823.2,1.2218,17872,1.249,17420.4,1.2991,16704.4,1.3408,15991.2,1.3119,16583.6,1.3014,19123.5,1.3201,17838.7,1.2938,17209.4,1.2694,18586.5,1.2165,16258.1,1.2037,15141.6,1.2292,19202.1,1.2256,17746.5,1.2015,19090.1,1.1786,18040.3,1.1856,17515.5,1.2103,17751.8,1.1938,21072.4,1.202,17170,1.2271,19439.5,1.277,19795.4,1.265,17574.9,1.2684,16165.4,1.2811,19464.6,1.2727,19932.1,1.2611,19961.2,1.2881,17343.4,1.3213,18924.2,1.2999,18574.1,1.3074,21350.6,1.3242,18594.6,1.3516,19823.1,1.3511,20844.4,1.3419,19640.2,1.3716,17735.4,1.3622,19813.6,1.3896,22160,1.4227,20664.3,1.4684,17877.4,1.457,20906.5,1.4718,21164.1,1.4748,21374.4,1.5527,22952.3,1.5751,21343.5,1.5557,23899.3,1.5553,22392.9,1.577,18274.1,1.4975,22786.7,1.437,22321.5,1.3322,17842.2,1.2732,16373.5,1.3449,15993.8,1.3239,16446.1,1.2785,17729,1.305,16643,1.319,16196.7,1.365,18252.1,1.4016,17570.4,1.4088,15836.8,1.4268),dim=c(2,60),dimnames=list(c('UITV','EUDO'),1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('UITV','EUDO'),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 = '2'
> #'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
EUDO UITV M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 1.2218 17823.2 1 0 0 0 0 0 0 0 0 0 0
2 1.2490 17872.0 0 1 0 0 0 0 0 0 0 0 0
3 1.2991 17420.4 0 0 1 0 0 0 0 0 0 0 0
4 1.3408 16704.4 0 0 0 1 0 0 0 0 0 0 0
5 1.3119 15991.2 0 0 0 0 1 0 0 0 0 0 0
6 1.3014 16583.6 0 0 0 0 0 1 0 0 0 0 0
7 1.3201 19123.5 0 0 0 0 0 0 1 0 0 0 0
8 1.2938 17838.7 0 0 0 0 0 0 0 1 0 0 0
9 1.2694 17209.4 0 0 0 0 0 0 0 0 1 0 0
10 1.2165 18586.5 0 0 0 0 0 0 0 0 0 1 0
11 1.2037 16258.1 0 0 0 0 0 0 0 0 0 0 1
12 1.2292 15141.6 0 0 0 0 0 0 0 0 0 0 0
13 1.2256 19202.1 1 0 0 0 0 0 0 0 0 0 0
14 1.2015 17746.5 0 1 0 0 0 0 0 0 0 0 0
15 1.1786 19090.1 0 0 1 0 0 0 0 0 0 0 0
16 1.1856 18040.3 0 0 0 1 0 0 0 0 0 0 0
17 1.2103 17515.5 0 0 0 0 1 0 0 0 0 0 0
18 1.1938 17751.8 0 0 0 0 0 1 0 0 0 0 0
19 1.2020 21072.4 0 0 0 0 0 0 1 0 0 0 0
20 1.2271 17170.0 0 0 0 0 0 0 0 1 0 0 0
21 1.2770 19439.5 0 0 0 0 0 0 0 0 1 0 0
22 1.2650 19795.4 0 0 0 0 0 0 0 0 0 1 0
23 1.2684 17574.9 0 0 0 0 0 0 0 0 0 0 1
24 1.2811 16165.4 0 0 0 0 0 0 0 0 0 0 0
25 1.2727 19464.6 1 0 0 0 0 0 0 0 0 0 0
26 1.2611 19932.1 0 1 0 0 0 0 0 0 0 0 0
27 1.2881 19961.2 0 0 1 0 0 0 0 0 0 0 0
28 1.3213 17343.4 0 0 0 1 0 0 0 0 0 0 0
29 1.2999 18924.2 0 0 0 0 1 0 0 0 0 0 0
30 1.3074 18574.1 0 0 0 0 0 1 0 0 0 0 0
31 1.3242 21350.6 0 0 0 0 0 0 1 0 0 0 0
32 1.3516 18594.6 0 0 0 0 0 0 0 1 0 0 0
33 1.3511 19823.1 0 0 0 0 0 0 0 0 1 0 0
34 1.3419 20844.4 0 0 0 0 0 0 0 0 0 1 0
35 1.3716 19640.2 0 0 0 0 0 0 0 0 0 0 1
36 1.3622 17735.4 0 0 0 0 0 0 0 0 0 0 0
37 1.3896 19813.6 1 0 0 0 0 0 0 0 0 0 0
38 1.4227 22160.0 0 1 0 0 0 0 0 0 0 0 0
39 1.4684 20664.3 0 0 1 0 0 0 0 0 0 0 0
40 1.4570 17877.4 0 0 0 1 0 0 0 0 0 0 0
41 1.4718 20906.5 0 0 0 0 1 0 0 0 0 0 0
42 1.4748 21164.1 0 0 0 0 0 1 0 0 0 0 0
43 1.5527 21374.4 0 0 0 0 0 0 1 0 0 0 0
44 1.5751 22952.3 0 0 0 0 0 0 0 1 0 0 0
45 1.5557 21343.5 0 0 0 0 0 0 0 0 1 0 0
46 1.5553 23899.3 0 0 0 0 0 0 0 0 0 1 0
47 1.5770 22392.9 0 0 0 0 0 0 0 0 0 0 1
48 1.4975 18274.1 0 0 0 0 0 0 0 0 0 0 0
49 1.4370 22786.7 1 0 0 0 0 0 0 0 0 0 0
50 1.3322 22321.5 0 1 0 0 0 0 0 0 0 0 0
51 1.2732 17842.2 0 0 1 0 0 0 0 0 0 0 0
52 1.3449 16373.5 0 0 0 1 0 0 0 0 0 0 0
53 1.3239 15993.8 0 0 0 0 1 0 0 0 0 0 0
54 1.2785 16446.1 0 0 0 0 0 1 0 0 0 0 0
55 1.3050 17729.0 0 0 0 0 0 0 1 0 0 0 0
56 1.3190 16643.0 0 0 0 0 0 0 0 1 0 0 0
57 1.3650 16196.7 0 0 0 0 0 0 0 0 1 0 0
58 1.4016 18252.1 0 0 0 0 0 0 0 0 0 1 0
59 1.4088 17570.4 0 0 0 0 0 0 0 0 0 0 1
60 1.4268 15836.8 0 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) UITV M1 M2 M3 M4
6.882e-01 4.036e-05 -1.787e-01 -2.023e-01 -1.533e-01 -5.515e-02
M5 M6 M7 M8 M9 M10
-8.567e-02 -1.076e-01 -1.598e-01 -8.712e-02 -8.337e-02 -1.504e-01
M11
-7.646e-02
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.1768347 -0.0556203 -0.0002460 0.0550178 0.1616770
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.882e-01 1.097e-01 6.273 1.04e-07 ***
UITV 4.036e-05 6.207e-06 6.502 4.66e-08 ***
M1 -1.787e-01 5.610e-02 -3.184 0.002574 **
M2 -2.023e-01 5.653e-02 -3.579 0.000814 ***
M3 -1.533e-01 5.451e-02 -2.813 0.007147 **
M4 -5.515e-02 5.265e-02 -1.048 0.300173
M5 -8.567e-02 5.306e-02 -1.615 0.113081
M6 -1.076e-01 5.329e-02 -2.020 0.049113 *
M7 -1.598e-01 5.681e-02 -2.812 0.007152 **
M8 -8.712e-02 5.396e-02 -1.615 0.113091
M9 -8.337e-02 5.420e-02 -1.538 0.130719
M10 -1.504e-01 5.717e-02 -2.631 0.011477 *
M11 -7.646e-02 5.403e-02 -1.415 0.163590
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.08301 on 47 degrees of freedom
Multiple R-squared: 0.5035, Adjusted R-squared: 0.3767
F-statistic: 3.971 on 12 and 47 DF, p-value: 0.0003095
> 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.22837053 0.45674106 0.771629469
[2,] 0.11435568 0.22871135 0.885644323
[3,] 0.06162334 0.12324668 0.938376662
[4,] 0.04282383 0.08564767 0.957176165
[5,] 0.07878401 0.15756802 0.921215992
[6,] 0.15991609 0.31983218 0.840083912
[7,] 0.19221741 0.38443483 0.807782585
[8,] 0.23389508 0.46779016 0.766104919
[9,] 0.23964712 0.47929424 0.760352880
[10,] 0.22794558 0.45589116 0.772054418
[11,] 0.20662556 0.41325112 0.793374439
[12,] 0.21307913 0.42615826 0.786920869
[13,] 0.19416212 0.38832424 0.805837881
[14,] 0.19950226 0.39900453 0.800497737
[15,] 0.17736730 0.35473459 0.822632705
[16,] 0.28778295 0.57556590 0.712217052
[17,] 0.28379691 0.56759383 0.716203087
[18,] 0.42246959 0.84493918 0.577530411
[19,] 0.57317455 0.85365090 0.426825449
[20,] 0.69988525 0.60022950 0.300114749
[21,] 0.85432359 0.29135282 0.145676411
[22,] 0.88546008 0.22907984 0.114539920
[23,] 0.92755108 0.14489783 0.072448916
[24,] 0.96416795 0.07166409 0.035832047
[25,] 0.96668945 0.06662109 0.033310547
[26,] 0.94507507 0.10984986 0.054924931
[27,] 0.89348067 0.21303866 0.106519328
[28,] 0.99078025 0.01843949 0.009219747
[29,] 0.98745989 0.02508022 0.012540108
> postscript(file="/var/www/html/rcomp/tmp/11x7f1258917583.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/2fk2b1258917583.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/3d4791258917583.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/4bexj1258917583.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/5z1ti1258917583.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.0070311644 0.0418420810 0.0611943910 0.0336180030 0.0640138821
6 7 8 9 10
0.0515787070 0.0199200662 -0.0271920000 -0.0299472269 -0.0713928068
11 12 13 14 15
-0.0641610278 -0.0700637082 -0.0588815591 -0.0005929219 -0.1266922680
16 17 18 19 20
-0.1754969743 -0.0991046447 -0.1031681427 -0.1768346979 -0.0669042424
21 22 23 24 25
-0.1123508133 -0.0716822492 -0.0526051570 -0.0594827844 -0.0223756766
26 27 28 29 30
-0.0292005531 -0.0523485946 -0.0116711058 -0.0663577240 -0.0227549725
31 32 33 34 35
-0.0658624445 0.0001009775 -0.0537323503 -0.0371183607 -0.0327576555
36 37 38 39 40
-0.0417456966 0.0804391919 0.0424846493 0.0995753139 0.1024774323
41 42 43 44 45
0.0255395366 0.0401164016 0.1616770221 0.0477305554 0.0895065212
46 47 48 49 50
0.0529903268 0.0615473837 0.0718131564 0.0078492082 -0.0545332553
51 52 53 54 55
0.0182711576 0.0510726449 0.0759089499 0.0342280066 0.0611000542
56 57 58 59 60
0.0462647096 0.1065238694 0.1272030899 0.0879764565 0.0994790328
> postscript(file="/var/www/html/rcomp/tmp/6b42g1258917583.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.0070311644 NA
1 0.0418420810 -0.0070311644
2 0.0611943910 0.0418420810
3 0.0336180030 0.0611943910
4 0.0640138821 0.0336180030
5 0.0515787070 0.0640138821
6 0.0199200662 0.0515787070
7 -0.0271920000 0.0199200662
8 -0.0299472269 -0.0271920000
9 -0.0713928068 -0.0299472269
10 -0.0641610278 -0.0713928068
11 -0.0700637082 -0.0641610278
12 -0.0588815591 -0.0700637082
13 -0.0005929219 -0.0588815591
14 -0.1266922680 -0.0005929219
15 -0.1754969743 -0.1266922680
16 -0.0991046447 -0.1754969743
17 -0.1031681427 -0.0991046447
18 -0.1768346979 -0.1031681427
19 -0.0669042424 -0.1768346979
20 -0.1123508133 -0.0669042424
21 -0.0716822492 -0.1123508133
22 -0.0526051570 -0.0716822492
23 -0.0594827844 -0.0526051570
24 -0.0223756766 -0.0594827844
25 -0.0292005531 -0.0223756766
26 -0.0523485946 -0.0292005531
27 -0.0116711058 -0.0523485946
28 -0.0663577240 -0.0116711058
29 -0.0227549725 -0.0663577240
30 -0.0658624445 -0.0227549725
31 0.0001009775 -0.0658624445
32 -0.0537323503 0.0001009775
33 -0.0371183607 -0.0537323503
34 -0.0327576555 -0.0371183607
35 -0.0417456966 -0.0327576555
36 0.0804391919 -0.0417456966
37 0.0424846493 0.0804391919
38 0.0995753139 0.0424846493
39 0.1024774323 0.0995753139
40 0.0255395366 0.1024774323
41 0.0401164016 0.0255395366
42 0.1616770221 0.0401164016
43 0.0477305554 0.1616770221
44 0.0895065212 0.0477305554
45 0.0529903268 0.0895065212
46 0.0615473837 0.0529903268
47 0.0718131564 0.0615473837
48 0.0078492082 0.0718131564
49 -0.0545332553 0.0078492082
50 0.0182711576 -0.0545332553
51 0.0510726449 0.0182711576
52 0.0759089499 0.0510726449
53 0.0342280066 0.0759089499
54 0.0611000542 0.0342280066
55 0.0462647096 0.0611000542
56 0.1065238694 0.0462647096
57 0.1272030899 0.1065238694
58 0.0879764565 0.1272030899
59 0.0994790328 0.0879764565
60 NA 0.0994790328
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 0.0418420810 -0.0070311644
[2,] 0.0611943910 0.0418420810
[3,] 0.0336180030 0.0611943910
[4,] 0.0640138821 0.0336180030
[5,] 0.0515787070 0.0640138821
[6,] 0.0199200662 0.0515787070
[7,] -0.0271920000 0.0199200662
[8,] -0.0299472269 -0.0271920000
[9,] -0.0713928068 -0.0299472269
[10,] -0.0641610278 -0.0713928068
[11,] -0.0700637082 -0.0641610278
[12,] -0.0588815591 -0.0700637082
[13,] -0.0005929219 -0.0588815591
[14,] -0.1266922680 -0.0005929219
[15,] -0.1754969743 -0.1266922680
[16,] -0.0991046447 -0.1754969743
[17,] -0.1031681427 -0.0991046447
[18,] -0.1768346979 -0.1031681427
[19,] -0.0669042424 -0.1768346979
[20,] -0.1123508133 -0.0669042424
[21,] -0.0716822492 -0.1123508133
[22,] -0.0526051570 -0.0716822492
[23,] -0.0594827844 -0.0526051570
[24,] -0.0223756766 -0.0594827844
[25,] -0.0292005531 -0.0223756766
[26,] -0.0523485946 -0.0292005531
[27,] -0.0116711058 -0.0523485946
[28,] -0.0663577240 -0.0116711058
[29,] -0.0227549725 -0.0663577240
[30,] -0.0658624445 -0.0227549725
[31,] 0.0001009775 -0.0658624445
[32,] -0.0537323503 0.0001009775
[33,] -0.0371183607 -0.0537323503
[34,] -0.0327576555 -0.0371183607
[35,] -0.0417456966 -0.0327576555
[36,] 0.0804391919 -0.0417456966
[37,] 0.0424846493 0.0804391919
[38,] 0.0995753139 0.0424846493
[39,] 0.1024774323 0.0995753139
[40,] 0.0255395366 0.1024774323
[41,] 0.0401164016 0.0255395366
[42,] 0.1616770221 0.0401164016
[43,] 0.0477305554 0.1616770221
[44,] 0.0895065212 0.0477305554
[45,] 0.0529903268 0.0895065212
[46,] 0.0615473837 0.0529903268
[47,] 0.0718131564 0.0615473837
[48,] 0.0078492082 0.0718131564
[49,] -0.0545332553 0.0078492082
[50,] 0.0182711576 -0.0545332553
[51,] 0.0510726449 0.0182711576
[52,] 0.0759089499 0.0510726449
[53,] 0.0342280066 0.0759089499
[54,] 0.0611000542 0.0342280066
[55,] 0.0462647096 0.0611000542
[56,] 0.1065238694 0.0462647096
[57,] 0.1272030899 0.1065238694
[58,] 0.0879764565 0.1272030899
[59,] 0.0994790328 0.0879764565
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 0.0418420810 -0.0070311644
2 0.0611943910 0.0418420810
3 0.0336180030 0.0611943910
4 0.0640138821 0.0336180030
5 0.0515787070 0.0640138821
6 0.0199200662 0.0515787070
7 -0.0271920000 0.0199200662
8 -0.0299472269 -0.0271920000
9 -0.0713928068 -0.0299472269
10 -0.0641610278 -0.0713928068
11 -0.0700637082 -0.0641610278
12 -0.0588815591 -0.0700637082
13 -0.0005929219 -0.0588815591
14 -0.1266922680 -0.0005929219
15 -0.1754969743 -0.1266922680
16 -0.0991046447 -0.1754969743
17 -0.1031681427 -0.0991046447
18 -0.1768346979 -0.1031681427
19 -0.0669042424 -0.1768346979
20 -0.1123508133 -0.0669042424
21 -0.0716822492 -0.1123508133
22 -0.0526051570 -0.0716822492
23 -0.0594827844 -0.0526051570
24 -0.0223756766 -0.0594827844
25 -0.0292005531 -0.0223756766
26 -0.0523485946 -0.0292005531
27 -0.0116711058 -0.0523485946
28 -0.0663577240 -0.0116711058
29 -0.0227549725 -0.0663577240
30 -0.0658624445 -0.0227549725
31 0.0001009775 -0.0658624445
32 -0.0537323503 0.0001009775
33 -0.0371183607 -0.0537323503
34 -0.0327576555 -0.0371183607
35 -0.0417456966 -0.0327576555
36 0.0804391919 -0.0417456966
37 0.0424846493 0.0804391919
38 0.0995753139 0.0424846493
39 0.1024774323 0.0995753139
40 0.0255395366 0.1024774323
41 0.0401164016 0.0255395366
42 0.1616770221 0.0401164016
43 0.0477305554 0.1616770221
44 0.0895065212 0.0477305554
45 0.0529903268 0.0895065212
46 0.0615473837 0.0529903268
47 0.0718131564 0.0615473837
48 0.0078492082 0.0718131564
49 -0.0545332553 0.0078492082
50 0.0182711576 -0.0545332553
51 0.0510726449 0.0182711576
52 0.0759089499 0.0510726449
53 0.0342280066 0.0759089499
54 0.0611000542 0.0342280066
55 0.0462647096 0.0611000542
56 0.1065238694 0.0462647096
57 0.1272030899 0.1065238694
58 0.0879764565 0.1272030899
59 0.0994790328 0.0879764565
> 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/7izci1258917583.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/81mdi1258917583.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/9rmji1258917583.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/10ka671258917583.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/11olkq1258917583.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/12shuz1258917583.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/13plxn1258917583.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/14l1vh1258917584.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/15qn051258917584.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/16oyek1258917584.tab")
+ }
>
> system("convert tmp/11x7f1258917583.ps tmp/11x7f1258917583.png")
> system("convert tmp/2fk2b1258917583.ps tmp/2fk2b1258917583.png")
> system("convert tmp/3d4791258917583.ps tmp/3d4791258917583.png")
> system("convert tmp/4bexj1258917583.ps tmp/4bexj1258917583.png")
> system("convert tmp/5z1ti1258917583.ps tmp/5z1ti1258917583.png")
> system("convert tmp/6b42g1258917583.ps tmp/6b42g1258917583.png")
> system("convert tmp/7izci1258917583.ps tmp/7izci1258917583.png")
> system("convert tmp/81mdi1258917583.ps tmp/81mdi1258917583.png")
> system("convert tmp/9rmji1258917583.ps tmp/9rmji1258917583.png")
> system("convert tmp/10ka671258917583.ps tmp/10ka671258917583.png")
>
>
> proc.time()
user system elapsed
2.454 1.585 2.886