R version 2.15.2 (2012-10-26) -- "Trick or Treat"
Copyright (C) 2012 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: i686-pc-linux-gnu (32-bit)
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(6.80
+ ,225.00
+ ,0.44
+ ,0.67
+ ,9.20
+ ,6.30
+ ,180.00
+ ,0.44
+ ,0.80
+ ,11.70
+ ,6.40
+ ,190.00
+ ,0.46
+ ,0.76
+ ,15.80
+ ,6.20
+ ,180.00
+ ,0.42
+ ,0.65
+ ,8.60
+ ,6.90
+ ,205.00
+ ,0.45
+ ,0.90
+ ,23.20
+ ,6.40
+ ,225.00
+ ,0.43
+ ,0.78
+ ,27.40
+ ,6.30
+ ,185.00
+ ,0.49
+ ,0.77
+ ,9.30
+ ,6.80
+ ,235.00
+ ,0.47
+ ,0.75
+ ,16.00
+ ,6.90
+ ,235.00
+ ,0.44
+ ,0.82
+ ,4.70
+ ,6.70
+ ,210.00
+ ,0.48
+ ,0.83
+ ,12.50
+ ,6.90
+ ,245.00
+ ,0.52
+ ,0.63
+ ,20.10
+ ,6.90
+ ,245.00
+ ,0.49
+ ,0.76
+ ,9.10
+ ,6.30
+ ,185.00
+ ,0.37
+ ,0.71
+ ,8.10
+ ,6.10
+ ,185.00
+ ,0.42
+ ,0.78
+ ,8.60
+ ,6.20
+ ,180.00
+ ,0.44
+ ,0.78
+ ,20.30
+ ,6.80
+ ,220.00
+ ,0.50
+ ,0.88
+ ,25.00
+ ,6.50
+ ,194.00
+ ,0.50
+ ,0.83
+ ,19.20
+ ,7.60
+ ,225.00
+ ,0.43
+ ,0.57
+ ,3.30
+ ,6.30
+ ,210.00
+ ,0.37
+ ,0.82
+ ,11.20
+ ,7.10
+ ,240.00
+ ,0.50
+ ,0.71
+ ,10.50
+ ,6.80
+ ,225.00
+ ,0.40
+ ,0.77
+ ,10.10
+ ,7.30
+ ,263.00
+ ,0.48
+ ,0.66
+ ,7.20
+ ,6.40
+ ,210.00
+ ,0.48
+ ,0.24
+ ,13.60
+ ,6.80
+ ,235.00
+ ,0.43
+ ,0.73
+ ,9.00
+ ,7.20
+ ,230.00
+ ,0.56
+ ,0.72
+ ,24.60
+ ,6.40
+ ,190.00
+ ,0.44
+ ,0.76
+ ,12.60
+ ,6.60
+ ,220.00
+ ,0.49
+ ,0.75
+ ,5.60
+ ,6.80
+ ,210.00
+ ,0.40
+ ,0.74
+ ,8.70
+ ,6.10
+ ,180.00
+ ,0.42
+ ,0.71
+ ,7.70
+ ,6.50
+ ,235.00
+ ,0.49
+ ,0.74
+ ,24.10
+ ,6.40
+ ,185.00
+ ,0.48
+ ,0.86
+ ,11.70
+ ,6.00
+ ,175.00
+ ,0.39
+ ,0.72
+ ,7.70
+ ,6.00
+ ,192.00
+ ,0.44
+ ,0.79
+ ,9.60
+ ,7.30
+ ,263.00
+ ,0.48
+ ,0.66
+ ,7.20
+ ,6.10
+ ,180.00
+ ,0.34
+ ,0.82
+ ,12.30
+ ,6.70
+ ,240.00
+ ,0.52
+ ,0.73
+ ,8.90
+ ,6.40
+ ,210.00
+ ,0.48
+ ,0.85
+ ,13.60
+ ,5.80
+ ,160.00
+ ,0.41
+ ,0.81
+ ,11.20
+ ,6.90
+ ,230.00
+ ,0.41
+ ,0.60
+ ,2.80
+ ,7.00
+ ,245.00
+ ,0.41
+ ,0.57
+ ,3.20
+ ,7.30
+ ,228.00
+ ,0.45
+ ,0.73
+ ,9.40
+ ,5.90
+ ,155.00
+ ,0.29
+ ,0.71
+ ,11.90
+ ,6.20
+ ,200.00
+ ,0.45
+ ,0.80
+ ,15.40
+ ,6.80
+ ,235.00
+ ,0.55
+ ,0.78
+ ,7.40
+ ,7.00
+ ,235.00
+ ,0.48
+ ,0.74
+ ,18.90
+ ,5.90
+ ,105.00
+ ,0.36
+ ,0.84
+ ,7.90
+ ,6.10
+ ,180.00
+ ,0.53
+ ,0.79
+ ,12.20
+ ,5.70
+ ,185.00
+ ,0.35
+ ,0.70
+ ,11.00
+ ,7.10
+ ,245.00
+ ,0.41
+ ,0.78
+ ,2.80
+ ,5.80
+ ,180.00
+ ,0.43
+ ,0.87
+ ,11.80
+ ,7.40
+ ,240.00
+ ,0.60
+ ,0.71
+ ,17.10
+ ,6.80
+ ,225.00
+ ,0.48
+ ,0.70
+ ,11.60
+ ,6.80
+ ,215.00
+ ,0.46
+ ,0.73
+ ,5.80
+ ,7.00
+ ,230.00
+ ,0.44
+ ,0.76
+ ,8.30)
+ ,dim=c(5
+ ,54)
+ ,dimnames=list(c('X1'
+ ,'X2'
+ ,'X3'
+ ,'X4'
+ ,'X5
')
+ ,1:54))
> y <- array(NA,dim=c(5,54),dimnames=list(c('X1','X2','X3','X4','X5
'),1:54))
> 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 = '3'
> par3 <- 'No Linear Trend'
> par2 <- 'Do not include Seasonal Dummies'
> par1 <- '3'
> #'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, 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
X3 X1 X2 X4 X5\r
1 0.44 6.8 225 0.67 9.2
2 0.44 6.3 180 0.80 11.7
3 0.46 6.4 190 0.76 15.8
4 0.42 6.2 180 0.65 8.6
5 0.45 6.9 205 0.90 23.2
6 0.43 6.4 225 0.78 27.4
7 0.49 6.3 185 0.77 9.3
8 0.47 6.8 235 0.75 16.0
9 0.44 6.9 235 0.82 4.7
10 0.48 6.7 210 0.83 12.5
11 0.52 6.9 245 0.63 20.1
12 0.49 6.9 245 0.76 9.1
13 0.37 6.3 185 0.71 8.1
14 0.42 6.1 185 0.78 8.6
15 0.44 6.2 180 0.78 20.3
16 0.50 6.8 220 0.88 25.0
17 0.50 6.5 194 0.83 19.2
18 0.43 7.6 225 0.57 3.3
19 0.37 6.3 210 0.82 11.2
20 0.50 7.1 240 0.71 10.5
21 0.40 6.8 225 0.77 10.1
22 0.48 7.3 263 0.66 7.2
23 0.48 6.4 210 0.24 13.6
24 0.43 6.8 235 0.73 9.0
25 0.56 7.2 230 0.72 24.6
26 0.44 6.4 190 0.76 12.6
27 0.49 6.6 220 0.75 5.6
28 0.40 6.8 210 0.74 8.7
29 0.42 6.1 180 0.71 7.7
30 0.49 6.5 235 0.74 24.1
31 0.48 6.4 185 0.86 11.7
32 0.39 6.0 175 0.72 7.7
33 0.44 6.0 192 0.79 9.6
34 0.48 7.3 263 0.66 7.2
35 0.34 6.1 180 0.82 12.3
36 0.52 6.7 240 0.73 8.9
37 0.48 6.4 210 0.85 13.6
38 0.41 5.8 160 0.81 11.2
39 0.41 6.9 230 0.60 2.8
40 0.41 7.0 245 0.57 3.2
41 0.45 7.3 228 0.73 9.4
42 0.29 5.9 155 0.71 11.9
43 0.45 6.2 200 0.80 15.4
44 0.55 6.8 235 0.78 7.4
45 0.48 7.0 235 0.74 18.9
46 0.36 5.9 105 0.84 7.9
47 0.53 6.1 180 0.79 12.2
48 0.35 5.7 185 0.70 11.0
49 0.41 7.1 245 0.78 2.8
50 0.43 5.8 180 0.87 11.8
51 0.60 7.4 240 0.71 17.1
52 0.48 6.8 225 0.70 11.6
53 0.46 6.8 215 0.73 5.8
54 0.44 7.0 230 0.76 8.3
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X1 X2 X4 `X5\\r`
0.046321 0.035851 0.000533 0.022167 0.003308
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.105568 -0.028352 0.002931 0.021291 0.111171
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.0463208 0.1254681 0.369 0.71358
X1 0.0358511 0.0250111 1.433 0.15809
X2 0.0005330 0.0003834 1.390 0.17074
X4 0.0221667 0.0673110 0.329 0.74332
`X5\\r` 0.0033081 0.0011061 2.991 0.00434 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.04578 on 49 degrees of freedom
Multiple R-squared: 0.3997, Adjusted R-squared: 0.3507
F-statistic: 8.158 on 4 and 49 DF, p-value: 4.005e-05
> 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.15938344 0.31876688 0.8406166
[2,] 0.18800481 0.37600961 0.8119952
[3,] 0.12291940 0.24583881 0.8770806
[4,] 0.21832702 0.43665403 0.7816730
[5,] 0.15465360 0.30930719 0.8453464
[6,] 0.28583647 0.57167295 0.7141635
[7,] 0.19645053 0.39290106 0.8035495
[8,] 0.13048114 0.26096229 0.8695189
[9,] 0.09765604 0.19531208 0.9023440
[10,] 0.09444565 0.18889130 0.9055544
[11,] 0.07589635 0.15179271 0.9241036
[12,] 0.16571232 0.33142465 0.8342877
[13,] 0.13363246 0.26726491 0.8663675
[14,] 0.16609126 0.33218253 0.8339087
[15,] 0.11757368 0.23514736 0.8824263
[16,] 0.13788881 0.27577762 0.8621112
[17,] 0.11151876 0.22303753 0.8884812
[18,] 0.10509710 0.21019419 0.8949029
[19,] 0.07047445 0.14094891 0.9295255
[20,] 0.09266230 0.18532460 0.9073377
[21,] 0.09909209 0.19818419 0.9009079
[22,] 0.07643067 0.15286133 0.9235693
[23,] 0.05424707 0.10849415 0.9457529
[24,] 0.04941829 0.09883658 0.9505817
[25,] 0.03387734 0.06775468 0.9661227
[26,] 0.02391850 0.04783700 0.9760815
[27,] 0.01412512 0.02825024 0.9858749
[28,] 0.05386368 0.10772736 0.9461363
[29,] 0.06890906 0.13781811 0.9310909
[30,] 0.04955641 0.09911283 0.9504436
[31,] 0.03074471 0.06148942 0.9692553
[32,] 0.02163145 0.04326290 0.9783686
[33,] 0.01694852 0.03389704 0.9830515
[34,] 0.01165600 0.02331200 0.9883440
[35,] 0.04480188 0.08960376 0.9551981
[36,] 0.02712187 0.05424373 0.9728781
[37,] 0.08441510 0.16883021 0.9155849
[38,] 0.22226068 0.44452135 0.7777393
[39,] 0.70452858 0.59094284 0.2954714
> postscript(file="/var/fisher/rcomp/tmp/1ltc61355171196.ps",horizontal=F,onefile=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/fisher/rcomp/tmp/2swiw1355171196.ps",horizontal=F,onefile=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/fisher/rcomp/tmp/3xsqg1355171196.ps",horizontal=F,onefile=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/fisher/rcomp/tmp/44cdo1355171196.ps",horizontal=F,onefile=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/fisher/rcomp/tmp/53zc11355171196.ps",horizontal=F,onefile=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 = 54
Frequency = 1
1 2 3 4 5 6
-0.015327417 0.015432791 0.013840830 0.012597997 -0.049663472 -0.073632558
7 8 9 10 11 12
0.071372041 -0.014926116 -0.012681470 0.021789834 0.015255260 0.018762590
13 14 15 16 17 18
-0.043328246 0.010636267 -0.008988345 -0.009585099 0.035324398 -0.032273923
19 20 21 22 23 24
-0.069347522 0.020734540 -0.060521371 -0.006670438 0.041984643 -0.031326143
25 26 27 28 29 30
0.035614013 0.004426722 0.054643774 -0.047229529 0.017830390 -0.010744644
31 32 33 34 35 36
0.047852502 -0.006140991 0.026960384 -0.006670438 -0.079825169 0.059924610
37 38 39 40 41 42
0.028462943 0.015451426 -0.028854247 -0.041093109 -0.026843714 -0.105567510
43 44 45 46 47 48
0.006117285 0.092858467 -0.021468143 0.001434888 0.111170642 -0.051189356
49 50 51 52 53 54
-0.048009999 0.021475885 0.088145793 0.016068163 0.019920432 -0.024180539
> postscript(file="/var/fisher/rcomp/tmp/6vkkq1355171196.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> dum <- cbind(lag(myerror,k=1),myerror)
> dum
Time Series:
Start = 0
End = 54
Frequency = 1
lag(myerror, k = 1) myerror
0 -0.015327417 NA
1 0.015432791 -0.015327417
2 0.013840830 0.015432791
3 0.012597997 0.013840830
4 -0.049663472 0.012597997
5 -0.073632558 -0.049663472
6 0.071372041 -0.073632558
7 -0.014926116 0.071372041
8 -0.012681470 -0.014926116
9 0.021789834 -0.012681470
10 0.015255260 0.021789834
11 0.018762590 0.015255260
12 -0.043328246 0.018762590
13 0.010636267 -0.043328246
14 -0.008988345 0.010636267
15 -0.009585099 -0.008988345
16 0.035324398 -0.009585099
17 -0.032273923 0.035324398
18 -0.069347522 -0.032273923
19 0.020734540 -0.069347522
20 -0.060521371 0.020734540
21 -0.006670438 -0.060521371
22 0.041984643 -0.006670438
23 -0.031326143 0.041984643
24 0.035614013 -0.031326143
25 0.004426722 0.035614013
26 0.054643774 0.004426722
27 -0.047229529 0.054643774
28 0.017830390 -0.047229529
29 -0.010744644 0.017830390
30 0.047852502 -0.010744644
31 -0.006140991 0.047852502
32 0.026960384 -0.006140991
33 -0.006670438 0.026960384
34 -0.079825169 -0.006670438
35 0.059924610 -0.079825169
36 0.028462943 0.059924610
37 0.015451426 0.028462943
38 -0.028854247 0.015451426
39 -0.041093109 -0.028854247
40 -0.026843714 -0.041093109
41 -0.105567510 -0.026843714
42 0.006117285 -0.105567510
43 0.092858467 0.006117285
44 -0.021468143 0.092858467
45 0.001434888 -0.021468143
46 0.111170642 0.001434888
47 -0.051189356 0.111170642
48 -0.048009999 -0.051189356
49 0.021475885 -0.048009999
50 0.088145793 0.021475885
51 0.016068163 0.088145793
52 0.019920432 0.016068163
53 -0.024180539 0.019920432
54 NA -0.024180539
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 0.015432791 -0.015327417
[2,] 0.013840830 0.015432791
[3,] 0.012597997 0.013840830
[4,] -0.049663472 0.012597997
[5,] -0.073632558 -0.049663472
[6,] 0.071372041 -0.073632558
[7,] -0.014926116 0.071372041
[8,] -0.012681470 -0.014926116
[9,] 0.021789834 -0.012681470
[10,] 0.015255260 0.021789834
[11,] 0.018762590 0.015255260
[12,] -0.043328246 0.018762590
[13,] 0.010636267 -0.043328246
[14,] -0.008988345 0.010636267
[15,] -0.009585099 -0.008988345
[16,] 0.035324398 -0.009585099
[17,] -0.032273923 0.035324398
[18,] -0.069347522 -0.032273923
[19,] 0.020734540 -0.069347522
[20,] -0.060521371 0.020734540
[21,] -0.006670438 -0.060521371
[22,] 0.041984643 -0.006670438
[23,] -0.031326143 0.041984643
[24,] 0.035614013 -0.031326143
[25,] 0.004426722 0.035614013
[26,] 0.054643774 0.004426722
[27,] -0.047229529 0.054643774
[28,] 0.017830390 -0.047229529
[29,] -0.010744644 0.017830390
[30,] 0.047852502 -0.010744644
[31,] -0.006140991 0.047852502
[32,] 0.026960384 -0.006140991
[33,] -0.006670438 0.026960384
[34,] -0.079825169 -0.006670438
[35,] 0.059924610 -0.079825169
[36,] 0.028462943 0.059924610
[37,] 0.015451426 0.028462943
[38,] -0.028854247 0.015451426
[39,] -0.041093109 -0.028854247
[40,] -0.026843714 -0.041093109
[41,] -0.105567510 -0.026843714
[42,] 0.006117285 -0.105567510
[43,] 0.092858467 0.006117285
[44,] -0.021468143 0.092858467
[45,] 0.001434888 -0.021468143
[46,] 0.111170642 0.001434888
[47,] -0.051189356 0.111170642
[48,] -0.048009999 -0.051189356
[49,] 0.021475885 -0.048009999
[50,] 0.088145793 0.021475885
[51,] 0.016068163 0.088145793
[52,] 0.019920432 0.016068163
[53,] -0.024180539 0.019920432
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 0.015432791 -0.015327417
2 0.013840830 0.015432791
3 0.012597997 0.013840830
4 -0.049663472 0.012597997
5 -0.073632558 -0.049663472
6 0.071372041 -0.073632558
7 -0.014926116 0.071372041
8 -0.012681470 -0.014926116
9 0.021789834 -0.012681470
10 0.015255260 0.021789834
11 0.018762590 0.015255260
12 -0.043328246 0.018762590
13 0.010636267 -0.043328246
14 -0.008988345 0.010636267
15 -0.009585099 -0.008988345
16 0.035324398 -0.009585099
17 -0.032273923 0.035324398
18 -0.069347522 -0.032273923
19 0.020734540 -0.069347522
20 -0.060521371 0.020734540
21 -0.006670438 -0.060521371
22 0.041984643 -0.006670438
23 -0.031326143 0.041984643
24 0.035614013 -0.031326143
25 0.004426722 0.035614013
26 0.054643774 0.004426722
27 -0.047229529 0.054643774
28 0.017830390 -0.047229529
29 -0.010744644 0.017830390
30 0.047852502 -0.010744644
31 -0.006140991 0.047852502
32 0.026960384 -0.006140991
33 -0.006670438 0.026960384
34 -0.079825169 -0.006670438
35 0.059924610 -0.079825169
36 0.028462943 0.059924610
37 0.015451426 0.028462943
38 -0.028854247 0.015451426
39 -0.041093109 -0.028854247
40 -0.026843714 -0.041093109
41 -0.105567510 -0.026843714
42 0.006117285 -0.105567510
43 0.092858467 0.006117285
44 -0.021468143 0.092858467
45 0.001434888 -0.021468143
46 0.111170642 0.001434888
47 -0.051189356 0.111170642
48 -0.048009999 -0.051189356
49 0.021475885 -0.048009999
50 0.088145793 0.021475885
51 0.016068163 0.088145793
52 0.019920432 0.016068163
53 -0.024180539 0.019920432
> 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/fisher/rcomp/tmp/7k9jw1355171196.ps",horizontal=F,onefile=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/fisher/rcomp/tmp/8hhls1355171196.ps",horizontal=F,onefile=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/fisher/rcomp/tmp/9kjej1355171196.ps",horizontal=F,onefile=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/fisher/rcomp/tmp/10iivq1355171196.ps",horizontal=F,onefile=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/fisher/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/fisher/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/fisher/rcomp/tmp/119k571355171196.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/fisher/rcomp/tmp/122v0x1355171196.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/fisher/rcomp/tmp/13p55x1355171196.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/fisher/rcomp/tmp/14l3rh1355171196.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/fisher/rcomp/tmp/15d7821355171196.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/fisher/rcomp/tmp/16dmch1355171196.tab")
+ }
>
> try(system("convert tmp/1ltc61355171196.ps tmp/1ltc61355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/2swiw1355171196.ps tmp/2swiw1355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/3xsqg1355171196.ps tmp/3xsqg1355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/44cdo1355171196.ps tmp/44cdo1355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/53zc11355171196.ps tmp/53zc11355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/6vkkq1355171196.ps tmp/6vkkq1355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/7k9jw1355171196.ps tmp/7k9jw1355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/8hhls1355171196.ps tmp/8hhls1355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/9kjej1355171196.ps tmp/9kjej1355171196.png",intern=TRUE))
character(0)
> try(system("convert tmp/10iivq1355171196.ps tmp/10iivq1355171196.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
5.795 1.508 7.310