R version 2.7.0 (2008-04-22)
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(101.5,1,100.7,1,110.6,1,96.8,1,100.0,1,104.8,1,86.8,1,92.0,1,100.2,1,106.6,1,102.1,1,93.7,1,97.6,1,96.9,1,105.6,1,102.8,1,101.7,1,104.2,1,92.7,1,91.9,1,106.5,1,112.3,1,102.8,1,96.5,1,101.0,0,98.9,0,105.1,0,103.0,0,99.0,0,104.3,0,94.6,0,90.4,0,108.9,0,111.4,0,100.8,0,102.5,0,98.2,0,98.7,0,113.3,0,104.6,0,99.3,0,111.8,0,97.3,0,97.7,0,115.6,0,111.9,0,107.0,0,107.1,0,100.6,0,99.2,0,108.4,0,103.0,0,99.8,0,115.0,0,90.8,0,95.9,0,114.4,0,108.2,0,112.6,0,109.1,0,105.0,0,105.0,0,118.5,0,103.7,0,112.5,0,116.6,0,96.6,0,101.9,0,116.5,0,119.3,0,115.4,0,108.5,0,111.5,0,108.8,0,121.8,0,109.6,0,112.2,0,119.6,0,104.1,0,105.3,0,115.0,0,124.1,0,116.8,0,107.5,0,115.6,0,116.2,0,116.3,0,119.0,0,111.9,0,118.6,0,106.9,0,103.2,0),dim=c(2,92),dimnames=list(c('Y','X'),1:92))
> y <- array(NA,dim=c(2,92),dimnames=list(c('Y','X'),1:92))
> 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
Y X M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11
1 101.5 1 1 0 0 0 0 0 0 0 0 0 0
2 100.7 1 0 1 0 0 0 0 0 0 0 0 0
3 110.6 1 0 0 1 0 0 0 0 0 0 0 0
4 96.8 1 0 0 0 1 0 0 0 0 0 0 0
5 100.0 1 0 0 0 0 1 0 0 0 0 0 0
6 104.8 1 0 0 0 0 0 1 0 0 0 0 0
7 86.8 1 0 0 0 0 0 0 1 0 0 0 0
8 92.0 1 0 0 0 0 0 0 0 1 0 0 0
9 100.2 1 0 0 0 0 0 0 0 0 1 0 0
10 106.6 1 0 0 0 0 0 0 0 0 0 1 0
11 102.1 1 0 0 0 0 0 0 0 0 0 0 1
12 93.7 1 0 0 0 0 0 0 0 0 0 0 0
13 97.6 1 1 0 0 0 0 0 0 0 0 0 0
14 96.9 1 0 1 0 0 0 0 0 0 0 0 0
15 105.6 1 0 0 1 0 0 0 0 0 0 0 0
16 102.8 1 0 0 0 1 0 0 0 0 0 0 0
17 101.7 1 0 0 0 0 1 0 0 0 0 0 0
18 104.2 1 0 0 0 0 0 1 0 0 0 0 0
19 92.7 1 0 0 0 0 0 0 1 0 0 0 0
20 91.9 1 0 0 0 0 0 0 0 1 0 0 0
21 106.5 1 0 0 0 0 0 0 0 0 1 0 0
22 112.3 1 0 0 0 0 0 0 0 0 0 1 0
23 102.8 1 0 0 0 0 0 0 0 0 0 0 1
24 96.5 1 0 0 0 0 0 0 0 0 0 0 0
25 101.0 0 1 0 0 0 0 0 0 0 0 0 0
26 98.9 0 0 1 0 0 0 0 0 0 0 0 0
27 105.1 0 0 0 1 0 0 0 0 0 0 0 0
28 103.0 0 0 0 0 1 0 0 0 0 0 0 0
29 99.0 0 0 0 0 0 1 0 0 0 0 0 0
30 104.3 0 0 0 0 0 0 1 0 0 0 0 0
31 94.6 0 0 0 0 0 0 0 1 0 0 0 0
32 90.4 0 0 0 0 0 0 0 0 1 0 0 0
33 108.9 0 0 0 0 0 0 0 0 0 1 0 0
34 111.4 0 0 0 0 0 0 0 0 0 0 1 0
35 100.8 0 0 0 0 0 0 0 0 0 0 0 1
36 102.5 0 0 0 0 0 0 0 0 0 0 0 0
37 98.2 0 1 0 0 0 0 0 0 0 0 0 0
38 98.7 0 0 1 0 0 0 0 0 0 0 0 0
39 113.3 0 0 0 1 0 0 0 0 0 0 0 0
40 104.6 0 0 0 0 1 0 0 0 0 0 0 0
41 99.3 0 0 0 0 0 1 0 0 0 0 0 0
42 111.8 0 0 0 0 0 0 1 0 0 0 0 0
43 97.3 0 0 0 0 0 0 0 1 0 0 0 0
44 97.7 0 0 0 0 0 0 0 0 1 0 0 0
45 115.6 0 0 0 0 0 0 0 0 0 1 0 0
46 111.9 0 0 0 0 0 0 0 0 0 0 1 0
47 107.0 0 0 0 0 0 0 0 0 0 0 0 1
48 107.1 0 0 0 0 0 0 0 0 0 0 0 0
49 100.6 0 1 0 0 0 0 0 0 0 0 0 0
50 99.2 0 0 1 0 0 0 0 0 0 0 0 0
51 108.4 0 0 0 1 0 0 0 0 0 0 0 0
52 103.0 0 0 0 0 1 0 0 0 0 0 0 0
53 99.8 0 0 0 0 0 1 0 0 0 0 0 0
54 115.0 0 0 0 0 0 0 1 0 0 0 0 0
55 90.8 0 0 0 0 0 0 0 1 0 0 0 0
56 95.9 0 0 0 0 0 0 0 0 1 0 0 0
57 114.4 0 0 0 0 0 0 0 0 0 1 0 0
58 108.2 0 0 0 0 0 0 0 0 0 0 1 0
59 112.6 0 0 0 0 0 0 0 0 0 0 0 1
60 109.1 0 0 0 0 0 0 0 0 0 0 0 0
61 105.0 0 1 0 0 0 0 0 0 0 0 0 0
62 105.0 0 0 1 0 0 0 0 0 0 0 0 0
63 118.5 0 0 0 1 0 0 0 0 0 0 0 0
64 103.7 0 0 0 0 1 0 0 0 0 0 0 0
65 112.5 0 0 0 0 0 1 0 0 0 0 0 0
66 116.6 0 0 0 0 0 0 1 0 0 0 0 0
67 96.6 0 0 0 0 0 0 0 1 0 0 0 0
68 101.9 0 0 0 0 0 0 0 0 1 0 0 0
69 116.5 0 0 0 0 0 0 0 0 0 1 0 0
70 119.3 0 0 0 0 0 0 0 0 0 0 1 0
71 115.4 0 0 0 0 0 0 0 0 0 0 0 1
72 108.5 0 0 0 0 0 0 0 0 0 0 0 0
73 111.5 0 1 0 0 0 0 0 0 0 0 0 0
74 108.8 0 0 1 0 0 0 0 0 0 0 0 0
75 121.8 0 0 0 1 0 0 0 0 0 0 0 0
76 109.6 0 0 0 0 1 0 0 0 0 0 0 0
77 112.2 0 0 0 0 0 1 0 0 0 0 0 0
78 119.6 0 0 0 0 0 0 1 0 0 0 0 0
79 104.1 0 0 0 0 0 0 0 1 0 0 0 0
80 105.3 0 0 0 0 0 0 0 0 1 0 0 0
81 115.0 0 0 0 0 0 0 0 0 0 1 0 0
82 124.1 0 0 0 0 0 0 0 0 0 0 1 0
83 116.8 0 0 0 0 0 0 0 0 0 0 0 1
84 107.5 0 0 0 0 0 0 0 0 0 0 0 0
85 115.6 0 1 0 0 0 0 0 0 0 0 0 0
86 116.2 0 0 1 0 0 0 0 0 0 0 0 0
87 116.3 0 0 0 1 0 0 0 0 0 0 0 0
88 119.0 0 0 0 0 1 0 0 0 0 0 0 0
89 111.9 0 0 0 0 0 1 0 0 0 0 0 0
90 118.6 0 0 0 0 0 0 1 0 0 0 0 0
91 106.9 0 0 0 0 0 0 0 1 0 0 0 0
92 103.2 0 0 0 0 0 0 0 0 1 0 0 0
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X M1 M2 M3 M4
105.7232 -7.5812 0.0471 -0.7779 8.6221 1.4846
M5 M6 M7 M8 M9 M10
0.7221 8.0346 -7.6029 -6.5404 7.4571 9.8429
M11
4.6571
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-9.58036 -3.68432 0.02779 3.34615 11.79219
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 105.7232 2.0327 52.012 < 2e-16 ***
X -7.5812 1.2576 -6.028 5.01e-08 ***
M1 0.0471 2.7399 0.017 0.986329
M2 -0.7779 2.7399 -0.284 0.777217
M3 8.6221 2.7399 3.147 0.002328 **
M4 1.4846 2.7399 0.542 0.589449
M5 0.7221 2.7399 0.264 0.792813
M6 8.0346 2.7399 2.932 0.004398 **
M7 -7.6029 2.7399 -2.775 0.006890 **
M8 -6.5404 2.7399 -2.387 0.019372 *
M9 7.4571 2.8294 2.636 0.010105 *
M10 9.8429 2.8294 3.479 0.000822 ***
M11 4.6571 2.8294 1.646 0.103736
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.293 on 79 degrees of freedom
Multiple R-squared: 0.6291, Adjusted R-squared: 0.5727
F-statistic: 11.17 on 12 and 79 DF, p-value: 1.185e-12
> 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.2886102865 0.5772205730 0.71138971
[2,] 0.1579794776 0.3159589552 0.84202052
[3,] 0.0769677286 0.1539354573 0.92303227
[4,] 0.0724429442 0.1448858883 0.92755706
[5,] 0.0346804160 0.0693608320 0.96531958
[6,] 0.0359765745 0.0719531490 0.96402343
[7,] 0.0325585079 0.0651170157 0.96744149
[8,] 0.0164216372 0.0328432744 0.98357836
[9,] 0.0091285484 0.0182570969 0.99087145
[10,] 0.0044208573 0.0088417146 0.99557914
[11,] 0.0022147499 0.0044294997 0.99778525
[12,] 0.0017065454 0.0034130909 0.99829345
[13,] 0.0011134954 0.0022269908 0.99888650
[14,] 0.0006695745 0.0013391491 0.99933043
[15,] 0.0004509884 0.0009019768 0.99954901
[16,] 0.0004320202 0.0008640405 0.99956798
[17,] 0.0003253878 0.0006507756 0.99967461
[18,] 0.0003757795 0.0007515590 0.99962422
[19,] 0.0001882344 0.0003764687 0.99981177
[20,] 0.0001992175 0.0003984350 0.99980078
[21,] 0.0003673537 0.0007347074 0.99963265
[22,] 0.0003589139 0.0007178278 0.99964109
[23,] 0.0002613192 0.0005226384 0.99973868
[24,] 0.0003286050 0.0006572100 0.99967139
[25,] 0.0002223525 0.0004447050 0.99977765
[26,] 0.0002187600 0.0004375199 0.99978124
[27,] 0.0004142069 0.0008284139 0.99958579
[28,] 0.0003856334 0.0007712668 0.99961437
[29,] 0.0003916668 0.0007833336 0.99960833
[30,] 0.0010934006 0.0021868013 0.99890660
[31,] 0.0007018396 0.0014036792 0.99929816
[32,] 0.0007023840 0.0014047681 0.99929762
[33,] 0.0011091780 0.0022183560 0.99889082
[34,] 0.0012384984 0.0024769968 0.99876150
[35,] 0.0016829628 0.0033659257 0.99831704
[36,] 0.0024358091 0.0048716181 0.99756419
[37,] 0.0023547278 0.0047094557 0.99764527
[38,] 0.0060355341 0.0120710682 0.99396447
[39,] 0.0088495948 0.0176991896 0.99115041
[40,] 0.0247975836 0.0495951673 0.97520242
[41,] 0.0296270115 0.0592540230 0.97037299
[42,] 0.0263462630 0.0526925260 0.97365374
[43,] 0.1135148898 0.2270297796 0.88648511
[44,] 0.1431960055 0.2863920110 0.85680399
[45,] 0.1506663865 0.3013327730 0.84933361
[46,] 0.2245969721 0.4491939442 0.77540303
[47,] 0.3030492337 0.6060984674 0.69695077
[48,] 0.3143908920 0.6287817841 0.68560911
[49,] 0.5656839241 0.8686321517 0.43431608
[50,] 0.6062696028 0.7874607943 0.39373040
[51,] 0.5835714602 0.8328570796 0.41642854
[52,] 0.7761358640 0.4477282720 0.22386414
[53,] 0.7454357789 0.5091284422 0.25456422
[54,] 0.6850508826 0.6298982349 0.31494912
[55,] 0.6929305736 0.6141388528 0.30706943
[56,] 0.6372406406 0.7255187188 0.36275936
[57,] 0.5393048130 0.9213903739 0.46069519
[58,] 0.5163351773 0.9673296454 0.48366482
[59,] 0.6055009761 0.7889980479 0.39449902
[60,] 0.6195764229 0.7608471542 0.38042358
[61,] 0.9507059847 0.0985880306 0.04929402
> postscript(file="/var/www/html/rcomp/tmp/1clow1229040817.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/2ywit1229040817.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/35fkm1229040818.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/4qirp1229040818.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/5shwo1229040818.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 = 92
Frequency = 1
1 2 3 4 5
3.3109375000 3.3359375000 3.8359375000 -2.8265625000 1.1359375000
6 7 8 9 10
-1.3765625000 -3.7390625000 0.3984375000 -5.3991071429 -1.3848214286
11 12 13 14 15
-0.6991071429 -4.4419642857 -0.5890625000 -0.4640625000 -1.1640625000
16 17 18 19 20
3.1734375000 2.8359375000 -1.9765625000 2.1609375000 0.2984375000
21 22 23 24 25
0.9008928571 4.3151785714 0.0008928571 -1.6419642857 -4.7703125000
26 27 28 29 30
-6.0453125000 -9.2453125000 -4.2078125000 -7.4453125000 -9.4578125000
31 32 33 34 35
-3.5203125000 -8.7828125000 -4.2803571429 -4.1660714286 -9.5803571429
36 37 38 39 40
-3.2232142857 -7.5703125000 -6.2453125000 -1.0453125000 -2.6078125000
41 42 43 44 45
-7.1453125000 -1.9578125000 -0.8203125000 -1.4828125000 2.4196428571
46 47 48 49 50
-3.6660714286 -3.3803571429 1.3767857143 -5.1703125000 -5.7453125000
51 52 53 54 55
-5.9453125000 -4.2078125000 -6.6453125000 1.2421875000 -7.3203125000
56 57 58 59 60
-3.2828125000 1.2196428571 -7.3660714286 2.2196428571 3.3767857143
61 62 63 64 65
-0.7703125000 0.0546875000 4.1546875000 -3.5078125000 6.0546875000
66 67 68 69 70
2.8421875000 -1.5203125000 2.7171875000 3.3196428571 3.7339285714
71 72 73 74 75
5.0196428571 2.7767857143 5.7296875000 3.8546875000 7.4546875000
76 77 78 79 80
2.3921875000 5.7546875000 5.8421875000 5.9796875000 6.1171875000
81 82 83 84 85
1.8196428571 8.5339285714 6.4196428571 1.7767857143 9.8296875000
86 87 88 89 90
11.2546875000 1.9546875000 11.7921875000 5.4546875000 4.8421875000
91 92
8.7796875000 4.0171875000
> postscript(file="/var/www/html/rcomp/tmp/6kevc1229040818.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 = 92
Frequency = 1
lag(myerror, k = 1) myerror
0 3.3109375000 NA
1 3.3359375000 3.3109375000
2 3.8359375000 3.3359375000
3 -2.8265625000 3.8359375000
4 1.1359375000 -2.8265625000
5 -1.3765625000 1.1359375000
6 -3.7390625000 -1.3765625000
7 0.3984375000 -3.7390625000
8 -5.3991071429 0.3984375000
9 -1.3848214286 -5.3991071429
10 -0.6991071429 -1.3848214286
11 -4.4419642857 -0.6991071429
12 -0.5890625000 -4.4419642857
13 -0.4640625000 -0.5890625000
14 -1.1640625000 -0.4640625000
15 3.1734375000 -1.1640625000
16 2.8359375000 3.1734375000
17 -1.9765625000 2.8359375000
18 2.1609375000 -1.9765625000
19 0.2984375000 2.1609375000
20 0.9008928571 0.2984375000
21 4.3151785714 0.9008928571
22 0.0008928571 4.3151785714
23 -1.6419642857 0.0008928571
24 -4.7703125000 -1.6419642857
25 -6.0453125000 -4.7703125000
26 -9.2453125000 -6.0453125000
27 -4.2078125000 -9.2453125000
28 -7.4453125000 -4.2078125000
29 -9.4578125000 -7.4453125000
30 -3.5203125000 -9.4578125000
31 -8.7828125000 -3.5203125000
32 -4.2803571429 -8.7828125000
33 -4.1660714286 -4.2803571429
34 -9.5803571429 -4.1660714286
35 -3.2232142857 -9.5803571429
36 -7.5703125000 -3.2232142857
37 -6.2453125000 -7.5703125000
38 -1.0453125000 -6.2453125000
39 -2.6078125000 -1.0453125000
40 -7.1453125000 -2.6078125000
41 -1.9578125000 -7.1453125000
42 -0.8203125000 -1.9578125000
43 -1.4828125000 -0.8203125000
44 2.4196428571 -1.4828125000
45 -3.6660714286 2.4196428571
46 -3.3803571429 -3.6660714286
47 1.3767857143 -3.3803571429
48 -5.1703125000 1.3767857143
49 -5.7453125000 -5.1703125000
50 -5.9453125000 -5.7453125000
51 -4.2078125000 -5.9453125000
52 -6.6453125000 -4.2078125000
53 1.2421875000 -6.6453125000
54 -7.3203125000 1.2421875000
55 -3.2828125000 -7.3203125000
56 1.2196428571 -3.2828125000
57 -7.3660714286 1.2196428571
58 2.2196428571 -7.3660714286
59 3.3767857143 2.2196428571
60 -0.7703125000 3.3767857143
61 0.0546875000 -0.7703125000
62 4.1546875000 0.0546875000
63 -3.5078125000 4.1546875000
64 6.0546875000 -3.5078125000
65 2.8421875000 6.0546875000
66 -1.5203125000 2.8421875000
67 2.7171875000 -1.5203125000
68 3.3196428571 2.7171875000
69 3.7339285714 3.3196428571
70 5.0196428571 3.7339285714
71 2.7767857143 5.0196428571
72 5.7296875000 2.7767857143
73 3.8546875000 5.7296875000
74 7.4546875000 3.8546875000
75 2.3921875000 7.4546875000
76 5.7546875000 2.3921875000
77 5.8421875000 5.7546875000
78 5.9796875000 5.8421875000
79 6.1171875000 5.9796875000
80 1.8196428571 6.1171875000
81 8.5339285714 1.8196428571
82 6.4196428571 8.5339285714
83 1.7767857143 6.4196428571
84 9.8296875000 1.7767857143
85 11.2546875000 9.8296875000
86 1.9546875000 11.2546875000
87 11.7921875000 1.9546875000
88 5.4546875000 11.7921875000
89 4.8421875000 5.4546875000
90 8.7796875000 4.8421875000
91 4.0171875000 8.7796875000
92 NA 4.0171875000
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 3.3359375000 3.3109375000
[2,] 3.8359375000 3.3359375000
[3,] -2.8265625000 3.8359375000
[4,] 1.1359375000 -2.8265625000
[5,] -1.3765625000 1.1359375000
[6,] -3.7390625000 -1.3765625000
[7,] 0.3984375000 -3.7390625000
[8,] -5.3991071429 0.3984375000
[9,] -1.3848214286 -5.3991071429
[10,] -0.6991071429 -1.3848214286
[11,] -4.4419642857 -0.6991071429
[12,] -0.5890625000 -4.4419642857
[13,] -0.4640625000 -0.5890625000
[14,] -1.1640625000 -0.4640625000
[15,] 3.1734375000 -1.1640625000
[16,] 2.8359375000 3.1734375000
[17,] -1.9765625000 2.8359375000
[18,] 2.1609375000 -1.9765625000
[19,] 0.2984375000 2.1609375000
[20,] 0.9008928571 0.2984375000
[21,] 4.3151785714 0.9008928571
[22,] 0.0008928571 4.3151785714
[23,] -1.6419642857 0.0008928571
[24,] -4.7703125000 -1.6419642857
[25,] -6.0453125000 -4.7703125000
[26,] -9.2453125000 -6.0453125000
[27,] -4.2078125000 -9.2453125000
[28,] -7.4453125000 -4.2078125000
[29,] -9.4578125000 -7.4453125000
[30,] -3.5203125000 -9.4578125000
[31,] -8.7828125000 -3.5203125000
[32,] -4.2803571429 -8.7828125000
[33,] -4.1660714286 -4.2803571429
[34,] -9.5803571429 -4.1660714286
[35,] -3.2232142857 -9.5803571429
[36,] -7.5703125000 -3.2232142857
[37,] -6.2453125000 -7.5703125000
[38,] -1.0453125000 -6.2453125000
[39,] -2.6078125000 -1.0453125000
[40,] -7.1453125000 -2.6078125000
[41,] -1.9578125000 -7.1453125000
[42,] -0.8203125000 -1.9578125000
[43,] -1.4828125000 -0.8203125000
[44,] 2.4196428571 -1.4828125000
[45,] -3.6660714286 2.4196428571
[46,] -3.3803571429 -3.6660714286
[47,] 1.3767857143 -3.3803571429
[48,] -5.1703125000 1.3767857143
[49,] -5.7453125000 -5.1703125000
[50,] -5.9453125000 -5.7453125000
[51,] -4.2078125000 -5.9453125000
[52,] -6.6453125000 -4.2078125000
[53,] 1.2421875000 -6.6453125000
[54,] -7.3203125000 1.2421875000
[55,] -3.2828125000 -7.3203125000
[56,] 1.2196428571 -3.2828125000
[57,] -7.3660714286 1.2196428571
[58,] 2.2196428571 -7.3660714286
[59,] 3.3767857143 2.2196428571
[60,] -0.7703125000 3.3767857143
[61,] 0.0546875000 -0.7703125000
[62,] 4.1546875000 0.0546875000
[63,] -3.5078125000 4.1546875000
[64,] 6.0546875000 -3.5078125000
[65,] 2.8421875000 6.0546875000
[66,] -1.5203125000 2.8421875000
[67,] 2.7171875000 -1.5203125000
[68,] 3.3196428571 2.7171875000
[69,] 3.7339285714 3.3196428571
[70,] 5.0196428571 3.7339285714
[71,] 2.7767857143 5.0196428571
[72,] 5.7296875000 2.7767857143
[73,] 3.8546875000 5.7296875000
[74,] 7.4546875000 3.8546875000
[75,] 2.3921875000 7.4546875000
[76,] 5.7546875000 2.3921875000
[77,] 5.8421875000 5.7546875000
[78,] 5.9796875000 5.8421875000
[79,] 6.1171875000 5.9796875000
[80,] 1.8196428571 6.1171875000
[81,] 8.5339285714 1.8196428571
[82,] 6.4196428571 8.5339285714
[83,] 1.7767857143 6.4196428571
[84,] 9.8296875000 1.7767857143
[85,] 11.2546875000 9.8296875000
[86,] 1.9546875000 11.2546875000
[87,] 11.7921875000 1.9546875000
[88,] 5.4546875000 11.7921875000
[89,] 4.8421875000 5.4546875000
[90,] 8.7796875000 4.8421875000
[91,] 4.0171875000 8.7796875000
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 3.3359375000 3.3109375000
2 3.8359375000 3.3359375000
3 -2.8265625000 3.8359375000
4 1.1359375000 -2.8265625000
5 -1.3765625000 1.1359375000
6 -3.7390625000 -1.3765625000
7 0.3984375000 -3.7390625000
8 -5.3991071429 0.3984375000
9 -1.3848214286 -5.3991071429
10 -0.6991071429 -1.3848214286
11 -4.4419642857 -0.6991071429
12 -0.5890625000 -4.4419642857
13 -0.4640625000 -0.5890625000
14 -1.1640625000 -0.4640625000
15 3.1734375000 -1.1640625000
16 2.8359375000 3.1734375000
17 -1.9765625000 2.8359375000
18 2.1609375000 -1.9765625000
19 0.2984375000 2.1609375000
20 0.9008928571 0.2984375000
21 4.3151785714 0.9008928571
22 0.0008928571 4.3151785714
23 -1.6419642857 0.0008928571
24 -4.7703125000 -1.6419642857
25 -6.0453125000 -4.7703125000
26 -9.2453125000 -6.0453125000
27 -4.2078125000 -9.2453125000
28 -7.4453125000 -4.2078125000
29 -9.4578125000 -7.4453125000
30 -3.5203125000 -9.4578125000
31 -8.7828125000 -3.5203125000
32 -4.2803571429 -8.7828125000
33 -4.1660714286 -4.2803571429
34 -9.5803571429 -4.1660714286
35 -3.2232142857 -9.5803571429
36 -7.5703125000 -3.2232142857
37 -6.2453125000 -7.5703125000
38 -1.0453125000 -6.2453125000
39 -2.6078125000 -1.0453125000
40 -7.1453125000 -2.6078125000
41 -1.9578125000 -7.1453125000
42 -0.8203125000 -1.9578125000
43 -1.4828125000 -0.8203125000
44 2.4196428571 -1.4828125000
45 -3.6660714286 2.4196428571
46 -3.3803571429 -3.6660714286
47 1.3767857143 -3.3803571429
48 -5.1703125000 1.3767857143
49 -5.7453125000 -5.1703125000
50 -5.9453125000 -5.7453125000
51 -4.2078125000 -5.9453125000
52 -6.6453125000 -4.2078125000
53 1.2421875000 -6.6453125000
54 -7.3203125000 1.2421875000
55 -3.2828125000 -7.3203125000
56 1.2196428571 -3.2828125000
57 -7.3660714286 1.2196428571
58 2.2196428571 -7.3660714286
59 3.3767857143 2.2196428571
60 -0.7703125000 3.3767857143
61 0.0546875000 -0.7703125000
62 4.1546875000 0.0546875000
63 -3.5078125000 4.1546875000
64 6.0546875000 -3.5078125000
65 2.8421875000 6.0546875000
66 -1.5203125000 2.8421875000
67 2.7171875000 -1.5203125000
68 3.3196428571 2.7171875000
69 3.7339285714 3.3196428571
70 5.0196428571 3.7339285714
71 2.7767857143 5.0196428571
72 5.7296875000 2.7767857143
73 3.8546875000 5.7296875000
74 7.4546875000 3.8546875000
75 2.3921875000 7.4546875000
76 5.7546875000 2.3921875000
77 5.8421875000 5.7546875000
78 5.9796875000 5.8421875000
79 6.1171875000 5.9796875000
80 1.8196428571 6.1171875000
81 8.5339285714 1.8196428571
82 6.4196428571 8.5339285714
83 1.7767857143 6.4196428571
84 9.8296875000 1.7767857143
85 11.2546875000 9.8296875000
86 1.9546875000 11.2546875000
87 11.7921875000 1.9546875000
88 5.4546875000 11.7921875000
89 4.8421875000 5.4546875000
90 8.7796875000 4.8421875000
91 4.0171875000 8.7796875000
> 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/7vzp61229040818.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/8ds1h1229040818.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/98a2p1229040818.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/10atjs1229040818.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/11tynp1229040818.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/125lhf1229040818.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/13vm6k1229040818.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/14nhtr1229040818.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/1503zr1229040818.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/16fcru1229040818.tab")
+ }
>
> system("convert tmp/1clow1229040817.ps tmp/1clow1229040817.png")
> system("convert tmp/2ywit1229040817.ps tmp/2ywit1229040817.png")
> system("convert tmp/35fkm1229040818.ps tmp/35fkm1229040818.png")
> system("convert tmp/4qirp1229040818.ps tmp/4qirp1229040818.png")
> system("convert tmp/5shwo1229040818.ps tmp/5shwo1229040818.png")
> system("convert tmp/6kevc1229040818.ps tmp/6kevc1229040818.png")
> system("convert tmp/7vzp61229040818.ps tmp/7vzp61229040818.png")
> system("convert tmp/8ds1h1229040818.ps tmp/8ds1h1229040818.png")
> system("convert tmp/98a2p1229040818.ps tmp/98a2p1229040818.png")
> system("convert tmp/10atjs1229040818.ps tmp/10atjs1229040818.png")
>
>
> proc.time()
user system elapsed
5.658 2.759 6.063