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(0,7.3,0,7.1,0,6.9,0,6.8,0,7.5,0,7.6,0,7.8,0,8,0,8.1,0,8.2,0,8.3,0,8.2,0,8,0,7.9,0,7.6,0,7.6,0,8.2,0,8.3,0,8.4,0,8.4,0,8.4,0,8.6,0,8.9,0,8.8,0,8.3,0,7.5,0,7.2,0,7.5,0,8.8,0,9.3,0,9.3,0,8.7,0,8.2,0,8.3,0,8.5,0,8.6,0,8.6,0,8.2,0,8.1,0,8,0,8.6,0,8.7,0,8.8,0,8.5,0,8.4,0,8.5,0,8.7,0,8.7,0,8.6,0,8.5,0,8.3,0,8.1,0,8.2,0,8.1,0,8.1,0,7.9,0,7.9,0,7.9,0,8,0,8,0,7.9,0,8,1,7.7,1,7.2,1,7.5,1,7.3,1,7,1,7,1,7,1,7.2,1,7.3,1,7.1,1,6.8,1,6.6,1,6.2,1,6.2,1,6.8,1,6.9,1,6.8),dim=c(2,79),dimnames=list(c('y','x'),1:79))
> y <- array(NA,dim=c(2,79),dimnames=list(c('y','x'),1:79))
> 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 0 7.3
2 0 7.1
3 0 6.9
4 0 6.8
5 0 7.5
6 0 7.6
7 0 7.8
8 0 8.0
9 0 8.1
10 0 8.2
11 0 8.3
12 0 8.2
13 0 8.0
14 0 7.9
15 0 7.6
16 0 7.6
17 0 8.2
18 0 8.3
19 0 8.4
20 0 8.4
21 0 8.4
22 0 8.6
23 0 8.9
24 0 8.8
25 0 8.3
26 0 7.5
27 0 7.2
28 0 7.5
29 0 8.8
30 0 9.3
31 0 9.3
32 0 8.7
33 0 8.2
34 0 8.3
35 0 8.5
36 0 8.6
37 0 8.6
38 0 8.2
39 0 8.1
40 0 8.0
41 0 8.6
42 0 8.7
43 0 8.8
44 0 8.5
45 0 8.4
46 0 8.5
47 0 8.7
48 0 8.7
49 0 8.6
50 0 8.5
51 0 8.3
52 0 8.1
53 0 8.2
54 0 8.1
55 0 8.1
56 0 7.9
57 0 7.9
58 0 7.9
59 0 8.0
60 0 8.0
61 0 7.9
62 0 8.0
63 1 7.7
64 1 7.2
65 1 7.5
66 1 7.3
67 1 7.0
68 1 7.0
69 1 7.0
70 1 7.2
71 1 7.3
72 1 7.1
73 1 6.8
74 1 6.6
75 1 6.2
76 1 6.2
77 1 6.8
78 1 6.9
79 1 6.8
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) x
3.5312 -0.4191
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.68102 -0.17805 -0.01040 0.15726 0.69621
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.53117 0.37986 9.296 3.23e-14 ***
x -0.41914 0.04783 -8.763 3.44e-13 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.2945 on 77 degrees of freedom
Multiple R-squared: 0.4993, Adjusted R-squared: 0.4928
F-statistic: 76.79 on 1 and 77 DF, p-value: 3.437e-13
> 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 0.000000e+00 1.000000e+00
[2,] 0 0.000000e+00 1.000000e+00
[3,] 0 0.000000e+00 1.000000e+00
[4,] 0 0.000000e+00 1.000000e+00
[5,] 0 0.000000e+00 1.000000e+00
[6,] 0 0.000000e+00 1.000000e+00
[7,] 0 0.000000e+00 1.000000e+00
[8,] 0 0.000000e+00 1.000000e+00
[9,] 0 0.000000e+00 1.000000e+00
[10,] 0 0.000000e+00 1.000000e+00
[11,] 0 0.000000e+00 1.000000e+00
[12,] 0 0.000000e+00 1.000000e+00
[13,] 0 0.000000e+00 1.000000e+00
[14,] 0 0.000000e+00 1.000000e+00
[15,] 0 0.000000e+00 1.000000e+00
[16,] 0 0.000000e+00 1.000000e+00
[17,] 0 0.000000e+00 1.000000e+00
[18,] 0 0.000000e+00 1.000000e+00
[19,] 0 0.000000e+00 1.000000e+00
[20,] 0 0.000000e+00 1.000000e+00
[21,] 0 0.000000e+00 1.000000e+00
[22,] 0 0.000000e+00 1.000000e+00
[23,] 0 0.000000e+00 1.000000e+00
[24,] 0 0.000000e+00 1.000000e+00
[25,] 0 0.000000e+00 1.000000e+00
[26,] 0 0.000000e+00 1.000000e+00
[27,] 0 0.000000e+00 1.000000e+00
[28,] 0 0.000000e+00 1.000000e+00
[29,] 0 0.000000e+00 1.000000e+00
[30,] 0 0.000000e+00 1.000000e+00
[31,] 0 0.000000e+00 1.000000e+00
[32,] 0 0.000000e+00 1.000000e+00
[33,] 0 0.000000e+00 1.000000e+00
[34,] 0 0.000000e+00 1.000000e+00
[35,] 0 0.000000e+00 1.000000e+00
[36,] 0 0.000000e+00 1.000000e+00
[37,] 0 0.000000e+00 1.000000e+00
[38,] 0 0.000000e+00 1.000000e+00
[39,] 0 0.000000e+00 1.000000e+00
[40,] 0 0.000000e+00 1.000000e+00
[41,] 0 0.000000e+00 1.000000e+00
[42,] 0 0.000000e+00 1.000000e+00
[43,] 0 0.000000e+00 1.000000e+00
[44,] 0 0.000000e+00 1.000000e+00
[45,] 0 0.000000e+00 1.000000e+00
[46,] 0 0.000000e+00 1.000000e+00
[47,] 0 0.000000e+00 1.000000e+00
[48,] 0 0.000000e+00 1.000000e+00
[49,] 0 0.000000e+00 1.000000e+00
[50,] 0 0.000000e+00 1.000000e+00
[51,] 0 0.000000e+00 1.000000e+00
[52,] 0 0.000000e+00 1.000000e+00
[53,] 0 0.000000e+00 1.000000e+00
[54,] 0 0.000000e+00 1.000000e+00
[55,] 0 0.000000e+00 1.000000e+00
[56,] 0 0.000000e+00 1.000000e+00
[57,] 0 0.000000e+00 1.000000e+00
[58,] 0 0.000000e+00 1.000000e+00
[59,] 1 0.000000e+00 0.000000e+00
[60,] 1 4.113792e-192 2.056896e-192
[61,] 1 0.000000e+00 0.000000e+00
[62,] 1 0.000000e+00 0.000000e+00
[63,] 1 3.862707e-146 1.931353e-146
[64,] 1 1.565069e-134 7.825343e-135
[65,] 1 0.000000e+00 0.000000e+00
[66,] 1 0.000000e+00 0.000000e+00
[67,] 1 1.279461e-90 6.397305e-91
[68,] 1 0.000000e+00 0.000000e+00
[69,] 1 2.088170e-60 1.044085e-60
[70,] 1 9.207747e-46 4.603874e-46
> postscript(file="/var/www/html/rcomp/tmp/127zj1227548463.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/257f91227548463.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/34cit1227548463.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/4q7ct1227548463.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/54dsm1227548463.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 = 79
Frequency = 1
1 2 3 4 5 6
-0.47144870 -0.55527663 -0.63910457 -0.68101854 -0.38762076 -0.34570679
7 8 9 10 11 12
-0.26187885 -0.17805091 -0.13613694 -0.09422298 -0.05230901 -0.09422298
13 14 15 16 17 18
-0.17805091 -0.21996488 -0.34570679 -0.34570679 -0.09422298 -0.05230901
19 20 21 22 23 24
-0.01039504 -0.01039504 -0.01039504 0.07343290 0.19917481 0.15726084
25 26 27 28 29 30
-0.05230901 -0.38762076 -0.51336266 -0.38762076 0.15726084 0.36683068
31 32 33 34 35 36
0.36683068 0.11534687 -0.09422298 -0.05230901 0.03151893 0.07343290
37 38 39 40 41 42
0.07343290 -0.09422298 -0.13613694 -0.17805091 0.07343290 0.11534687
43 44 45 46 47 48
0.15726084 0.03151893 -0.01039504 0.03151893 0.11534687 0.11534687
49 50 51 52 53 54
0.07343290 0.03151893 -0.05230901 -0.13613694 -0.09422298 -0.13613694
55 56 57 58 59 60
-0.13613694 -0.21996488 -0.21996488 -0.21996488 -0.17805091 -0.17805091
61 62 63 64 65 66
-0.21996488 -0.17805091 0.69620718 0.48663734 0.61237924 0.52855130
67 68 69 70 71 72
0.40280940 0.40280940 0.40280940 0.48663734 0.52855130 0.44472337
73 74 75 76 77 78
0.31898146 0.23515352 0.06749765 0.06749765 0.31898146 0.36089543
79
0.31898146
> postscript(file="/var/www/html/rcomp/tmp/6caco1227548463.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 = 79
Frequency = 1
lag(myerror, k = 1) myerror
0 -0.47144870 NA
1 -0.55527663 -0.47144870
2 -0.63910457 -0.55527663
3 -0.68101854 -0.63910457
4 -0.38762076 -0.68101854
5 -0.34570679 -0.38762076
6 -0.26187885 -0.34570679
7 -0.17805091 -0.26187885
8 -0.13613694 -0.17805091
9 -0.09422298 -0.13613694
10 -0.05230901 -0.09422298
11 -0.09422298 -0.05230901
12 -0.17805091 -0.09422298
13 -0.21996488 -0.17805091
14 -0.34570679 -0.21996488
15 -0.34570679 -0.34570679
16 -0.09422298 -0.34570679
17 -0.05230901 -0.09422298
18 -0.01039504 -0.05230901
19 -0.01039504 -0.01039504
20 -0.01039504 -0.01039504
21 0.07343290 -0.01039504
22 0.19917481 0.07343290
23 0.15726084 0.19917481
24 -0.05230901 0.15726084
25 -0.38762076 -0.05230901
26 -0.51336266 -0.38762076
27 -0.38762076 -0.51336266
28 0.15726084 -0.38762076
29 0.36683068 0.15726084
30 0.36683068 0.36683068
31 0.11534687 0.36683068
32 -0.09422298 0.11534687
33 -0.05230901 -0.09422298
34 0.03151893 -0.05230901
35 0.07343290 0.03151893
36 0.07343290 0.07343290
37 -0.09422298 0.07343290
38 -0.13613694 -0.09422298
39 -0.17805091 -0.13613694
40 0.07343290 -0.17805091
41 0.11534687 0.07343290
42 0.15726084 0.11534687
43 0.03151893 0.15726084
44 -0.01039504 0.03151893
45 0.03151893 -0.01039504
46 0.11534687 0.03151893
47 0.11534687 0.11534687
48 0.07343290 0.11534687
49 0.03151893 0.07343290
50 -0.05230901 0.03151893
51 -0.13613694 -0.05230901
52 -0.09422298 -0.13613694
53 -0.13613694 -0.09422298
54 -0.13613694 -0.13613694
55 -0.21996488 -0.13613694
56 -0.21996488 -0.21996488
57 -0.21996488 -0.21996488
58 -0.17805091 -0.21996488
59 -0.17805091 -0.17805091
60 -0.21996488 -0.17805091
61 -0.17805091 -0.21996488
62 0.69620718 -0.17805091
63 0.48663734 0.69620718
64 0.61237924 0.48663734
65 0.52855130 0.61237924
66 0.40280940 0.52855130
67 0.40280940 0.40280940
68 0.40280940 0.40280940
69 0.48663734 0.40280940
70 0.52855130 0.48663734
71 0.44472337 0.52855130
72 0.31898146 0.44472337
73 0.23515352 0.31898146
74 0.06749765 0.23515352
75 0.06749765 0.06749765
76 0.31898146 0.06749765
77 0.36089543 0.31898146
78 0.31898146 0.36089543
79 NA 0.31898146
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -0.55527663 -0.47144870
[2,] -0.63910457 -0.55527663
[3,] -0.68101854 -0.63910457
[4,] -0.38762076 -0.68101854
[5,] -0.34570679 -0.38762076
[6,] -0.26187885 -0.34570679
[7,] -0.17805091 -0.26187885
[8,] -0.13613694 -0.17805091
[9,] -0.09422298 -0.13613694
[10,] -0.05230901 -0.09422298
[11,] -0.09422298 -0.05230901
[12,] -0.17805091 -0.09422298
[13,] -0.21996488 -0.17805091
[14,] -0.34570679 -0.21996488
[15,] -0.34570679 -0.34570679
[16,] -0.09422298 -0.34570679
[17,] -0.05230901 -0.09422298
[18,] -0.01039504 -0.05230901
[19,] -0.01039504 -0.01039504
[20,] -0.01039504 -0.01039504
[21,] 0.07343290 -0.01039504
[22,] 0.19917481 0.07343290
[23,] 0.15726084 0.19917481
[24,] -0.05230901 0.15726084
[25,] -0.38762076 -0.05230901
[26,] -0.51336266 -0.38762076
[27,] -0.38762076 -0.51336266
[28,] 0.15726084 -0.38762076
[29,] 0.36683068 0.15726084
[30,] 0.36683068 0.36683068
[31,] 0.11534687 0.36683068
[32,] -0.09422298 0.11534687
[33,] -0.05230901 -0.09422298
[34,] 0.03151893 -0.05230901
[35,] 0.07343290 0.03151893
[36,] 0.07343290 0.07343290
[37,] -0.09422298 0.07343290
[38,] -0.13613694 -0.09422298
[39,] -0.17805091 -0.13613694
[40,] 0.07343290 -0.17805091
[41,] 0.11534687 0.07343290
[42,] 0.15726084 0.11534687
[43,] 0.03151893 0.15726084
[44,] -0.01039504 0.03151893
[45,] 0.03151893 -0.01039504
[46,] 0.11534687 0.03151893
[47,] 0.11534687 0.11534687
[48,] 0.07343290 0.11534687
[49,] 0.03151893 0.07343290
[50,] -0.05230901 0.03151893
[51,] -0.13613694 -0.05230901
[52,] -0.09422298 -0.13613694
[53,] -0.13613694 -0.09422298
[54,] -0.13613694 -0.13613694
[55,] -0.21996488 -0.13613694
[56,] -0.21996488 -0.21996488
[57,] -0.21996488 -0.21996488
[58,] -0.17805091 -0.21996488
[59,] -0.17805091 -0.17805091
[60,] -0.21996488 -0.17805091
[61,] -0.17805091 -0.21996488
[62,] 0.69620718 -0.17805091
[63,] 0.48663734 0.69620718
[64,] 0.61237924 0.48663734
[65,] 0.52855130 0.61237924
[66,] 0.40280940 0.52855130
[67,] 0.40280940 0.40280940
[68,] 0.40280940 0.40280940
[69,] 0.48663734 0.40280940
[70,] 0.52855130 0.48663734
[71,] 0.44472337 0.52855130
[72,] 0.31898146 0.44472337
[73,] 0.23515352 0.31898146
[74,] 0.06749765 0.23515352
[75,] 0.06749765 0.06749765
[76,] 0.31898146 0.06749765
[77,] 0.36089543 0.31898146
[78,] 0.31898146 0.36089543
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -0.55527663 -0.47144870
2 -0.63910457 -0.55527663
3 -0.68101854 -0.63910457
4 -0.38762076 -0.68101854
5 -0.34570679 -0.38762076
6 -0.26187885 -0.34570679
7 -0.17805091 -0.26187885
8 -0.13613694 -0.17805091
9 -0.09422298 -0.13613694
10 -0.05230901 -0.09422298
11 -0.09422298 -0.05230901
12 -0.17805091 -0.09422298
13 -0.21996488 -0.17805091
14 -0.34570679 -0.21996488
15 -0.34570679 -0.34570679
16 -0.09422298 -0.34570679
17 -0.05230901 -0.09422298
18 -0.01039504 -0.05230901
19 -0.01039504 -0.01039504
20 -0.01039504 -0.01039504
21 0.07343290 -0.01039504
22 0.19917481 0.07343290
23 0.15726084 0.19917481
24 -0.05230901 0.15726084
25 -0.38762076 -0.05230901
26 -0.51336266 -0.38762076
27 -0.38762076 -0.51336266
28 0.15726084 -0.38762076
29 0.36683068 0.15726084
30 0.36683068 0.36683068
31 0.11534687 0.36683068
32 -0.09422298 0.11534687
33 -0.05230901 -0.09422298
34 0.03151893 -0.05230901
35 0.07343290 0.03151893
36 0.07343290 0.07343290
37 -0.09422298 0.07343290
38 -0.13613694 -0.09422298
39 -0.17805091 -0.13613694
40 0.07343290 -0.17805091
41 0.11534687 0.07343290
42 0.15726084 0.11534687
43 0.03151893 0.15726084
44 -0.01039504 0.03151893
45 0.03151893 -0.01039504
46 0.11534687 0.03151893
47 0.11534687 0.11534687
48 0.07343290 0.11534687
49 0.03151893 0.07343290
50 -0.05230901 0.03151893
51 -0.13613694 -0.05230901
52 -0.09422298 -0.13613694
53 -0.13613694 -0.09422298
54 -0.13613694 -0.13613694
55 -0.21996488 -0.13613694
56 -0.21996488 -0.21996488
57 -0.21996488 -0.21996488
58 -0.17805091 -0.21996488
59 -0.17805091 -0.17805091
60 -0.21996488 -0.17805091
61 -0.17805091 -0.21996488
62 0.69620718 -0.17805091
63 0.48663734 0.69620718
64 0.61237924 0.48663734
65 0.52855130 0.61237924
66 0.40280940 0.52855130
67 0.40280940 0.40280940
68 0.40280940 0.40280940
69 0.48663734 0.40280940
70 0.52855130 0.48663734
71 0.44472337 0.52855130
72 0.31898146 0.44472337
73 0.23515352 0.31898146
74 0.06749765 0.23515352
75 0.06749765 0.06749765
76 0.31898146 0.06749765
77 0.36089543 0.31898146
78 0.31898146 0.36089543
> 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/7xwtp1227548463.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/8anks1227548463.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/9ntmw1227548463.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/10zc891227548463.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/11vmdj1227548463.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/12y9hr1227548463.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/13m4g81227548463.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/14k7631227548463.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/15gvh91227548463.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/16lqg61227548463.tab")
+ }
>
> system("convert tmp/127zj1227548463.ps tmp/127zj1227548463.png")
> system("convert tmp/257f91227548463.ps tmp/257f91227548463.png")
> system("convert tmp/34cit1227548463.ps tmp/34cit1227548463.png")
> system("convert tmp/4q7ct1227548463.ps tmp/4q7ct1227548463.png")
> system("convert tmp/54dsm1227548463.ps tmp/54dsm1227548463.png")
> system("convert tmp/6caco1227548463.ps tmp/6caco1227548463.png")
> system("convert tmp/7xwtp1227548463.ps tmp/7xwtp1227548463.png")
> system("convert tmp/8anks1227548463.ps tmp/8anks1227548463.png")
> system("convert tmp/9ntmw1227548463.ps tmp/9ntmw1227548463.png")
> system("convert tmp/10zc891227548463.ps tmp/10zc891227548463.png")
>
>
> proc.time()
user system elapsed
2.82 1.66 3.37