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(124,0,113,0,109,0,109,0,106,0,101,0,98,0,93,0,91,0,122,0,139,0,140,0,132,0,117,0,114,0,113,0,110,0,107,0,103,0,98,0,98,0,137,0,148,0,147,0,139,0,130,0,128,0,127,0,123,0,118,0,114,0,108,0,111,0,151,0,159,0,158,0,148,0,138,0,137,0,136,0,133,0,126,0,120,0,114,0,116,0,153,0,162,0,161,0,149,0,139,0,135,0,130,0,127,0,122,0,117,0,112,0,113,0,149,0,157,0,157,0,147,0,137,0,132,0,125,0,123,0,117,0,114,0,111,0,112,0,144,0,150,0,149,0,134,0,123,0,116,0,117,1,111,1,105,1,102,1,95,1,93,1,124,1,130,1,124,1,115,1,106,1,105,1,105,1,101,1,95,1,93,1,84,1,87,1,116,1,120,1,117,1,109,1),dim=c(2,97),dimnames=list(c('Y','X'),1:97))
> y <- array(NA,dim=c(2,97),dimnames=list(c('Y','X'),1:97))
> 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 124 0
2 113 0
3 109 0
4 109 0
5 106 0
6 101 0
7 98 0
8 93 0
9 91 0
10 122 0
11 139 0
12 140 0
13 132 0
14 117 0
15 114 0
16 113 0
17 110 0
18 107 0
19 103 0
20 98 0
21 98 0
22 137 0
23 148 0
24 147 0
25 139 0
26 130 0
27 128 0
28 127 0
29 123 0
30 118 0
31 114 0
32 108 0
33 111 0
34 151 0
35 159 0
36 158 0
37 148 0
38 138 0
39 137 0
40 136 0
41 133 0
42 126 0
43 120 0
44 114 0
45 116 0
46 153 0
47 162 0
48 161 0
49 149 0
50 139 0
51 135 0
52 130 0
53 127 0
54 122 0
55 117 0
56 112 0
57 113 0
58 149 0
59 157 0
60 157 0
61 147 0
62 137 0
63 132 0
64 125 0
65 123 0
66 117 0
67 114 0
68 111 0
69 112 0
70 144 0
71 150 0
72 149 0
73 134 0
74 123 0
75 116 0
76 117 1
77 111 1
78 105 1
79 102 1
80 95 1
81 93 1
82 124 1
83 130 1
84 124 1
85 115 1
86 106 1
87 105 1
88 105 1
89 101 1
90 95 1
91 93 1
92 84 1
93 87 1
94 116 1
95 120 1
96 117 1
97 109 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X
126.93 -19.93
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-35.933 -12.933 -1.933 12.067 35.067
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 126.933 1.961 64.739 < 2e-16 ***
X -19.933 4.117 -4.842 4.98e-06 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 16.98 on 95 degrees of freedom
Multiple R-squared: 0.1979, Adjusted R-squared: 0.1895
F-statistic: 23.44 on 1 and 95 DF, p-value: 4.983e-06
> 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.1241086 0.24821729 0.875891353
[2,] 0.1014631 0.20292616 0.898536921
[3,] 0.0963011 0.19260220 0.903698901
[4,] 0.1269552 0.25391033 0.873044835
[5,] 0.1618351 0.32367010 0.838164949
[6,] 0.1898509 0.37970180 0.810149101
[7,] 0.4767564 0.95351278 0.523243610
[8,] 0.6578591 0.68428185 0.342140923
[9,] 0.6678589 0.66428216 0.332141081
[10,] 0.5923115 0.81537702 0.407688509
[11,] 0.5168971 0.96620584 0.483102922
[12,] 0.4454145 0.89082906 0.554585468
[13,] 0.3872024 0.77440484 0.612797582
[14,] 0.3480624 0.69612488 0.651937559
[15,] 0.3398595 0.67971906 0.660140469
[16,] 0.3832530 0.76650603 0.616746987
[17,] 0.4331257 0.86625148 0.566874258
[18,] 0.5260402 0.94791968 0.473959840
[19,] 0.7289940 0.54201201 0.271006004
[20,] 0.8358982 0.32820366 0.164101830
[21,] 0.8533729 0.29325420 0.146627098
[22,] 0.8316604 0.33667925 0.168339627
[23,] 0.8020507 0.39589856 0.197949280
[24,] 0.7669502 0.46609967 0.233049834
[25,] 0.7245845 0.55083096 0.275415480
[26,] 0.6844817 0.63103653 0.315518263
[27,] 0.6565446 0.68691073 0.343455365
[28,] 0.6647204 0.67055923 0.335279613
[29,] 0.6584023 0.68319549 0.341597746
[30,] 0.7734369 0.45312615 0.226563077
[31,] 0.9040300 0.19193991 0.095969955
[32,] 0.9596458 0.08070846 0.040354231
[33,] 0.9686397 0.06272067 0.031360335
[34,] 0.9634245 0.07315110 0.036575550
[35,] 0.9561534 0.08769312 0.043846561
[36,] 0.9463912 0.10721761 0.053608804
[37,] 0.9318640 0.13627190 0.068135951
[38,] 0.9122630 0.17547395 0.087736975
[39,] 0.8947777 0.21044465 0.105222326
[40,] 0.8899713 0.22005732 0.110028660
[41,] 0.8810344 0.23793116 0.118965578
[42,] 0.9143002 0.17139951 0.085699757
[43,] 0.9657678 0.06846447 0.034232237
[44,] 0.9874672 0.02506554 0.012532769
[45,] 0.9898526 0.02029483 0.010147417
[46,] 0.9872515 0.02549708 0.012748542
[47,] 0.9825434 0.03491315 0.017456574
[48,] 0.9751463 0.04970731 0.024853656
[49,] 0.9652688 0.06946232 0.034731158
[50,] 0.9546731 0.09065375 0.045326873
[51,] 0.9478853 0.10422948 0.052114742
[52,] 0.9507840 0.09843205 0.049216026
[53,] 0.9535288 0.09294242 0.046471212
[54,] 0.9582845 0.08343095 0.041715476
[55,] 0.9780759 0.04384824 0.021924122
[56,] 0.9905425 0.01891490 0.009457452
[57,] 0.9922892 0.01542160 0.007710800
[58,] 0.9899356 0.02012875 0.010064373
[59,] 0.9853468 0.02930639 0.014653193
[60,] 0.9781760 0.04364804 0.021824021
[61,] 0.9687165 0.06256700 0.031283499
[62,] 0.9613810 0.07723805 0.038619026
[63,] 0.9590134 0.08197313 0.040986563
[64,] 0.9660969 0.06780626 0.033903128
[65,] 0.9757498 0.04850045 0.024250225
[66,] 0.9700699 0.05986014 0.029930069
[67,] 0.9750301 0.04993984 0.024969920
[68,] 0.9839610 0.03207804 0.016039018
[69,] 0.9804531 0.03909389 0.019546947
[70,] 0.9708248 0.05835049 0.029175247
[71,] 0.9566200 0.08675998 0.043379989
[72,] 0.9445121 0.11097579 0.055487897
[73,] 0.9210204 0.15795911 0.078979557
[74,] 0.8880057 0.22398859 0.111994293
[75,] 0.8485496 0.30290085 0.151450423
[76,] 0.8259457 0.34810869 0.174054343
[77,] 0.8156300 0.36874009 0.184370047
[78,] 0.8206150 0.35877003 0.179385014
[79,] 0.8888793 0.22224144 0.111120720
[80,] 0.9143458 0.17130839 0.085654193
[81,] 0.8951400 0.20971995 0.104859974
[82,] 0.8398938 0.32021239 0.160106197
[83,] 0.7637278 0.47254435 0.236272176
[84,] 0.6665199 0.66696017 0.333480086
[85,] 0.5499383 0.90012334 0.450061672
[86,] 0.4502022 0.90040440 0.549797799
[87,] 0.3702443 0.74048862 0.629755690
[88,] 0.4940994 0.98819873 0.505900634
> postscript(file="/var/www/html/rcomp/tmp/1h99z1227611906.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/2rb711227611906.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/396v91227611906.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/4fh8w1227611906.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/545lp1227611906.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 = 97
Frequency = 1
1 2 3 4 5 6
-2.93333333 -13.93333333 -17.93333333 -17.93333333 -20.93333333 -25.93333333
7 8 9 10 11 12
-28.93333333 -33.93333333 -35.93333333 -4.93333333 12.06666667 13.06666667
13 14 15 16 17 18
5.06666667 -9.93333333 -12.93333333 -13.93333333 -16.93333333 -19.93333333
19 20 21 22 23 24
-23.93333333 -28.93333333 -28.93333333 10.06666667 21.06666667 20.06666667
25 26 27 28 29 30
12.06666667 3.06666667 1.06666667 0.06666667 -3.93333333 -8.93333333
31 32 33 34 35 36
-12.93333333 -18.93333333 -15.93333333 24.06666667 32.06666667 31.06666667
37 38 39 40 41 42
21.06666667 11.06666667 10.06666667 9.06666667 6.06666667 -0.93333333
43 44 45 46 47 48
-6.93333333 -12.93333333 -10.93333333 26.06666667 35.06666667 34.06666667
49 50 51 52 53 54
22.06666667 12.06666667 8.06666667 3.06666667 0.06666667 -4.93333333
55 56 57 58 59 60
-9.93333333 -14.93333333 -13.93333333 22.06666667 30.06666667 30.06666667
61 62 63 64 65 66
20.06666667 10.06666667 5.06666667 -1.93333333 -3.93333333 -9.93333333
67 68 69 70 71 72
-12.93333333 -15.93333333 -14.93333333 17.06666667 23.06666667 22.06666667
73 74 75 76 77 78
7.06666667 -3.93333333 -10.93333333 10.00000000 4.00000000 -2.00000000
79 80 81 82 83 84
-5.00000000 -12.00000000 -14.00000000 17.00000000 23.00000000 17.00000000
85 86 87 88 89 90
8.00000000 -1.00000000 -2.00000000 -2.00000000 -6.00000000 -12.00000000
91 92 93 94 95 96
-14.00000000 -23.00000000 -20.00000000 9.00000000 13.00000000 10.00000000
97
2.00000000
> postscript(file="/var/www/html/rcomp/tmp/6wwtl1227611906.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 = 97
Frequency = 1
lag(myerror, k = 1) myerror
0 -2.93333333 NA
1 -13.93333333 -2.93333333
2 -17.93333333 -13.93333333
3 -17.93333333 -17.93333333
4 -20.93333333 -17.93333333
5 -25.93333333 -20.93333333
6 -28.93333333 -25.93333333
7 -33.93333333 -28.93333333
8 -35.93333333 -33.93333333
9 -4.93333333 -35.93333333
10 12.06666667 -4.93333333
11 13.06666667 12.06666667
12 5.06666667 13.06666667
13 -9.93333333 5.06666667
14 -12.93333333 -9.93333333
15 -13.93333333 -12.93333333
16 -16.93333333 -13.93333333
17 -19.93333333 -16.93333333
18 -23.93333333 -19.93333333
19 -28.93333333 -23.93333333
20 -28.93333333 -28.93333333
21 10.06666667 -28.93333333
22 21.06666667 10.06666667
23 20.06666667 21.06666667
24 12.06666667 20.06666667
25 3.06666667 12.06666667
26 1.06666667 3.06666667
27 0.06666667 1.06666667
28 -3.93333333 0.06666667
29 -8.93333333 -3.93333333
30 -12.93333333 -8.93333333
31 -18.93333333 -12.93333333
32 -15.93333333 -18.93333333
33 24.06666667 -15.93333333
34 32.06666667 24.06666667
35 31.06666667 32.06666667
36 21.06666667 31.06666667
37 11.06666667 21.06666667
38 10.06666667 11.06666667
39 9.06666667 10.06666667
40 6.06666667 9.06666667
41 -0.93333333 6.06666667
42 -6.93333333 -0.93333333
43 -12.93333333 -6.93333333
44 -10.93333333 -12.93333333
45 26.06666667 -10.93333333
46 35.06666667 26.06666667
47 34.06666667 35.06666667
48 22.06666667 34.06666667
49 12.06666667 22.06666667
50 8.06666667 12.06666667
51 3.06666667 8.06666667
52 0.06666667 3.06666667
53 -4.93333333 0.06666667
54 -9.93333333 -4.93333333
55 -14.93333333 -9.93333333
56 -13.93333333 -14.93333333
57 22.06666667 -13.93333333
58 30.06666667 22.06666667
59 30.06666667 30.06666667
60 20.06666667 30.06666667
61 10.06666667 20.06666667
62 5.06666667 10.06666667
63 -1.93333333 5.06666667
64 -3.93333333 -1.93333333
65 -9.93333333 -3.93333333
66 -12.93333333 -9.93333333
67 -15.93333333 -12.93333333
68 -14.93333333 -15.93333333
69 17.06666667 -14.93333333
70 23.06666667 17.06666667
71 22.06666667 23.06666667
72 7.06666667 22.06666667
73 -3.93333333 7.06666667
74 -10.93333333 -3.93333333
75 10.00000000 -10.93333333
76 4.00000000 10.00000000
77 -2.00000000 4.00000000
78 -5.00000000 -2.00000000
79 -12.00000000 -5.00000000
80 -14.00000000 -12.00000000
81 17.00000000 -14.00000000
82 23.00000000 17.00000000
83 17.00000000 23.00000000
84 8.00000000 17.00000000
85 -1.00000000 8.00000000
86 -2.00000000 -1.00000000
87 -2.00000000 -2.00000000
88 -6.00000000 -2.00000000
89 -12.00000000 -6.00000000
90 -14.00000000 -12.00000000
91 -23.00000000 -14.00000000
92 -20.00000000 -23.00000000
93 9.00000000 -20.00000000
94 13.00000000 9.00000000
95 10.00000000 13.00000000
96 2.00000000 10.00000000
97 NA 2.00000000
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -13.93333333 -2.93333333
[2,] -17.93333333 -13.93333333
[3,] -17.93333333 -17.93333333
[4,] -20.93333333 -17.93333333
[5,] -25.93333333 -20.93333333
[6,] -28.93333333 -25.93333333
[7,] -33.93333333 -28.93333333
[8,] -35.93333333 -33.93333333
[9,] -4.93333333 -35.93333333
[10,] 12.06666667 -4.93333333
[11,] 13.06666667 12.06666667
[12,] 5.06666667 13.06666667
[13,] -9.93333333 5.06666667
[14,] -12.93333333 -9.93333333
[15,] -13.93333333 -12.93333333
[16,] -16.93333333 -13.93333333
[17,] -19.93333333 -16.93333333
[18,] -23.93333333 -19.93333333
[19,] -28.93333333 -23.93333333
[20,] -28.93333333 -28.93333333
[21,] 10.06666667 -28.93333333
[22,] 21.06666667 10.06666667
[23,] 20.06666667 21.06666667
[24,] 12.06666667 20.06666667
[25,] 3.06666667 12.06666667
[26,] 1.06666667 3.06666667
[27,] 0.06666667 1.06666667
[28,] -3.93333333 0.06666667
[29,] -8.93333333 -3.93333333
[30,] -12.93333333 -8.93333333
[31,] -18.93333333 -12.93333333
[32,] -15.93333333 -18.93333333
[33,] 24.06666667 -15.93333333
[34,] 32.06666667 24.06666667
[35,] 31.06666667 32.06666667
[36,] 21.06666667 31.06666667
[37,] 11.06666667 21.06666667
[38,] 10.06666667 11.06666667
[39,] 9.06666667 10.06666667
[40,] 6.06666667 9.06666667
[41,] -0.93333333 6.06666667
[42,] -6.93333333 -0.93333333
[43,] -12.93333333 -6.93333333
[44,] -10.93333333 -12.93333333
[45,] 26.06666667 -10.93333333
[46,] 35.06666667 26.06666667
[47,] 34.06666667 35.06666667
[48,] 22.06666667 34.06666667
[49,] 12.06666667 22.06666667
[50,] 8.06666667 12.06666667
[51,] 3.06666667 8.06666667
[52,] 0.06666667 3.06666667
[53,] -4.93333333 0.06666667
[54,] -9.93333333 -4.93333333
[55,] -14.93333333 -9.93333333
[56,] -13.93333333 -14.93333333
[57,] 22.06666667 -13.93333333
[58,] 30.06666667 22.06666667
[59,] 30.06666667 30.06666667
[60,] 20.06666667 30.06666667
[61,] 10.06666667 20.06666667
[62,] 5.06666667 10.06666667
[63,] -1.93333333 5.06666667
[64,] -3.93333333 -1.93333333
[65,] -9.93333333 -3.93333333
[66,] -12.93333333 -9.93333333
[67,] -15.93333333 -12.93333333
[68,] -14.93333333 -15.93333333
[69,] 17.06666667 -14.93333333
[70,] 23.06666667 17.06666667
[71,] 22.06666667 23.06666667
[72,] 7.06666667 22.06666667
[73,] -3.93333333 7.06666667
[74,] -10.93333333 -3.93333333
[75,] 10.00000000 -10.93333333
[76,] 4.00000000 10.00000000
[77,] -2.00000000 4.00000000
[78,] -5.00000000 -2.00000000
[79,] -12.00000000 -5.00000000
[80,] -14.00000000 -12.00000000
[81,] 17.00000000 -14.00000000
[82,] 23.00000000 17.00000000
[83,] 17.00000000 23.00000000
[84,] 8.00000000 17.00000000
[85,] -1.00000000 8.00000000
[86,] -2.00000000 -1.00000000
[87,] -2.00000000 -2.00000000
[88,] -6.00000000 -2.00000000
[89,] -12.00000000 -6.00000000
[90,] -14.00000000 -12.00000000
[91,] -23.00000000 -14.00000000
[92,] -20.00000000 -23.00000000
[93,] 9.00000000 -20.00000000
[94,] 13.00000000 9.00000000
[95,] 10.00000000 13.00000000
[96,] 2.00000000 10.00000000
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -13.93333333 -2.93333333
2 -17.93333333 -13.93333333
3 -17.93333333 -17.93333333
4 -20.93333333 -17.93333333
5 -25.93333333 -20.93333333
6 -28.93333333 -25.93333333
7 -33.93333333 -28.93333333
8 -35.93333333 -33.93333333
9 -4.93333333 -35.93333333
10 12.06666667 -4.93333333
11 13.06666667 12.06666667
12 5.06666667 13.06666667
13 -9.93333333 5.06666667
14 -12.93333333 -9.93333333
15 -13.93333333 -12.93333333
16 -16.93333333 -13.93333333
17 -19.93333333 -16.93333333
18 -23.93333333 -19.93333333
19 -28.93333333 -23.93333333
20 -28.93333333 -28.93333333
21 10.06666667 -28.93333333
22 21.06666667 10.06666667
23 20.06666667 21.06666667
24 12.06666667 20.06666667
25 3.06666667 12.06666667
26 1.06666667 3.06666667
27 0.06666667 1.06666667
28 -3.93333333 0.06666667
29 -8.93333333 -3.93333333
30 -12.93333333 -8.93333333
31 -18.93333333 -12.93333333
32 -15.93333333 -18.93333333
33 24.06666667 -15.93333333
34 32.06666667 24.06666667
35 31.06666667 32.06666667
36 21.06666667 31.06666667
37 11.06666667 21.06666667
38 10.06666667 11.06666667
39 9.06666667 10.06666667
40 6.06666667 9.06666667
41 -0.93333333 6.06666667
42 -6.93333333 -0.93333333
43 -12.93333333 -6.93333333
44 -10.93333333 -12.93333333
45 26.06666667 -10.93333333
46 35.06666667 26.06666667
47 34.06666667 35.06666667
48 22.06666667 34.06666667
49 12.06666667 22.06666667
50 8.06666667 12.06666667
51 3.06666667 8.06666667
52 0.06666667 3.06666667
53 -4.93333333 0.06666667
54 -9.93333333 -4.93333333
55 -14.93333333 -9.93333333
56 -13.93333333 -14.93333333
57 22.06666667 -13.93333333
58 30.06666667 22.06666667
59 30.06666667 30.06666667
60 20.06666667 30.06666667
61 10.06666667 20.06666667
62 5.06666667 10.06666667
63 -1.93333333 5.06666667
64 -3.93333333 -1.93333333
65 -9.93333333 -3.93333333
66 -12.93333333 -9.93333333
67 -15.93333333 -12.93333333
68 -14.93333333 -15.93333333
69 17.06666667 -14.93333333
70 23.06666667 17.06666667
71 22.06666667 23.06666667
72 7.06666667 22.06666667
73 -3.93333333 7.06666667
74 -10.93333333 -3.93333333
75 10.00000000 -10.93333333
76 4.00000000 10.00000000
77 -2.00000000 4.00000000
78 -5.00000000 -2.00000000
79 -12.00000000 -5.00000000
80 -14.00000000 -12.00000000
81 17.00000000 -14.00000000
82 23.00000000 17.00000000
83 17.00000000 23.00000000
84 8.00000000 17.00000000
85 -1.00000000 8.00000000
86 -2.00000000 -1.00000000
87 -2.00000000 -2.00000000
88 -6.00000000 -2.00000000
89 -12.00000000 -6.00000000
90 -14.00000000 -12.00000000
91 -23.00000000 -14.00000000
92 -20.00000000 -23.00000000
93 9.00000000 -20.00000000
94 13.00000000 9.00000000
95 10.00000000 13.00000000
96 2.00000000 10.00000000
> 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/7kjhj1227611906.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/8ihoi1227611906.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/9zgfr1227611906.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/10mph91227611906.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/11o23l1227611906.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/12k4f41227611906.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/13niqd1227611906.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/1477gq1227611906.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/15ujhn1227611906.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/16e4lt1227611906.tab")
+ }
>
> system("convert tmp/1h99z1227611906.ps tmp/1h99z1227611906.png")
> system("convert tmp/2rb711227611906.ps tmp/2rb711227611906.png")
> system("convert tmp/396v91227611906.ps tmp/396v91227611906.png")
> system("convert tmp/4fh8w1227611906.ps tmp/4fh8w1227611906.png")
> system("convert tmp/545lp1227611906.ps tmp/545lp1227611906.png")
> system("convert tmp/6wwtl1227611906.ps tmp/6wwtl1227611906.png")
> system("convert tmp/7kjhj1227611906.ps tmp/7kjhj1227611906.png")
> system("convert tmp/8ihoi1227611906.ps tmp/8ihoi1227611906.png")
> system("convert tmp/9zgfr1227611906.ps tmp/9zgfr1227611906.png")
> system("convert tmp/10mph91227611906.ps tmp/10mph91227611906.png")
>
>
> proc.time()
user system elapsed
2.849 1.597 3.796