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(0.7905
+ ,0.313
+ ,0.7744
+ ,0.779
+ ,0.7719
+ ,0.364
+ ,0.7905
+ ,0.7744
+ ,0.7811
+ ,0.363
+ ,0.7719
+ ,0.7905
+ ,0.7557
+ ,-0.155
+ ,0.7811
+ ,0.7719
+ ,0.7637
+ ,0.052
+ ,0.7557
+ ,0.7811
+ ,0.7595
+ ,0.568
+ ,0.7637
+ ,0.7557
+ ,0.7471
+ ,0.668
+ ,0.7595
+ ,0.7637
+ ,0.7615
+ ,1.378
+ ,0.7471
+ ,0.7595
+ ,0.7487
+ ,0.252
+ ,0.7615
+ ,0.7471
+ ,0.7389
+ ,-0.402
+ ,0.7487
+ ,0.7615
+ ,0.7337
+ ,-0.05
+ ,0.7389
+ ,0.7487
+ ,0.751
+ ,0.555
+ ,0.7337
+ ,0.7389
+ ,0.7382
+ ,0.05
+ ,0.751
+ ,0.7337
+ ,0.7159
+ ,0.15
+ ,0.7382
+ ,0.751
+ ,0.7542
+ ,0.45
+ ,0.7159
+ ,0.7382
+ ,0.7636
+ ,0.299
+ ,0.7542
+ ,0.7159
+ ,0.7433
+ ,0.199
+ ,0.7636
+ ,0.7542
+ ,0.7658
+ ,0.496
+ ,0.7433
+ ,0.7636
+ ,0.7627
+ ,0.444
+ ,0.7658
+ ,0.7433
+ ,0.748
+ ,-0.393
+ ,0.7627
+ ,0.7658
+ ,0.7692
+ ,-0.444
+ ,0.748
+ ,0.7627
+ ,0.785
+ ,0.198
+ ,0.7692
+ ,0.748
+ ,0.7913
+ ,0.494
+ ,0.785
+ ,0.7692
+ ,0.772
+ ,0.133
+ ,0.7913
+ ,0.785
+ ,0.788
+ ,0.388
+ ,0.772
+ ,0.7913
+ ,0.807
+ ,0.484
+ ,0.788
+ ,0.772
+ ,0.8268
+ ,0.278
+ ,0.807
+ ,0.788
+ ,0.8244
+ ,0.369
+ ,0.8268
+ ,0.807
+ ,0.8487
+ ,0.165
+ ,0.8244
+ ,0.8268
+ ,0.8572
+ ,0.155
+ ,0.8487
+ ,0.8244
+ ,0.8214
+ ,0.087
+ ,0.8572
+ ,0.8487
+ ,0.8827
+ ,0.414
+ ,0.8214
+ ,0.8572
+ ,0.9216
+ ,0.36
+ ,0.8827
+ ,0.8214
+ ,0.8865
+ ,0.975
+ ,0.9216
+ ,0.8827
+ ,0.8816
+ ,0.27
+ ,0.8865
+ ,0.9216
+ ,0.8884
+ ,0.359
+ ,0.8816
+ ,0.8865
+ ,0.9466
+ ,0.169
+ ,0.8884
+ ,0.8816
+ ,0.918
+ ,0.381
+ ,0.9466
+ ,0.8884
+ ,0.9337
+ ,0.154
+ ,0.918
+ ,0.9466
+ ,0.9559
+ ,0.486
+ ,0.9337
+ ,0.918
+ ,0.9626
+ ,0.925
+ ,0.9559
+ ,0.9337
+ ,0.9434
+ ,0.728
+ ,0.9626
+ ,0.9559
+ ,0.8639
+ ,-0.014
+ ,0.9434
+ ,0.9626
+ ,0.7996
+ ,0.046
+ ,0.8639
+ ,0.9434
+ ,0.668
+ ,-0.819
+ ,0.7996
+ ,0.8639
+ ,0.6572
+ ,-1.674
+ ,0.668
+ ,0.7996
+ ,0.6928
+ ,-0.788
+ ,0.6572
+ ,0.668
+ ,0.6438
+ ,0.279
+ ,0.6928
+ ,0.6572
+ ,0.6454
+ ,0.396
+ ,0.6438
+ ,0.6928
+ ,0.6873
+ ,-0.141
+ ,0.6454
+ ,0.6438
+ ,0.7265
+ ,-0.019
+ ,0.6873
+ ,0.6454
+ ,0.7912
+ ,0.099
+ ,0.7265
+ ,0.6873
+ ,0.8114
+ ,0.742
+ ,0.7912
+ ,0.7265
+ ,0.8281
+ ,0.005
+ ,0.8114
+ ,0.7912
+ ,0.8393
+ ,0.448
+ ,0.8281
+ ,0.8114)
+ ,dim=c(4
+ ,55)
+ ,dimnames=list(c('USDOLLAR'
+ ,'Amerikaanse_inflatie'
+ ,'Y[t-1]'
+ ,'Y[t-2]')
+ ,1:55))
> y <- array(NA,dim=c(4,55),dimnames=list(c('USDOLLAR','Amerikaanse_inflatie','Y[t-1]','Y[t-2]'),1:55))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = '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
USDOLLAR Amerikaanse_inflatie Y[t-1] Y[t-2] t
1 0.7905 0.313 0.7744 0.7790 1
2 0.7719 0.364 0.7905 0.7744 2
3 0.7811 0.363 0.7719 0.7905 3
4 0.7557 -0.155 0.7811 0.7719 4
5 0.7637 0.052 0.7557 0.7811 5
6 0.7595 0.568 0.7637 0.7557 6
7 0.7471 0.668 0.7595 0.7637 7
8 0.7615 1.378 0.7471 0.7595 8
9 0.7487 0.252 0.7615 0.7471 9
10 0.7389 -0.402 0.7487 0.7615 10
11 0.7337 -0.050 0.7389 0.7487 11
12 0.7510 0.555 0.7337 0.7389 12
13 0.7382 0.050 0.7510 0.7337 13
14 0.7159 0.150 0.7382 0.7510 14
15 0.7542 0.450 0.7159 0.7382 15
16 0.7636 0.299 0.7542 0.7159 16
17 0.7433 0.199 0.7636 0.7542 17
18 0.7658 0.496 0.7433 0.7636 18
19 0.7627 0.444 0.7658 0.7433 19
20 0.7480 -0.393 0.7627 0.7658 20
21 0.7692 -0.444 0.7480 0.7627 21
22 0.7850 0.198 0.7692 0.7480 22
23 0.7913 0.494 0.7850 0.7692 23
24 0.7720 0.133 0.7913 0.7850 24
25 0.7880 0.388 0.7720 0.7913 25
26 0.8070 0.484 0.7880 0.7720 26
27 0.8268 0.278 0.8070 0.7880 27
28 0.8244 0.369 0.8268 0.8070 28
29 0.8487 0.165 0.8244 0.8268 29
30 0.8572 0.155 0.8487 0.8244 30
31 0.8214 0.087 0.8572 0.8487 31
32 0.8827 0.414 0.8214 0.8572 32
33 0.9216 0.360 0.8827 0.8214 33
34 0.8865 0.975 0.9216 0.8827 34
35 0.8816 0.270 0.8865 0.9216 35
36 0.8884 0.359 0.8816 0.8865 36
37 0.9466 0.169 0.8884 0.8816 37
38 0.9180 0.381 0.9466 0.8884 38
39 0.9337 0.154 0.9180 0.9466 39
40 0.9559 0.486 0.9337 0.9180 40
41 0.9626 0.925 0.9559 0.9337 41
42 0.9434 0.728 0.9626 0.9559 42
43 0.8639 -0.014 0.9434 0.9626 43
44 0.7996 0.046 0.8639 0.9434 44
45 0.6680 -0.819 0.7996 0.8639 45
46 0.6572 -1.674 0.6680 0.7996 46
47 0.6928 -0.788 0.6572 0.6680 47
48 0.6438 0.279 0.6928 0.6572 48
49 0.6454 0.396 0.6438 0.6928 49
50 0.6873 -0.141 0.6454 0.6438 50
51 0.7265 -0.019 0.6873 0.6454 51
52 0.7912 0.099 0.7265 0.6873 52
53 0.8114 0.742 0.7912 0.7265 53
54 0.8281 0.005 0.8114 0.7912 54
55 0.8393 0.448 0.8281 0.8114 55
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Amerikaanse_inflatie `Y[t-1]`
0.1119357 0.0182913 1.0508570
`Y[t-2]` t
-0.2053135 0.0002817
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.104525 -0.017044 0.001070 0.015781 0.068574
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.1119357 0.0458104 2.443 0.0181 *
Amerikaanse_inflatie 0.0182913 0.0118453 1.544 0.1288
`Y[t-1]` 1.0508570 0.1611941 6.519 3.39e-08 ***
`Y[t-2]` -0.2053135 0.1498345 -1.370 0.1767
t 0.0002817 0.0002876 0.980 0.3320
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.03173 on 50 degrees of freedom
Multiple R-squared: 0.8543, Adjusted R-squared: 0.8426
F-statistic: 73.28 on 4 and 50 DF, p-value: < 2.2e-16
> 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.920149e-03 5.840297e-03 0.9970799
[2,] 6.015243e-03 1.203049e-02 0.9939848
[3,] 1.292593e-03 2.585187e-03 0.9987074
[4,] 2.600225e-04 5.200449e-04 0.9997400
[5,] 1.843540e-04 3.687080e-04 0.9998156
[6,] 1.141135e-04 2.282269e-04 0.9998859
[7,] 3.759486e-05 7.518973e-05 0.9999624
[8,] 6.596196e-05 1.319239e-04 0.9999340
[9,] 9.134999e-04 1.827000e-03 0.9990865
[10,] 5.439271e-04 1.087854e-03 0.9994561
[11,] 5.433196e-04 1.086639e-03 0.9994567
[12,] 2.512515e-04 5.025029e-04 0.9997487
[13,] 9.717049e-05 1.943410e-04 0.9999028
[14,] 1.219004e-04 2.438009e-04 0.9998781
[15,] 9.605315e-05 1.921063e-04 0.9999039
[16,] 3.696251e-05 7.392501e-05 0.9999630
[17,] 2.815469e-05 5.630937e-05 0.9999718
[18,] 1.020416e-05 2.040832e-05 0.9999898
[19,] 5.969686e-06 1.193937e-05 0.9999940
[20,] 4.000215e-06 8.000430e-06 0.9999960
[21,] 1.641065e-06 3.282130e-06 0.9999984
[22,] 8.175265e-07 1.635053e-06 0.9999992
[23,] 2.837041e-07 5.674082e-07 0.9999997
[24,] 1.865153e-06 3.730306e-06 0.9999981
[25,] 5.824156e-06 1.164831e-05 0.9999942
[26,] 2.313567e-05 4.627133e-05 0.9999769
[27,] 1.115412e-04 2.230824e-04 0.9998885
[28,] 6.058051e-05 1.211610e-04 0.9999394
[29,] 2.438405e-05 4.876810e-05 0.9999756
[30,] 2.600138e-04 5.200275e-04 0.9997400
[31,] 1.792759e-04 3.585518e-04 0.9998207
[32,] 1.977478e-04 3.954957e-04 0.9998023
[33,] 2.865698e-04 5.731396e-04 0.9997134
[34,] 4.940557e-04 9.881114e-04 0.9995059
[35,] 3.707754e-03 7.415507e-03 0.9962922
[36,] 2.643692e-02 5.287385e-02 0.9735631
[37,] 7.281323e-01 5.437355e-01 0.2718677
[38,] 7.492350e-01 5.015300e-01 0.2507650
[39,] 6.379981e-01 7.240038e-01 0.3620019
[40,] 8.172846e-01 3.654309e-01 0.1827154
> postscript(file="/var/www/html/rcomp/tmp/1nspx1260707247.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/2s2nu1260707247.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/3v3361260707247.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/4gwr01260707247.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/55lge1260707247.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 = 55
Frequency = 1
1 2 3 4 5
1.871305e-02 -1.896471e-02 1.282339e-02 -1.687009e-02 1.564258e-02
6 7 8 9 10
-1.189923e-02 -2.035393e-02 -7.054125e-03 -1.721800e-02 1.070334e-03
11 12 13 14 15
-3.179497e-03 6.224966e-03 -1.686705e-02 -2.427496e-02 2.906207e-02
16 17 18 19 20
-3.883925e-03 -2.465101e-02 1.539714e-02 -1.484553e-02 -6.640160e-03
21 22 23 24 25
3.002215e-02 8.501175e-03 -3.145621e-03 -1.950057e-02 1.312848e-02
26 27 28 29 30
9.314582e-03 1.591966e-02 -5.332537e-03 2.900448e-02 1.137715e-02
31 32 33 34 35
-2.740388e-02 6.699903e-02 3.483733e-02 -4.008612e-02 1.249937e-02
36 37 38 39 40
1.533246e-02 6.857428e-02 -2.394890e-02 3.762532e-02 3.110050e-02
41 42 43 44 45
9.383341e-03 -8.977723e-03 -5.363518e-02 -3.971322e-02 -1.045252e-01
46 47 48 49 50
2.512329e-02 2.856551e-02 -7.986090e-02 -2.188150e-02 1.781753e-02
51 52 53 54 55
1.080191e-02 4.047091e-02 -1.131424e-02 1.064127e-02 5.456293e-05
> postscript(file="/var/www/html/rcomp/tmp/6gkn11260707247.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 = 55
Frequency = 1
lag(myerror, k = 1) myerror
0 1.871305e-02 NA
1 -1.896471e-02 1.871305e-02
2 1.282339e-02 -1.896471e-02
3 -1.687009e-02 1.282339e-02
4 1.564258e-02 -1.687009e-02
5 -1.189923e-02 1.564258e-02
6 -2.035393e-02 -1.189923e-02
7 -7.054125e-03 -2.035393e-02
8 -1.721800e-02 -7.054125e-03
9 1.070334e-03 -1.721800e-02
10 -3.179497e-03 1.070334e-03
11 6.224966e-03 -3.179497e-03
12 -1.686705e-02 6.224966e-03
13 -2.427496e-02 -1.686705e-02
14 2.906207e-02 -2.427496e-02
15 -3.883925e-03 2.906207e-02
16 -2.465101e-02 -3.883925e-03
17 1.539714e-02 -2.465101e-02
18 -1.484553e-02 1.539714e-02
19 -6.640160e-03 -1.484553e-02
20 3.002215e-02 -6.640160e-03
21 8.501175e-03 3.002215e-02
22 -3.145621e-03 8.501175e-03
23 -1.950057e-02 -3.145621e-03
24 1.312848e-02 -1.950057e-02
25 9.314582e-03 1.312848e-02
26 1.591966e-02 9.314582e-03
27 -5.332537e-03 1.591966e-02
28 2.900448e-02 -5.332537e-03
29 1.137715e-02 2.900448e-02
30 -2.740388e-02 1.137715e-02
31 6.699903e-02 -2.740388e-02
32 3.483733e-02 6.699903e-02
33 -4.008612e-02 3.483733e-02
34 1.249937e-02 -4.008612e-02
35 1.533246e-02 1.249937e-02
36 6.857428e-02 1.533246e-02
37 -2.394890e-02 6.857428e-02
38 3.762532e-02 -2.394890e-02
39 3.110050e-02 3.762532e-02
40 9.383341e-03 3.110050e-02
41 -8.977723e-03 9.383341e-03
42 -5.363518e-02 -8.977723e-03
43 -3.971322e-02 -5.363518e-02
44 -1.045252e-01 -3.971322e-02
45 2.512329e-02 -1.045252e-01
46 2.856551e-02 2.512329e-02
47 -7.986090e-02 2.856551e-02
48 -2.188150e-02 -7.986090e-02
49 1.781753e-02 -2.188150e-02
50 1.080191e-02 1.781753e-02
51 4.047091e-02 1.080191e-02
52 -1.131424e-02 4.047091e-02
53 1.064127e-02 -1.131424e-02
54 5.456293e-05 1.064127e-02
55 NA 5.456293e-05
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -1.896471e-02 0.018713054
[2,] 1.282339e-02 -0.018964715
[3,] -1.687009e-02 0.012823392
[4,] 1.564258e-02 -0.016870093
[5,] -1.189923e-02 0.015642583
[6,] -2.035393e-02 -0.011899229
[7,] -7.054125e-03 -0.020353926
[8,] -1.721800e-02 -0.007054125
[9,] 1.070334e-03 -0.017218000
[10,] -3.179497e-03 0.001070334
[11,] 6.224966e-03 -0.003179497
[12,] -1.686705e-02 0.006224966
[13,] -2.427496e-02 -0.016867046
[14,] 2.906207e-02 -0.024274957
[15,] -3.883925e-03 0.029062072
[16,] -2.465101e-02 -0.003883925
[17,] 1.539714e-02 -0.024651012
[18,] -1.484553e-02 0.015397138
[19,] -6.640160e-03 -0.014845533
[20,] 3.002215e-02 -0.006640160
[21,] 8.501175e-03 0.030022151
[22,] -3.145621e-03 0.008501175
[23,] -1.950057e-02 -0.003145621
[24,] 1.312848e-02 -0.019500572
[25,] 9.314582e-03 0.013128484
[26,] 1.591966e-02 0.009314582
[27,] -5.332537e-03 0.015919656
[28,] 2.900448e-02 -0.005332537
[29,] 1.137715e-02 0.029004485
[30,] -2.740388e-02 0.011377149
[31,] 6.699903e-02 -0.027403879
[32,] 3.483733e-02 0.066999032
[33,] -4.008612e-02 0.034837334
[34,] 1.249937e-02 -0.040086116
[35,] 1.533246e-02 0.012499368
[36,] 6.857428e-02 0.015332462
[37,] -2.394890e-02 0.068574277
[38,] 3.762532e-02 -0.023948899
[39,] 3.110050e-02 0.037625316
[40,] 9.383341e-03 0.031100504
[41,] -8.977723e-03 0.009383341
[42,] -5.363518e-02 -0.008977723
[43,] -3.971322e-02 -0.053635182
[44,] -1.045252e-01 -0.039713224
[45,] 2.512329e-02 -0.104525228
[46,] 2.856551e-02 0.025123295
[47,] -7.986090e-02 0.028565506
[48,] -2.188150e-02 -0.079860898
[49,] 1.781753e-02 -0.021881501
[50,] 1.080191e-02 0.017817531
[51,] 4.047091e-02 0.010801912
[52,] -1.131424e-02 0.040470908
[53,] 1.064127e-02 -0.011314238
[54,] 5.456293e-05 0.010641268
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -1.896471e-02 0.018713054
2 1.282339e-02 -0.018964715
3 -1.687009e-02 0.012823392
4 1.564258e-02 -0.016870093
5 -1.189923e-02 0.015642583
6 -2.035393e-02 -0.011899229
7 -7.054125e-03 -0.020353926
8 -1.721800e-02 -0.007054125
9 1.070334e-03 -0.017218000
10 -3.179497e-03 0.001070334
11 6.224966e-03 -0.003179497
12 -1.686705e-02 0.006224966
13 -2.427496e-02 -0.016867046
14 2.906207e-02 -0.024274957
15 -3.883925e-03 0.029062072
16 -2.465101e-02 -0.003883925
17 1.539714e-02 -0.024651012
18 -1.484553e-02 0.015397138
19 -6.640160e-03 -0.014845533
20 3.002215e-02 -0.006640160
21 8.501175e-03 0.030022151
22 -3.145621e-03 0.008501175
23 -1.950057e-02 -0.003145621
24 1.312848e-02 -0.019500572
25 9.314582e-03 0.013128484
26 1.591966e-02 0.009314582
27 -5.332537e-03 0.015919656
28 2.900448e-02 -0.005332537
29 1.137715e-02 0.029004485
30 -2.740388e-02 0.011377149
31 6.699903e-02 -0.027403879
32 3.483733e-02 0.066999032
33 -4.008612e-02 0.034837334
34 1.249937e-02 -0.040086116
35 1.533246e-02 0.012499368
36 6.857428e-02 0.015332462
37 -2.394890e-02 0.068574277
38 3.762532e-02 -0.023948899
39 3.110050e-02 0.037625316
40 9.383341e-03 0.031100504
41 -8.977723e-03 0.009383341
42 -5.363518e-02 -0.008977723
43 -3.971322e-02 -0.053635182
44 -1.045252e-01 -0.039713224
45 2.512329e-02 -0.104525228
46 2.856551e-02 0.025123295
47 -7.986090e-02 0.028565506
48 -2.188150e-02 -0.079860898
49 1.781753e-02 -0.021881501
50 1.080191e-02 0.017817531
51 4.047091e-02 0.010801912
52 -1.131424e-02 0.040470908
53 1.064127e-02 -0.011314238
54 5.456293e-05 0.010641268
> 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/7fg5c1260707247.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/8f9ri1260707247.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/93qy41260707247.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/10ap5x1260707247.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/117bh71260707247.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/12vsm31260707247.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/13gk9i1260707247.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/1490301260707247.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/15419o1260707247.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/16w0eg1260707247.tab")
+ }
>
> try(system("convert tmp/1nspx1260707247.ps tmp/1nspx1260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/2s2nu1260707247.ps tmp/2s2nu1260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/3v3361260707247.ps tmp/3v3361260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/4gwr01260707247.ps tmp/4gwr01260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/55lge1260707247.ps tmp/55lge1260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/6gkn11260707247.ps tmp/6gkn11260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/7fg5c1260707247.ps tmp/7fg5c1260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/8f9ri1260707247.ps tmp/8f9ri1260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/93qy41260707247.ps tmp/93qy41260707247.png",intern=TRUE))
character(0)
> try(system("convert tmp/10ap5x1260707247.ps tmp/10ap5x1260707247.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.377 1.543 3.339