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(395.3,0,395.1,0,403.5,0,403.3,0,405.7,0,406.7,0,407.2,0,412.4,0,415.9,0,414.0,0,411.8,0,409.9,0,412.4,0,415.9,0,416.3,0,417.2,0,421.8,0,421.4,0,415.1,0,412.4,0,411.8,0,408.8,0,404.5,0,402.5,0,409.4,0,410.7,0,413.4,0,415.2,0,417.7,0,417.8,0,417.9,0,418.4,0,418.2,0,416.6,0,418.9,0,421.0,0,423.5,0,432.3,0,432.3,0,428.6,0,426.7,0,427.3,0,428.5,0,437.0,0,442.0,0,444.9,0,441.4,0,440.3,0,447.1,0,455.3,0,478.6,0,486.5,0,487.8,0,485.9,0,483.8,0,488.4,0,494.0,0,493.6,0,487.3,0,482.1,0,484.2,0,496.8,0,501.1,0,499.8,0,495.5,0,498.1,0,503.8,0,516.2,0,526.1,0,527.1,0,525.1,0,528.9,0,540.1,0,549.0,0,556.0,0,568.9,0,589.1,0,590.3,0,603.3,0,638.8,0,643.0,0,656.7,0,656.1,0,654.1,0,659.9,0,662.1,0,669.2,0,673.1,0,678.3,0,677.4,0,678.5,0,672.4,0,665.3,0,667.9,0,672.1,0,662.5,0,682.3,0,692.1,0,702.7,0,721.4,0,733.2,0,747.7,0,737.6,0,729.3,0,706.1,0,674.3,0,659.0,0,645.7,0,646.1,0,633.0,1,622.3,1,628.2,1,637.3,1,639.6,1,638.5,1,650.5,1,655.4,1),dim=c(2,117),dimnames=list(c('Y','X'),1:117))
> y <- array(NA,dim=c(2,117),dimnames=list(c('Y','X'),1:117))
> 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
Y X
1 395.3 0
2 395.1 0
3 403.5 0
4 403.3 0
5 405.7 0
6 406.7 0
7 407.2 0
8 412.4 0
9 415.9 0
10 414.0 0
11 411.8 0
12 409.9 0
13 412.4 0
14 415.9 0
15 416.3 0
16 417.2 0
17 421.8 0
18 421.4 0
19 415.1 0
20 412.4 0
21 411.8 0
22 408.8 0
23 404.5 0
24 402.5 0
25 409.4 0
26 410.7 0
27 413.4 0
28 415.2 0
29 417.7 0
30 417.8 0
31 417.9 0
32 418.4 0
33 418.2 0
34 416.6 0
35 418.9 0
36 421.0 0
37 423.5 0
38 432.3 0
39 432.3 0
40 428.6 0
41 426.7 0
42 427.3 0
43 428.5 0
44 437.0 0
45 442.0 0
46 444.9 0
47 441.4 0
48 440.3 0
49 447.1 0
50 455.3 0
51 478.6 0
52 486.5 0
53 487.8 0
54 485.9 0
55 483.8 0
56 488.4 0
57 494.0 0
58 493.6 0
59 487.3 0
60 482.1 0
61 484.2 0
62 496.8 0
63 501.1 0
64 499.8 0
65 495.5 0
66 498.1 0
67 503.8 0
68 516.2 0
69 526.1 0
70 527.1 0
71 525.1 0
72 528.9 0
73 540.1 0
74 549.0 0
75 556.0 0
76 568.9 0
77 589.1 0
78 590.3 0
79 603.3 0
80 638.8 0
81 643.0 0
82 656.7 0
83 656.1 0
84 654.1 0
85 659.9 0
86 662.1 0
87 669.2 0
88 673.1 0
89 678.3 0
90 677.4 0
91 678.5 0
92 672.4 0
93 665.3 0
94 667.9 0
95 672.1 0
96 662.5 0
97 682.3 0
98 692.1 0
99 702.7 0
100 721.4 0
101 733.2 0
102 747.7 0
103 737.6 0
104 729.3 0
105 706.1 0
106 674.3 0
107 659.0 0
108 645.7 0
109 646.1 0
110 633.0 1
111 622.3 1
112 628.2 1
113 637.3 1
114 639.6 1
115 638.5 1
116 650.5 1
117 655.4 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X
516.3 121.8
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-121.19 -98.39 -27.89 122.51 231.41
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 516.29 10.33 49.961 < 2e-16 ***
X 121.81 39.52 3.082 0.00257 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 107.9 on 115 degrees of freedom
Multiple R-squared: 0.07631, Adjusted R-squared: 0.06828
F-statistic: 9.501 on 1 and 115 DF, p-value: 0.002571
> 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,] 2.158547e-04 4.317094e-04 9.997841e-01
[2,] 1.648074e-05 3.296148e-05 9.999835e-01
[3,] 1.215917e-06 2.431834e-06 9.999988e-01
[4,] 2.136375e-07 4.272751e-07 9.999998e-01
[5,] 5.318057e-08 1.063611e-07 9.999999e-01
[6,] 6.898226e-09 1.379645e-08 1.000000e+00
[7,] 6.271078e-10 1.254216e-09 1.000000e+00
[8,] 4.681690e-11 9.363380e-11 1.000000e+00
[9,] 4.184377e-12 8.368753e-12 1.000000e+00
[10,] 5.754033e-13 1.150807e-12 1.000000e+00
[11,] 7.697728e-14 1.539546e-13 1.000000e+00
[12,] 1.094892e-14 2.189783e-14 1.000000e+00
[13,] 3.228724e-15 6.457449e-15 1.000000e+00
[14,] 7.189178e-16 1.437836e-15 1.000000e+00
[15,] 6.880145e-17 1.376029e-16 1.000000e+00
[16,] 5.677169e-18 1.135434e-17 1.000000e+00
[17,] 4.572141e-19 9.144283e-19 1.000000e+00
[18,] 3.764497e-20 7.528994e-20 1.000000e+00
[19,] 4.299665e-21 8.599329e-21 1.000000e+00
[20,] 6.321720e-22 1.264344e-21 1.000000e+00
[21,] 5.256743e-23 1.051349e-22 1.000000e+00
[22,] 4.369551e-24 8.739103e-24 1.000000e+00
[23,] 4.046856e-25 8.093712e-25 1.000000e+00
[24,] 4.386524e-26 8.773048e-26 1.000000e+00
[25,] 6.541177e-27 1.308235e-26 1.000000e+00
[26,] 9.802040e-28 1.960408e-27 1.000000e+00
[27,] 1.483394e-28 2.966789e-28 1.000000e+00
[28,] 2.413120e-29 4.826241e-29 1.000000e+00
[29,] 3.823441e-30 7.646881e-30 1.000000e+00
[30,] 5.064102e-31 1.012820e-30 1.000000e+00
[31,] 9.145619e-32 1.829124e-31 1.000000e+00
[32,] 2.359338e-32 4.718677e-32 1.000000e+00
[33,] 9.882926e-33 1.976585e-32 1.000000e+00
[34,] 4.122433e-32 8.244865e-32 1.000000e+00
[35,] 9.991677e-32 1.998335e-31 1.000000e+00
[36,] 7.643776e-32 1.528755e-31 1.000000e+00
[37,] 3.922008e-32 7.844016e-32 1.000000e+00
[38,] 2.212295e-32 4.424590e-32 1.000000e+00
[39,] 1.518639e-32 3.037277e-32 1.000000e+00
[40,] 5.463348e-32 1.092670e-31 1.000000e+00
[41,] 4.562870e-31 9.125740e-31 1.000000e+00
[42,] 4.517548e-30 9.035096e-30 1.000000e+00
[43,] 1.512014e-29 3.024028e-29 1.000000e+00
[44,] 3.741066e-29 7.482131e-29 1.000000e+00
[45,] 2.480052e-28 4.960103e-28 1.000000e+00
[46,] 5.156220e-27 1.031244e-26 1.000000e+00
[47,] 5.018383e-24 1.003677e-23 1.000000e+00
[48,] 2.400162e-21 4.800325e-21 1.000000e+00
[49,] 2.516796e-19 5.033593e-19 1.000000e+00
[50,] 7.686738e-18 1.537348e-17 1.000000e+00
[51,] 1.110318e-16 2.220635e-16 1.000000e+00
[52,] 1.596138e-15 3.192276e-15 1.000000e+00
[53,] 2.371694e-14 4.743389e-14 1.000000e+00
[54,] 2.462239e-13 4.924479e-13 1.000000e+00
[55,] 1.569525e-12 3.139051e-12 1.000000e+00
[56,] 8.083971e-12 1.616794e-11 1.000000e+00
[57,] 4.570555e-11 9.141111e-11 1.000000e+00
[58,] 3.689433e-10 7.378867e-10 1.000000e+00
[59,] 3.097765e-09 6.195529e-09 1.000000e+00
[60,] 2.355306e-08 4.710612e-08 1.000000e+00
[61,] 1.685434e-07 3.370867e-07 9.999998e-01
[62,] 1.303061e-06 2.606123e-06 9.999987e-01
[63,] 1.081609e-05 2.163219e-05 9.999892e-01
[64,] 9.241412e-05 1.848282e-04 9.999076e-01
[65,] 7.159135e-04 1.431827e-03 9.992841e-01
[66,] 4.560524e-03 9.121047e-03 9.954395e-01
[67,] 2.431358e-02 4.862716e-02 9.756864e-01
[68,] 1.030091e-01 2.060183e-01 8.969909e-01
[69,] 3.056272e-01 6.112544e-01 6.943728e-01
[70,] 6.154870e-01 7.690260e-01 3.845130e-01
[71,] 8.745635e-01 2.508731e-01 1.254365e-01
[72,] 9.773306e-01 4.533872e-02 2.266936e-02
[73,] 9.970279e-01 5.944271e-03 2.972135e-03
[74,] 9.997919e-01 4.161005e-04 2.080503e-04
[75,] 9.999881e-01 2.373332e-05 1.186666e-05
[76,] 9.999982e-01 3.698071e-06 1.849035e-06
[77,] 9.999996e-01 7.541188e-07 3.770594e-07
[78,] 9.999999e-01 2.455711e-07 1.227856e-07
[79,] 1.000000e+00 9.990526e-08 4.995263e-08
[80,] 1.000000e+00 4.636854e-08 2.318427e-08
[81,] 1.000000e+00 2.816925e-08 1.408462e-08
[82,] 1.000000e+00 2.048201e-08 1.024101e-08
[83,] 1.000000e+00 1.901898e-08 9.509488e-09
[84,] 1.000000e+00 2.092123e-08 1.046062e-08
[85,] 1.000000e+00 2.679882e-08 1.339941e-08
[86,] 1.000000e+00 3.788249e-08 1.894125e-08
[87,] 1.000000e+00 5.877517e-08 2.938758e-08
[88,] 1.000000e+00 9.165107e-08 4.582554e-08
[89,] 9.999999e-01 1.294759e-07 6.473796e-08
[90,] 9.999999e-01 1.973767e-07 9.868835e-08
[91,] 9.999998e-01 3.346492e-07 1.673246e-07
[92,] 9.999998e-01 4.314270e-07 2.157135e-07
[93,] 9.999996e-01 8.918947e-07 4.459474e-07
[94,] 9.999990e-01 2.078836e-06 1.039418e-06
[95,] 9.999975e-01 4.974064e-06 2.487032e-06
[96,] 9.999954e-01 9.249920e-06 4.624960e-06
[97,] 9.999942e-01 1.152752e-05 5.763759e-06
[98,] 9.999976e-01 4.855768e-06 2.427884e-06
[99,] 9.999992e-01 1.623528e-06 8.117638e-07
[100,] 9.999999e-01 1.606265e-07 8.031324e-08
[101,] 1.000000e+00 1.099450e-08 5.497250e-09
[102,] 1.000000e+00 1.801165e-08 9.005826e-09
[103,] 9.999999e-01 1.129894e-07 5.649472e-08
[104,] 9.999994e-01 1.221898e-06 6.109491e-07
[105,] 9.999937e-01 1.253874e-05 6.269369e-06
[106,] 9.999430e-01 1.139467e-04 5.697336e-05
[107,] 9.998069e-01 3.862517e-04 1.931259e-04
[108,] 9.991379e-01 1.724229e-03 8.621143e-04
> postscript(file="/var/www/html/rcomp/tmp/1y3fd1258467749.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/2pqwk1258467749.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/3ikwp1258467749.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/4m4ow1258467749.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/5vgb71258467749.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 = 117
Frequency = 1
1 2 3 4 5
-120.98990826 -121.18990826 -112.78990826 -112.98990826 -110.58990826
6 7 8 9 10
-109.58990826 -109.08990826 -103.88990826 -100.38990826 -102.28990826
11 12 13 14 15
-104.48990826 -106.38990826 -103.88990826 -100.38990826 -99.98990826
16 17 18 19 20
-99.08990826 -94.48990826 -94.88990826 -101.18990826 -103.88990826
21 22 23 24 25
-104.48990826 -107.48990826 -111.78990826 -113.78990826 -106.88990826
26 27 28 29 30
-105.58990826 -102.88990826 -101.08990826 -98.58990826 -98.48990826
31 32 33 34 35
-98.38990826 -97.88990826 -98.08990826 -99.68990826 -97.38990826
36 37 38 39 40
-95.28990826 -92.78990826 -83.98990826 -83.98990826 -87.68990826
41 42 43 44 45
-89.58990826 -88.98990826 -87.78990826 -79.28990826 -74.28990826
46 47 48 49 50
-71.38990826 -74.88990826 -75.98990826 -69.18990826 -60.98990826
51 52 53 54 55
-37.68990826 -29.78990826 -28.48990826 -30.38990826 -32.48990826
56 57 58 59 60
-27.88990826 -22.28990826 -22.68990826 -28.98990826 -34.18990826
61 62 63 64 65
-32.08990826 -19.48990826 -15.18990826 -16.48990826 -20.78990826
66 67 68 69 70
-18.18990826 -12.48990826 -0.08990826 9.81009174 10.81009174
71 72 73 74 75
8.81009174 12.61009174 23.81009174 32.71009174 39.71009174
76 77 78 79 80
52.61009174 72.81009174 74.01009174 87.01009174 122.51009174
81 82 83 84 85
126.71009174 140.41009174 139.81009174 137.81009174 143.61009174
86 87 88 89 90
145.81009174 152.91009174 156.81009174 162.01009174 161.11009174
91 92 93 94 95
162.21009174 156.11009174 149.01009174 151.61009174 155.81009174
96 97 98 99 100
146.21009174 166.01009174 175.81009174 186.41009174 205.11009174
101 102 103 104 105
216.91009174 231.41009174 221.31009174 213.01009174 189.81009174
106 107 108 109 110
158.01009174 142.71009174 129.41009174 129.81009174 -5.10000000
111 112 113 114 115
-15.80000000 -9.90000000 -0.80000000 1.50000000 0.40000000
116 117
12.40000000 17.30000000
> postscript(file="/var/www/html/rcomp/tmp/6se6p1258467749.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 = 117
Frequency = 1
lag(myerror, k = 1) myerror
0 -120.98990826 NA
1 -121.18990826 -120.98990826
2 -112.78990826 -121.18990826
3 -112.98990826 -112.78990826
4 -110.58990826 -112.98990826
5 -109.58990826 -110.58990826
6 -109.08990826 -109.58990826
7 -103.88990826 -109.08990826
8 -100.38990826 -103.88990826
9 -102.28990826 -100.38990826
10 -104.48990826 -102.28990826
11 -106.38990826 -104.48990826
12 -103.88990826 -106.38990826
13 -100.38990826 -103.88990826
14 -99.98990826 -100.38990826
15 -99.08990826 -99.98990826
16 -94.48990826 -99.08990826
17 -94.88990826 -94.48990826
18 -101.18990826 -94.88990826
19 -103.88990826 -101.18990826
20 -104.48990826 -103.88990826
21 -107.48990826 -104.48990826
22 -111.78990826 -107.48990826
23 -113.78990826 -111.78990826
24 -106.88990826 -113.78990826
25 -105.58990826 -106.88990826
26 -102.88990826 -105.58990826
27 -101.08990826 -102.88990826
28 -98.58990826 -101.08990826
29 -98.48990826 -98.58990826
30 -98.38990826 -98.48990826
31 -97.88990826 -98.38990826
32 -98.08990826 -97.88990826
33 -99.68990826 -98.08990826
34 -97.38990826 -99.68990826
35 -95.28990826 -97.38990826
36 -92.78990826 -95.28990826
37 -83.98990826 -92.78990826
38 -83.98990826 -83.98990826
39 -87.68990826 -83.98990826
40 -89.58990826 -87.68990826
41 -88.98990826 -89.58990826
42 -87.78990826 -88.98990826
43 -79.28990826 -87.78990826
44 -74.28990826 -79.28990826
45 -71.38990826 -74.28990826
46 -74.88990826 -71.38990826
47 -75.98990826 -74.88990826
48 -69.18990826 -75.98990826
49 -60.98990826 -69.18990826
50 -37.68990826 -60.98990826
51 -29.78990826 -37.68990826
52 -28.48990826 -29.78990826
53 -30.38990826 -28.48990826
54 -32.48990826 -30.38990826
55 -27.88990826 -32.48990826
56 -22.28990826 -27.88990826
57 -22.68990826 -22.28990826
58 -28.98990826 -22.68990826
59 -34.18990826 -28.98990826
60 -32.08990826 -34.18990826
61 -19.48990826 -32.08990826
62 -15.18990826 -19.48990826
63 -16.48990826 -15.18990826
64 -20.78990826 -16.48990826
65 -18.18990826 -20.78990826
66 -12.48990826 -18.18990826
67 -0.08990826 -12.48990826
68 9.81009174 -0.08990826
69 10.81009174 9.81009174
70 8.81009174 10.81009174
71 12.61009174 8.81009174
72 23.81009174 12.61009174
73 32.71009174 23.81009174
74 39.71009174 32.71009174
75 52.61009174 39.71009174
76 72.81009174 52.61009174
77 74.01009174 72.81009174
78 87.01009174 74.01009174
79 122.51009174 87.01009174
80 126.71009174 122.51009174
81 140.41009174 126.71009174
82 139.81009174 140.41009174
83 137.81009174 139.81009174
84 143.61009174 137.81009174
85 145.81009174 143.61009174
86 152.91009174 145.81009174
87 156.81009174 152.91009174
88 162.01009174 156.81009174
89 161.11009174 162.01009174
90 162.21009174 161.11009174
91 156.11009174 162.21009174
92 149.01009174 156.11009174
93 151.61009174 149.01009174
94 155.81009174 151.61009174
95 146.21009174 155.81009174
96 166.01009174 146.21009174
97 175.81009174 166.01009174
98 186.41009174 175.81009174
99 205.11009174 186.41009174
100 216.91009174 205.11009174
101 231.41009174 216.91009174
102 221.31009174 231.41009174
103 213.01009174 221.31009174
104 189.81009174 213.01009174
105 158.01009174 189.81009174
106 142.71009174 158.01009174
107 129.41009174 142.71009174
108 129.81009174 129.41009174
109 -5.10000000 129.81009174
110 -15.80000000 -5.10000000
111 -9.90000000 -15.80000000
112 -0.80000000 -9.90000000
113 1.50000000 -0.80000000
114 0.40000000 1.50000000
115 12.40000000 0.40000000
116 17.30000000 12.40000000
117 NA 17.30000000
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -121.18990826 -120.98990826
[2,] -112.78990826 -121.18990826
[3,] -112.98990826 -112.78990826
[4,] -110.58990826 -112.98990826
[5,] -109.58990826 -110.58990826
[6,] -109.08990826 -109.58990826
[7,] -103.88990826 -109.08990826
[8,] -100.38990826 -103.88990826
[9,] -102.28990826 -100.38990826
[10,] -104.48990826 -102.28990826
[11,] -106.38990826 -104.48990826
[12,] -103.88990826 -106.38990826
[13,] -100.38990826 -103.88990826
[14,] -99.98990826 -100.38990826
[15,] -99.08990826 -99.98990826
[16,] -94.48990826 -99.08990826
[17,] -94.88990826 -94.48990826
[18,] -101.18990826 -94.88990826
[19,] -103.88990826 -101.18990826
[20,] -104.48990826 -103.88990826
[21,] -107.48990826 -104.48990826
[22,] -111.78990826 -107.48990826
[23,] -113.78990826 -111.78990826
[24,] -106.88990826 -113.78990826
[25,] -105.58990826 -106.88990826
[26,] -102.88990826 -105.58990826
[27,] -101.08990826 -102.88990826
[28,] -98.58990826 -101.08990826
[29,] -98.48990826 -98.58990826
[30,] -98.38990826 -98.48990826
[31,] -97.88990826 -98.38990826
[32,] -98.08990826 -97.88990826
[33,] -99.68990826 -98.08990826
[34,] -97.38990826 -99.68990826
[35,] -95.28990826 -97.38990826
[36,] -92.78990826 -95.28990826
[37,] -83.98990826 -92.78990826
[38,] -83.98990826 -83.98990826
[39,] -87.68990826 -83.98990826
[40,] -89.58990826 -87.68990826
[41,] -88.98990826 -89.58990826
[42,] -87.78990826 -88.98990826
[43,] -79.28990826 -87.78990826
[44,] -74.28990826 -79.28990826
[45,] -71.38990826 -74.28990826
[46,] -74.88990826 -71.38990826
[47,] -75.98990826 -74.88990826
[48,] -69.18990826 -75.98990826
[49,] -60.98990826 -69.18990826
[50,] -37.68990826 -60.98990826
[51,] -29.78990826 -37.68990826
[52,] -28.48990826 -29.78990826
[53,] -30.38990826 -28.48990826
[54,] -32.48990826 -30.38990826
[55,] -27.88990826 -32.48990826
[56,] -22.28990826 -27.88990826
[57,] -22.68990826 -22.28990826
[58,] -28.98990826 -22.68990826
[59,] -34.18990826 -28.98990826
[60,] -32.08990826 -34.18990826
[61,] -19.48990826 -32.08990826
[62,] -15.18990826 -19.48990826
[63,] -16.48990826 -15.18990826
[64,] -20.78990826 -16.48990826
[65,] -18.18990826 -20.78990826
[66,] -12.48990826 -18.18990826
[67,] -0.08990826 -12.48990826
[68,] 9.81009174 -0.08990826
[69,] 10.81009174 9.81009174
[70,] 8.81009174 10.81009174
[71,] 12.61009174 8.81009174
[72,] 23.81009174 12.61009174
[73,] 32.71009174 23.81009174
[74,] 39.71009174 32.71009174
[75,] 52.61009174 39.71009174
[76,] 72.81009174 52.61009174
[77,] 74.01009174 72.81009174
[78,] 87.01009174 74.01009174
[79,] 122.51009174 87.01009174
[80,] 126.71009174 122.51009174
[81,] 140.41009174 126.71009174
[82,] 139.81009174 140.41009174
[83,] 137.81009174 139.81009174
[84,] 143.61009174 137.81009174
[85,] 145.81009174 143.61009174
[86,] 152.91009174 145.81009174
[87,] 156.81009174 152.91009174
[88,] 162.01009174 156.81009174
[89,] 161.11009174 162.01009174
[90,] 162.21009174 161.11009174
[91,] 156.11009174 162.21009174
[92,] 149.01009174 156.11009174
[93,] 151.61009174 149.01009174
[94,] 155.81009174 151.61009174
[95,] 146.21009174 155.81009174
[96,] 166.01009174 146.21009174
[97,] 175.81009174 166.01009174
[98,] 186.41009174 175.81009174
[99,] 205.11009174 186.41009174
[100,] 216.91009174 205.11009174
[101,] 231.41009174 216.91009174
[102,] 221.31009174 231.41009174
[103,] 213.01009174 221.31009174
[104,] 189.81009174 213.01009174
[105,] 158.01009174 189.81009174
[106,] 142.71009174 158.01009174
[107,] 129.41009174 142.71009174
[108,] 129.81009174 129.41009174
[109,] -5.10000000 129.81009174
[110,] -15.80000000 -5.10000000
[111,] -9.90000000 -15.80000000
[112,] -0.80000000 -9.90000000
[113,] 1.50000000 -0.80000000
[114,] 0.40000000 1.50000000
[115,] 12.40000000 0.40000000
[116,] 17.30000000 12.40000000
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -121.18990826 -120.98990826
2 -112.78990826 -121.18990826
3 -112.98990826 -112.78990826
4 -110.58990826 -112.98990826
5 -109.58990826 -110.58990826
6 -109.08990826 -109.58990826
7 -103.88990826 -109.08990826
8 -100.38990826 -103.88990826
9 -102.28990826 -100.38990826
10 -104.48990826 -102.28990826
11 -106.38990826 -104.48990826
12 -103.88990826 -106.38990826
13 -100.38990826 -103.88990826
14 -99.98990826 -100.38990826
15 -99.08990826 -99.98990826
16 -94.48990826 -99.08990826
17 -94.88990826 -94.48990826
18 -101.18990826 -94.88990826
19 -103.88990826 -101.18990826
20 -104.48990826 -103.88990826
21 -107.48990826 -104.48990826
22 -111.78990826 -107.48990826
23 -113.78990826 -111.78990826
24 -106.88990826 -113.78990826
25 -105.58990826 -106.88990826
26 -102.88990826 -105.58990826
27 -101.08990826 -102.88990826
28 -98.58990826 -101.08990826
29 -98.48990826 -98.58990826
30 -98.38990826 -98.48990826
31 -97.88990826 -98.38990826
32 -98.08990826 -97.88990826
33 -99.68990826 -98.08990826
34 -97.38990826 -99.68990826
35 -95.28990826 -97.38990826
36 -92.78990826 -95.28990826
37 -83.98990826 -92.78990826
38 -83.98990826 -83.98990826
39 -87.68990826 -83.98990826
40 -89.58990826 -87.68990826
41 -88.98990826 -89.58990826
42 -87.78990826 -88.98990826
43 -79.28990826 -87.78990826
44 -74.28990826 -79.28990826
45 -71.38990826 -74.28990826
46 -74.88990826 -71.38990826
47 -75.98990826 -74.88990826
48 -69.18990826 -75.98990826
49 -60.98990826 -69.18990826
50 -37.68990826 -60.98990826
51 -29.78990826 -37.68990826
52 -28.48990826 -29.78990826
53 -30.38990826 -28.48990826
54 -32.48990826 -30.38990826
55 -27.88990826 -32.48990826
56 -22.28990826 -27.88990826
57 -22.68990826 -22.28990826
58 -28.98990826 -22.68990826
59 -34.18990826 -28.98990826
60 -32.08990826 -34.18990826
61 -19.48990826 -32.08990826
62 -15.18990826 -19.48990826
63 -16.48990826 -15.18990826
64 -20.78990826 -16.48990826
65 -18.18990826 -20.78990826
66 -12.48990826 -18.18990826
67 -0.08990826 -12.48990826
68 9.81009174 -0.08990826
69 10.81009174 9.81009174
70 8.81009174 10.81009174
71 12.61009174 8.81009174
72 23.81009174 12.61009174
73 32.71009174 23.81009174
74 39.71009174 32.71009174
75 52.61009174 39.71009174
76 72.81009174 52.61009174
77 74.01009174 72.81009174
78 87.01009174 74.01009174
79 122.51009174 87.01009174
80 126.71009174 122.51009174
81 140.41009174 126.71009174
82 139.81009174 140.41009174
83 137.81009174 139.81009174
84 143.61009174 137.81009174
85 145.81009174 143.61009174
86 152.91009174 145.81009174
87 156.81009174 152.91009174
88 162.01009174 156.81009174
89 161.11009174 162.01009174
90 162.21009174 161.11009174
91 156.11009174 162.21009174
92 149.01009174 156.11009174
93 151.61009174 149.01009174
94 155.81009174 151.61009174
95 146.21009174 155.81009174
96 166.01009174 146.21009174
97 175.81009174 166.01009174
98 186.41009174 175.81009174
99 205.11009174 186.41009174
100 216.91009174 205.11009174
101 231.41009174 216.91009174
102 221.31009174 231.41009174
103 213.01009174 221.31009174
104 189.81009174 213.01009174
105 158.01009174 189.81009174
106 142.71009174 158.01009174
107 129.41009174 142.71009174
108 129.81009174 129.41009174
109 -5.10000000 129.81009174
110 -15.80000000 -5.10000000
111 -9.90000000 -15.80000000
112 -0.80000000 -9.90000000
113 1.50000000 -0.80000000
114 0.40000000 1.50000000
115 12.40000000 0.40000000
116 17.30000000 12.40000000
> 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/7egek1258467749.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/839k81258467749.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/9o9t61258467749.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/10tnxt1258467749.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/1100m81258467749.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/1236h11258467749.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/13jtd31258467749.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/14s3wc1258467749.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/15tt811258467749.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/16gp2o1258467749.tab")
+ }
>
> system("convert tmp/1y3fd1258467749.ps tmp/1y3fd1258467749.png")
> system("convert tmp/2pqwk1258467749.ps tmp/2pqwk1258467749.png")
> system("convert tmp/3ikwp1258467749.ps tmp/3ikwp1258467749.png")
> system("convert tmp/4m4ow1258467749.ps tmp/4m4ow1258467749.png")
> system("convert tmp/5vgb71258467749.ps tmp/5vgb71258467749.png")
> system("convert tmp/6se6p1258467749.ps tmp/6se6p1258467749.png")
> system("convert tmp/7egek1258467749.ps tmp/7egek1258467749.png")
> system("convert tmp/839k81258467749.ps tmp/839k81258467749.png")
> system("convert tmp/9o9t61258467749.ps tmp/9o9t61258467749.png")
> system("convert tmp/10tnxt1258467749.ps tmp/10tnxt1258467749.png")
>
>
> proc.time()
user system elapsed
3.169 1.673 4.919