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(216.234
+ ,627
+ ,1.59
+ ,213.586
+ ,696
+ ,1.26
+ ,209.465
+ ,825
+ ,1.13
+ ,204.045
+ ,677
+ ,1.92
+ ,200.237
+ ,656
+ ,2.61
+ ,203.666
+ ,785
+ ,2.26
+ ,241.476
+ ,412
+ ,2.41
+ ,260.307
+ ,352
+ ,2.26
+ ,243.324
+ ,839
+ ,2.03
+ ,244.460
+ ,729
+ ,2.86
+ ,233.575
+ ,696
+ ,2.55
+ ,237.217
+ ,641
+ ,2.27
+ ,235.243
+ ,695
+ ,2.26
+ ,230.354
+ ,638
+ ,2.57
+ ,227.184
+ ,762
+ ,3.07
+ ,221.678
+ ,635
+ ,2.76
+ ,217.142
+ ,721
+ ,2.51
+ ,219.452
+ ,854
+ ,2.87
+ ,256.446
+ ,418
+ ,3.14
+ ,265.845
+ ,367
+ ,3.11
+ ,248.624
+ ,824
+ ,3.16
+ ,241.114
+ ,687
+ ,2.47
+ ,229.245
+ ,601
+ ,2.57
+ ,231.805
+ ,676
+ ,2.89
+ ,219.277
+ ,740
+ ,2.63
+ ,219.313
+ ,691
+ ,2.38
+ ,212.610
+ ,683
+ ,1.69
+ ,214.771
+ ,594
+ ,1.96
+ ,211.142
+ ,729
+ ,2.19
+ ,211.457
+ ,731
+ ,1.87
+ ,240.048
+ ,386
+ ,1.6
+ ,240.636
+ ,331
+ ,1.63
+ ,230.580
+ ,707
+ ,1.22
+ ,208.795
+ ,715
+ ,1.21
+ ,197.922
+ ,657
+ ,1.49
+ ,194.596
+ ,653
+ ,1.64
+ ,194.581
+ ,642
+ ,1.66
+ ,185.686
+ ,643
+ ,1.77
+ ,178.106
+ ,718
+ ,1.82
+ ,172.608
+ ,654
+ ,1.78
+ ,167.302
+ ,632
+ ,1.28
+ ,168.053
+ ,731
+ ,1.29
+ ,202.300
+ ,392
+ ,1.37
+ ,202.388
+ ,344
+ ,1.12
+ ,182.516
+ ,792
+ ,1.51
+ ,173.476
+ ,852
+ ,2.24
+ ,166.444
+ ,649
+ ,2.94
+ ,171.297
+ ,629
+ ,3.09
+ ,169.701
+ ,685
+ ,3.46
+ ,164.182
+ ,617
+ ,3.64
+ ,161.914
+ ,715
+ ,4.39
+ ,159.612
+ ,715
+ ,4.15
+ ,151.001
+ ,629
+ ,5.21
+ ,158.114
+ ,916
+ ,5.8
+ ,186.530
+ ,531
+ ,5.91
+ ,187.069
+ ,357
+ ,5.39
+ ,174.330
+ ,917
+ ,5.46
+ ,169.362
+ ,828
+ ,4.72
+ ,166.827
+ ,708
+ ,3.14
+ ,178.037
+ ,858
+ ,2.63
+ ,186.413
+ ,775
+ ,2.32
+ ,189.226
+ ,785
+ ,1.93
+ ,191.563
+ ,1006
+ ,0.62
+ ,188.906
+ ,789
+ ,0.6
+ ,186.005
+ ,734
+ ,-0.37
+ ,195.309
+ ,906
+ ,-1.1
+ ,223.532
+ ,532
+ ,-1.68
+ ,226.899
+ ,387
+ ,-0.78
+ ,214.126
+ ,991
+ ,-1.19
+ ,206.903
+ ,841
+ ,-0.97
+ ,204.442
+ ,892
+ ,-0.12
+ ,220.375
+ ,782
+ ,0.26)
+ ,dim=c(3
+ ,72)
+ ,dimnames=list(c('werklozen'
+ ,'faillissementen'
+ ,'inflatie')
+ ,1:72))
> y <- array(NA,dim=c(3,72),dimnames=list(c('werklozen','faillissementen','inflatie'),1:72))
> 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
werklozen faillissementen inflatie
1 216.234 627 1.59
2 213.586 696 1.26
3 209.465 825 1.13
4 204.045 677 1.92
5 200.237 656 2.61
6 203.666 785 2.26
7 241.476 412 2.41
8 260.307 352 2.26
9 243.324 839 2.03
10 244.460 729 2.86
11 233.575 696 2.55
12 237.217 641 2.27
13 235.243 695 2.26
14 230.354 638 2.57
15 227.184 762 3.07
16 221.678 635 2.76
17 217.142 721 2.51
18 219.452 854 2.87
19 256.446 418 3.14
20 265.845 367 3.11
21 248.624 824 3.16
22 241.114 687 2.47
23 229.245 601 2.57
24 231.805 676 2.89
25 219.277 740 2.63
26 219.313 691 2.38
27 212.610 683 1.69
28 214.771 594 1.96
29 211.142 729 2.19
30 211.457 731 1.87
31 240.048 386 1.60
32 240.636 331 1.63
33 230.580 707 1.22
34 208.795 715 1.21
35 197.922 657 1.49
36 194.596 653 1.64
37 194.581 642 1.66
38 185.686 643 1.77
39 178.106 718 1.82
40 172.608 654 1.78
41 167.302 632 1.28
42 168.053 731 1.29
43 202.300 392 1.37
44 202.388 344 1.12
45 182.516 792 1.51
46 173.476 852 2.24
47 166.444 649 2.94
48 171.297 629 3.09
49 169.701 685 3.46
50 164.182 617 3.64
51 161.914 715 4.39
52 159.612 715 4.15
53 151.001 629 5.21
54 158.114 916 5.80
55 186.530 531 5.91
56 187.069 357 5.39
57 174.330 917 5.46
58 169.362 828 4.72
59 166.827 708 3.14
60 178.037 858 2.63
61 186.413 775 2.32
62 189.226 785 1.93
63 191.563 1006 0.62
64 188.906 789 0.60
65 186.005 734 -0.37
66 195.309 906 -1.10
67 223.532 532 -1.68
68 226.899 387 -0.78
69 214.126 991 -1.19
70 206.903 841 -0.97
71 204.442 892 -0.12
72 220.375 782 0.26
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) faillissementen inflatie
260.19095 -0.06473 -5.17544
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-45.3574 -18.9415 0.4205 17.7605 58.1218
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 260.19095 14.56527 17.864 < 2e-16 ***
faillissementen -0.06473 0.01944 -3.329 0.00140 **
inflatie -5.17544 1.95542 -2.647 0.01006 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 25.54 on 69 degrees of freedom
Multiple R-squared: 0.1941, Adjusted R-squared: 0.1708
F-statistic: 8.31 on 2 and 69 DF, p-value: 0.0005841
> 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.003253696 0.0065073925 9.967463e-01
[2,] 0.013748132 0.0274962649 9.862519e-01
[3,] 0.011150931 0.0223018611 9.888491e-01
[4,] 0.257676227 0.5153524536 7.423238e-01
[5,] 0.275593703 0.5511874063 7.244063e-01
[6,] 0.204315696 0.4086313913 7.956843e-01
[7,] 0.154047824 0.3080956486 8.459522e-01
[8,] 0.117680838 0.2353616755 8.823192e-01
[9,] 0.080115524 0.1602310480 9.198845e-01
[10,] 0.057268851 0.1145377029 9.427311e-01
[11,] 0.042959209 0.0859184171 9.570408e-01
[12,] 0.029399852 0.0587997041 9.706001e-01
[13,] 0.020605832 0.0412116647 9.793942e-01
[14,] 0.019052191 0.0381043813 9.809478e-01
[15,] 0.026261105 0.0525222091 9.737389e-01
[16,] 0.075652718 0.1513054356 9.243473e-01
[17,] 0.106129796 0.2122595913 8.938702e-01
[18,] 0.107171422 0.2143428439 8.928286e-01
[19,] 0.136049112 0.2720982234 8.639509e-01
[20,] 0.153536578 0.3070731569 8.464634e-01
[21,] 0.165644961 0.3312899216 8.343550e-01
[22,] 0.145418875 0.2908377495 8.545811e-01
[23,] 0.142803409 0.2856068182 8.571966e-01
[24,] 0.149417039 0.2988340778 8.505830e-01
[25,] 0.143494504 0.2869890072 8.565055e-01
[26,] 0.194704804 0.3894096081 8.052952e-01
[27,] 0.302388218 0.6047764356 6.976118e-01
[28,] 0.584001700 0.8319966009 4.159983e-01
[29,] 0.583205934 0.8335881318 4.167941e-01
[30,] 0.608639799 0.7827204029 3.913602e-01
[31,] 0.654375170 0.6912496601 3.456248e-01
[32,] 0.687096599 0.6258068012 3.129034e-01
[33,] 0.770151693 0.4596966139 2.298483e-01
[34,] 0.843838186 0.3123236283 1.561618e-01
[35,] 0.929701775 0.1405964503 7.029823e-02
[36,] 0.977969738 0.0440605236 2.203026e-02
[37,] 0.990352410 0.0192951799 9.647590e-03
[38,] 0.986223660 0.0275526807 1.377634e-02
[39,] 0.979161161 0.0416776781 2.083884e-02
[40,] 0.971201247 0.0575975055 2.879875e-02
[41,] 0.975049067 0.0499018657 2.495093e-02
[42,] 0.995480865 0.0090382699 4.519135e-03
[43,] 0.998273832 0.0034523368 1.726168e-03
[44,] 0.999030439 0.0019391214 9.695607e-04
[45,] 0.999634709 0.0007305815 3.652907e-04
[46,] 0.999711600 0.0005767990 2.883995e-04
[47,] 0.999797120 0.0004057594 2.028797e-04
[48,] 0.999946964 0.0001060716 5.303582e-05
[49,] 0.999870142 0.0002597155 1.298578e-04
[50,] 0.999777926 0.0004441477 2.220738e-04
[51,] 0.999520256 0.0009594878 4.797439e-04
[52,] 0.999442957 0.0011140856 5.570428e-04
[53,] 0.998802138 0.0023957232 1.197862e-03
[54,] 0.998312531 0.0033749377 1.687469e-03
[55,] 0.995791343 0.0084173146 4.208657e-03
[56,] 0.989623311 0.0207533776 1.037669e-02
[57,] 0.975747703 0.0485045936 2.425230e-02
[58,] 0.946319486 0.1073610281 5.368051e-02
[59,] 0.918134327 0.1637313470 8.186567e-02
[60,] 0.976838520 0.0463229590 2.316148e-02
[61,] 0.975129346 0.0497413077 2.487065e-02
> postscript(file="/var/www/html/rcomp/tmp/1ci601292940993.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/www/html/rcomp/tmp/2m95k1292940993.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/www/html/rcomp/tmp/3m95k1292940993.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/www/html/rcomp/tmp/4m95k1292940993.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/www/html/rcomp/tmp/5m95k1292940993.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 = 72
Frequency = 1
1 2 3 4 5 6
4.8553104 4.9655226 8.5213930 -2.3894864 -3.9856847 5.9815896
7 8 9 10 11 12
20.4250399 34.5961535 47.9444528 46.2561856 31.6308362 30.2637740
13 14 15 16 17 18
31.7332334 24.7592265 32.2029914 16.8723810 16.6089733 29.3907125
19 20 21 22 23 24
39.5614659 45.5041676 58.1218038 38.1732657 21.2553579 30.3259612
25 26 27 28 29 30
20.5948231 16.1653811 5.3735201 3.1712581 9.4706429 8.2589554
31 32 33 34 35 36
13.1220553 10.3053785 22.4644931 1.1455482 -12.0324479 -14.8410371
37 38 39 40 41 42
-15.4645164 -23.7254921 -26.1922568 -36.0397498 -45.3574442 -38.1467980
43 44 45 46 47 48
-25.4279381 -29.7406540 -18.5969049 -19.9752651 -36.5238738 -32.1890819
49 50 51 52 53 54
-28.2455041 -37.2343057 -29.2775624 -32.8216672 -41.5131556 -12.7702341
55 56 57 58 59 60
-8.7045154 -22.1190979 1.7508436 -12.8076098 -31.2869418 -13.0074876
61 62 63 64 65 66
-11.6081460 -10.1663046 -0.3046412 -17.1107310 -28.5918448 -11.9330107
67 68 69 70 71 72
-10.9193556 -12.2797584 11.9199253 -3.8734056 1.3657509 12.1455372
> postscript(file="/var/www/html/rcomp/tmp/6fjmn1292940993.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 = 72
Frequency = 1
lag(myerror, k = 1) myerror
0 4.8553104 NA
1 4.9655226 4.8553104
2 8.5213930 4.9655226
3 -2.3894864 8.5213930
4 -3.9856847 -2.3894864
5 5.9815896 -3.9856847
6 20.4250399 5.9815896
7 34.5961535 20.4250399
8 47.9444528 34.5961535
9 46.2561856 47.9444528
10 31.6308362 46.2561856
11 30.2637740 31.6308362
12 31.7332334 30.2637740
13 24.7592265 31.7332334
14 32.2029914 24.7592265
15 16.8723810 32.2029914
16 16.6089733 16.8723810
17 29.3907125 16.6089733
18 39.5614659 29.3907125
19 45.5041676 39.5614659
20 58.1218038 45.5041676
21 38.1732657 58.1218038
22 21.2553579 38.1732657
23 30.3259612 21.2553579
24 20.5948231 30.3259612
25 16.1653811 20.5948231
26 5.3735201 16.1653811
27 3.1712581 5.3735201
28 9.4706429 3.1712581
29 8.2589554 9.4706429
30 13.1220553 8.2589554
31 10.3053785 13.1220553
32 22.4644931 10.3053785
33 1.1455482 22.4644931
34 -12.0324479 1.1455482
35 -14.8410371 -12.0324479
36 -15.4645164 -14.8410371
37 -23.7254921 -15.4645164
38 -26.1922568 -23.7254921
39 -36.0397498 -26.1922568
40 -45.3574442 -36.0397498
41 -38.1467980 -45.3574442
42 -25.4279381 -38.1467980
43 -29.7406540 -25.4279381
44 -18.5969049 -29.7406540
45 -19.9752651 -18.5969049
46 -36.5238738 -19.9752651
47 -32.1890819 -36.5238738
48 -28.2455041 -32.1890819
49 -37.2343057 -28.2455041
50 -29.2775624 -37.2343057
51 -32.8216672 -29.2775624
52 -41.5131556 -32.8216672
53 -12.7702341 -41.5131556
54 -8.7045154 -12.7702341
55 -22.1190979 -8.7045154
56 1.7508436 -22.1190979
57 -12.8076098 1.7508436
58 -31.2869418 -12.8076098
59 -13.0074876 -31.2869418
60 -11.6081460 -13.0074876
61 -10.1663046 -11.6081460
62 -0.3046412 -10.1663046
63 -17.1107310 -0.3046412
64 -28.5918448 -17.1107310
65 -11.9330107 -28.5918448
66 -10.9193556 -11.9330107
67 -12.2797584 -10.9193556
68 11.9199253 -12.2797584
69 -3.8734056 11.9199253
70 1.3657509 -3.8734056
71 12.1455372 1.3657509
72 NA 12.1455372
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 4.9655226 4.8553104
[2,] 8.5213930 4.9655226
[3,] -2.3894864 8.5213930
[4,] -3.9856847 -2.3894864
[5,] 5.9815896 -3.9856847
[6,] 20.4250399 5.9815896
[7,] 34.5961535 20.4250399
[8,] 47.9444528 34.5961535
[9,] 46.2561856 47.9444528
[10,] 31.6308362 46.2561856
[11,] 30.2637740 31.6308362
[12,] 31.7332334 30.2637740
[13,] 24.7592265 31.7332334
[14,] 32.2029914 24.7592265
[15,] 16.8723810 32.2029914
[16,] 16.6089733 16.8723810
[17,] 29.3907125 16.6089733
[18,] 39.5614659 29.3907125
[19,] 45.5041676 39.5614659
[20,] 58.1218038 45.5041676
[21,] 38.1732657 58.1218038
[22,] 21.2553579 38.1732657
[23,] 30.3259612 21.2553579
[24,] 20.5948231 30.3259612
[25,] 16.1653811 20.5948231
[26,] 5.3735201 16.1653811
[27,] 3.1712581 5.3735201
[28,] 9.4706429 3.1712581
[29,] 8.2589554 9.4706429
[30,] 13.1220553 8.2589554
[31,] 10.3053785 13.1220553
[32,] 22.4644931 10.3053785
[33,] 1.1455482 22.4644931
[34,] -12.0324479 1.1455482
[35,] -14.8410371 -12.0324479
[36,] -15.4645164 -14.8410371
[37,] -23.7254921 -15.4645164
[38,] -26.1922568 -23.7254921
[39,] -36.0397498 -26.1922568
[40,] -45.3574442 -36.0397498
[41,] -38.1467980 -45.3574442
[42,] -25.4279381 -38.1467980
[43,] -29.7406540 -25.4279381
[44,] -18.5969049 -29.7406540
[45,] -19.9752651 -18.5969049
[46,] -36.5238738 -19.9752651
[47,] -32.1890819 -36.5238738
[48,] -28.2455041 -32.1890819
[49,] -37.2343057 -28.2455041
[50,] -29.2775624 -37.2343057
[51,] -32.8216672 -29.2775624
[52,] -41.5131556 -32.8216672
[53,] -12.7702341 -41.5131556
[54,] -8.7045154 -12.7702341
[55,] -22.1190979 -8.7045154
[56,] 1.7508436 -22.1190979
[57,] -12.8076098 1.7508436
[58,] -31.2869418 -12.8076098
[59,] -13.0074876 -31.2869418
[60,] -11.6081460 -13.0074876
[61,] -10.1663046 -11.6081460
[62,] -0.3046412 -10.1663046
[63,] -17.1107310 -0.3046412
[64,] -28.5918448 -17.1107310
[65,] -11.9330107 -28.5918448
[66,] -10.9193556 -11.9330107
[67,] -12.2797584 -10.9193556
[68,] 11.9199253 -12.2797584
[69,] -3.8734056 11.9199253
[70,] 1.3657509 -3.8734056
[71,] 12.1455372 1.3657509
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 4.9655226 4.8553104
2 8.5213930 4.9655226
3 -2.3894864 8.5213930
4 -3.9856847 -2.3894864
5 5.9815896 -3.9856847
6 20.4250399 5.9815896
7 34.5961535 20.4250399
8 47.9444528 34.5961535
9 46.2561856 47.9444528
10 31.6308362 46.2561856
11 30.2637740 31.6308362
12 31.7332334 30.2637740
13 24.7592265 31.7332334
14 32.2029914 24.7592265
15 16.8723810 32.2029914
16 16.6089733 16.8723810
17 29.3907125 16.6089733
18 39.5614659 29.3907125
19 45.5041676 39.5614659
20 58.1218038 45.5041676
21 38.1732657 58.1218038
22 21.2553579 38.1732657
23 30.3259612 21.2553579
24 20.5948231 30.3259612
25 16.1653811 20.5948231
26 5.3735201 16.1653811
27 3.1712581 5.3735201
28 9.4706429 3.1712581
29 8.2589554 9.4706429
30 13.1220553 8.2589554
31 10.3053785 13.1220553
32 22.4644931 10.3053785
33 1.1455482 22.4644931
34 -12.0324479 1.1455482
35 -14.8410371 -12.0324479
36 -15.4645164 -14.8410371
37 -23.7254921 -15.4645164
38 -26.1922568 -23.7254921
39 -36.0397498 -26.1922568
40 -45.3574442 -36.0397498
41 -38.1467980 -45.3574442
42 -25.4279381 -38.1467980
43 -29.7406540 -25.4279381
44 -18.5969049 -29.7406540
45 -19.9752651 -18.5969049
46 -36.5238738 -19.9752651
47 -32.1890819 -36.5238738
48 -28.2455041 -32.1890819
49 -37.2343057 -28.2455041
50 -29.2775624 -37.2343057
51 -32.8216672 -29.2775624
52 -41.5131556 -32.8216672
53 -12.7702341 -41.5131556
54 -8.7045154 -12.7702341
55 -22.1190979 -8.7045154
56 1.7508436 -22.1190979
57 -12.8076098 1.7508436
58 -31.2869418 -12.8076098
59 -13.0074876 -31.2869418
60 -11.6081460 -13.0074876
61 -10.1663046 -11.6081460
62 -0.3046412 -10.1663046
63 -17.1107310 -0.3046412
64 -28.5918448 -17.1107310
65 -11.9330107 -28.5918448
66 -10.9193556 -11.9330107
67 -12.2797584 -10.9193556
68 11.9199253 -12.2797584
69 -3.8734056 11.9199253
70 1.3657509 -3.8734056
71 12.1455372 1.3657509
> 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/7qslq1292940993.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/www/html/rcomp/tmp/8qslq1292940993.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/www/html/rcomp/tmp/9qslq1292940993.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/www/html/rcomp/tmp/1001kt1292940993.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/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/11m21h1292940993.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/1272zn1292940993.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/133ufw1292940993.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/14pcwk1292940993.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/15sduq1292940993.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/16vwte1292940993.tab")
+ }
>
> try(system("convert tmp/1ci601292940993.ps tmp/1ci601292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/2m95k1292940993.ps tmp/2m95k1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/3m95k1292940993.ps tmp/3m95k1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/4m95k1292940993.ps tmp/4m95k1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/5m95k1292940993.ps tmp/5m95k1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/6fjmn1292940993.ps tmp/6fjmn1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/7qslq1292940993.ps tmp/7qslq1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/8qslq1292940993.ps tmp/8qslq1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/9qslq1292940993.ps tmp/9qslq1292940993.png",intern=TRUE))
character(0)
> try(system("convert tmp/1001kt1292940993.ps tmp/1001kt1292940993.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.658 1.710 7.972