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(99.2
+ ,96.7
+ ,101.0
+ ,99.0
+ ,98.1
+ ,100.1
+ ,631
+ ,923
+ ,-12
+ ,-10.8
+ ,654
+ ,294
+ ,-13
+ ,-12.2
+ ,671
+ ,833
+ ,-16
+ ,-14.1
+ ,586
+ ,840
+ ,-10
+ ,-15.2
+ ,600
+ ,969
+ ,-4
+ ,-15.8
+ ,625
+ ,568
+ ,-9
+ ,-15.8
+ ,558
+ ,110
+ ,-8
+ ,-14.9
+ ,630
+ ,577
+ ,-9
+ ,-12.6
+ ,628
+ ,654
+ ,-3
+ ,-9.9
+ ,603
+ ,184
+ ,-13
+ ,-7.8
+ ,656
+ ,255
+ ,-3
+ ,-6
+ ,600
+ ,730
+ ,-1
+ ,-5
+ ,670
+ ,326
+ ,-2
+ ,-4.5
+ ,678
+ ,423
+ ,0
+ ,-3.9
+ ,641
+ ,502
+ ,0
+ ,-2.9
+ ,625
+ ,311
+ ,-3
+ ,-1.5
+ ,628
+ ,177
+ ,0
+ ,-0.5
+ ,589
+ ,767
+ ,5
+ ,0
+ ,582
+ ,471
+ ,3
+ ,0.5
+ ,636
+ ,248
+ ,4
+ ,0.9
+ ,599
+ ,885
+ ,3
+ ,0.8
+ ,621
+ ,694
+ ,1
+ ,0.1
+ ,637
+ ,406
+ ,-1
+ ,-1
+ ,595
+ ,994
+ ,0
+ ,-2
+ ,696
+ ,308
+ ,-2
+ ,-3
+ ,674
+ ,201
+ ,-1
+ ,-3.7
+ ,648
+ ,861
+ ,2
+ ,-4.7
+ ,649
+ ,605
+ ,0
+ ,-6.4
+ ,672
+ ,392
+ ,-6
+ ,-7.5
+ ,598
+ ,396
+ ,-7
+ ,-7.8
+ ,613
+ ,177
+ ,-6
+ ,-7.7
+ ,638
+ ,104
+ ,-4
+ ,-6.6
+ ,615
+ ,632
+ ,-9
+ ,-4.2
+ ,634
+ ,465
+ ,-2
+ ,-2
+ ,638
+ ,686
+ ,-3
+ ,-0.7
+ ,604
+ ,243
+ ,2
+ ,0.1
+ ,706
+ ,669
+ ,3
+ ,0.9
+ ,677
+ ,185
+ ,1
+ ,2.1
+ ,644
+ ,328
+ ,0
+ ,3.5
+ ,644
+ ,825
+ ,1
+ ,4.9
+ ,605
+ ,707
+ ,1
+ ,5.7
+ ,600
+ ,136
+ ,3
+ ,6.2
+ ,612
+ ,166
+ ,5
+ ,6.5
+ ,599
+ ,659
+ ,5
+ ,6.5
+ ,634
+ ,210
+ ,4
+ ,6.3
+ ,618
+ ,234
+ ,11
+ ,6.2
+ ,613
+ ,576
+ ,8
+ ,6.4
+ ,627
+ ,200
+ ,-1
+ ,6.3
+ ,668
+ ,973
+ ,4
+ ,5.8
+ ,651
+ ,479
+ ,4
+ ,5.1
+ ,619
+ ,661
+ ,4
+ ,5.1
+ ,644
+ ,260
+ ,6
+ ,5.8
+ ,579
+ ,936
+ ,6
+ ,6.7
+ ,601
+ ,752
+ ,6
+ ,7.1
+ ,595
+ ,376
+ ,6
+ ,6.7
+ ,588
+ ,902
+ ,4
+ ,5.5
+ ,634
+ ,341
+ ,1
+ ,4.2
+ ,594
+ ,305
+ ,6
+ ,3
+ ,606
+ ,200
+ ,0
+ ,2.2
+ ,610
+ ,926
+ ,2
+ ,2
+ ,633
+ ,685
+ ,-2
+ ,1.8
+ ,639
+ ,696
+ ,0
+ ,1.8
+ ,659
+ ,451
+ ,1
+ ,1.5
+ ,593
+ ,248
+ ,-3
+ ,0.4
+ ,606
+ ,677
+ ,-3
+ ,-0.9
+ ,599
+ ,434
+ ,-5
+ ,-1.7
+ ,569
+ ,578
+ ,-7
+ ,-2.6
+ ,629
+ ,873
+ ,-7
+ ,-4.4
+ ,613
+ ,438
+ ,-5
+ ,-8.3
+ ,604
+ ,172
+ ,-13
+ ,-14.4
+ ,658
+ ,328
+ ,-16
+ ,-21.3
+ ,612
+ ,633
+ ,-20
+ ,-26.5
+ ,707
+ ,372
+ ,-18
+ ,-29.2
+ ,739
+ ,770
+ ,-21
+ ,-30.8
+ ,777
+ ,535
+ ,-20
+ ,-30.9
+ ,685
+ ,030
+ ,-16
+ ,-29.5
+ ,730
+ ,234
+ ,-14
+ ,-27.1
+ ,714
+ ,154
+ ,-12
+ ,-24.4
+ ,630
+ ,872
+ ,-10
+ ,-21.9
+ ,719
+ ,492
+ ,-3
+ ,-19.3
+ ,677
+ ,023
+ ,-4
+ ,-17
+ ,679
+ ,272
+ ,-4
+ ,-13.8
+ ,718
+ ,317
+ ,-1
+ ,-9.9
+ ,645
+ ,672
+ ,-8
+ ,-7.9)
+ ,dim=c(3
+ ,86)
+ ,dimnames=list(c('Werkloosheid'
+ ,'Consumenten'
+ ,'Ondernemers')
+ ,1:86))
> y <- array(NA,dim=c(3,86),dimnames=list(c('Werkloosheid','Consumenten','Ondernemers'),1:86))
> 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
Werkloosheid Consumenten Ondernemers
1 99.2 96.7 101.0
2 99.0 98.1 100.1
3 631.0 923.0 -12.0
4 -10.8 654.0 294.0
5 -13.0 -12.2 671.0
6 833.0 -16.0 -14.1
7 586.0 840.0 -10.0
8 -15.2 600.0 969.0
9 -4.0 -15.8 625.0
10 568.0 -9.0 -15.8
11 558.0 110.0 -8.0
12 -14.9 630.0 577.0
13 -9.0 -12.6 628.0
14 654.0 -3.0 -9.9
15 603.0 184.0 -13.0
16 -7.8 656.0 255.0
17 -3.0 -6.0 600.0
18 730.0 -1.0 -5.0
19 670.0 326.0 -2.0
20 -4.5 678.0 423.0
21 0.0 -3.9 641.0
22 502.0 0.0 -2.9
23 625.0 311.0 -3.0
24 -1.5 628.0 177.0
25 0.0 -0.5 589.0
26 767.0 5.0 0.0
27 582.0 471.0 3.0
28 0.5 636.0 248.0
29 4.0 0.9 599.0
30 885.0 3.0 0.8
31 621.0 694.0 1.0
32 0.1 637.0 406.0
33 -1.0 -1.0 595.0
34 994.0 0.0 -2.0
35 696.0 308.0 -2.0
36 -3.0 674.0 201.0
37 -1.0 -3.7 648.0
38 861.0 2.0 -4.7
39 649.0 605.0 0.0
40 -6.4 672.0 392.0
41 -6.0 -7.5 598.0
42 396.0 -7.0 -7.8
43 613.0 177.0 -6.0
44 -7.7 638.0 104.0
45 -4.0 -6.6 615.0
46 632.0 -9.0 -4.2
47 634.0 465.0 -2.0
48 -2.0 638.0 686.0
49 -3.0 -0.7 604.0
50 243.0 2.0 0.1
51 706.0 669.0 3.0
52 0.9 677.0 185.0
53 1.0 2.1 644.0
54 328.0 0.0 3.5
55 644.0 825.0 1.0
56 4.9 605.0 707.0
57 1.0 5.7 600.0
58 136.0 3.0 6.2
59 612.0 166.0 5.0
60 6.5 599.0 659.0
61 5.0 6.5 634.0
62 210.0 4.0 6.3
63 618.0 234.0 11.0
64 6.2 613.0 576.0
65 8.0 6.4 627.0
66 200.0 -1.0 6.3
67 668.0 973.0 4.0
68 5.8 651.0 479.0
69 4.0 5.1 619.0
70 661.0 4.0 5.1
71 644.0 260.0 6.0
72 5.8 579.0 936.0
73 6.0 6.7 601.0
74 752.0 6.0 7.1
75 595.0 376.0 6.0
76 6.7 588.0 902.0
77 4.0 5.5 634.0
78 341.0 1.0 4.2
79 594.0 305.0 6.0
80 3.0 606.0 200.0
81 0.0 2.2 610.0
82 926.0 2.0 2.0
83 633.0 685.0 -2.0
84 1.8 639.0 696.0
85 0.0 1.8 659.0
86 451.0 1.0 1.5
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Consumenten Ondernemers
548.48570 -0.08245 -0.85219
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-414.956 -72.849 1.963 120.846 443.810
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 548.48570 33.92409 16.168 <2e-16 ***
Consumenten -0.08245 0.06847 -1.204 0.232
Ondernemers -0.85219 0.06910 -12.332 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 196.8 on 83 degrees of freedom
Multiple R-squared: 0.6488, Adjusted R-squared: 0.6403
F-statistic: 76.67 on 2 and 83 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,] 0.9748959 0.05020810 0.02510405
[2,] 0.9524710 0.09505806 0.04752903
[3,] 0.9562782 0.08744356 0.04372178
[4,] 0.9216364 0.15672729 0.07836364
[5,] 0.8912587 0.21748257 0.10874128
[6,] 0.8447777 0.31044462 0.15522231
[7,] 0.7891702 0.42165956 0.21082978
[8,] 0.7125939 0.57481211 0.28740606
[9,] 0.6790275 0.64194491 0.32097246
[10,] 0.6095335 0.78093306 0.39046653
[11,] 0.6938552 0.61228956 0.30614478
[12,] 0.6180209 0.76395826 0.38197913
[13,] 0.6171712 0.76565760 0.38282880
[14,] 0.5831204 0.83375921 0.41687961
[15,] 0.5391737 0.92165260 0.46082630
[16,] 0.4615116 0.92302313 0.53848843
[17,] 0.3898781 0.77975624 0.61012188
[18,] 0.3385372 0.67707443 0.66146279
[19,] 0.4666716 0.93334326 0.53332837
[20,] 0.3984326 0.79686511 0.60156744
[21,] 0.4113348 0.82266960 0.58866520
[22,] 0.3560183 0.71203661 0.64398170
[23,] 0.4045901 0.80918020 0.59540990
[24,] 0.3406284 0.68125687 0.65937156
[25,] 0.4482382 0.89647646 0.55176177
[26,] 0.4187098 0.83741961 0.58129019
[27,] 0.3813203 0.76264052 0.61867974
[28,] 0.3216405 0.64328100 0.67835950
[29,] 0.5496829 0.90063424 0.45031712
[30,] 0.5286752 0.94264958 0.47132479
[31,] 0.6258384 0.74832312 0.37416156
[32,] 0.5634106 0.87317878 0.43658939
[33,] 0.6485211 0.70295771 0.35147886
[34,] 0.6265158 0.74696838 0.37348419
[35,] 0.6110756 0.77784885 0.38892442
[36,] 0.5510558 0.89788834 0.44894417
[37,] 0.5515321 0.89693585 0.44846793
[38,] 0.5015783 0.99684348 0.49842174
[39,] 0.7325555 0.53488893 0.26744447
[40,] 0.6777997 0.64440059 0.32220030
[41,] 0.6427147 0.71457050 0.35728525
[42,] 0.6067088 0.78658241 0.39329121
[43,] 0.5797497 0.84050061 0.42025031
[44,] 0.5170990 0.96580207 0.48290104
[45,] 0.6061975 0.78760506 0.39380253
[46,] 0.6147813 0.77043732 0.38521866
[47,] 0.7688279 0.46234428 0.23117214
[48,] 0.7155881 0.56882377 0.28441188
[49,] 0.7220661 0.55586774 0.27793387
[50,] 0.6921395 0.61572099 0.30786049
[51,] 0.6505881 0.69882373 0.34941187
[52,] 0.5867668 0.82646633 0.41323317
[53,] 0.7711855 0.45762904 0.22881452
[54,] 0.7258024 0.54839518 0.27419759
[55,] 0.6722262 0.65554758 0.32777379
[56,] 0.6051001 0.78979987 0.39489994
[57,] 0.7208365 0.55832702 0.27916351
[58,] 0.6674853 0.66502941 0.33251470
[59,] 0.6118803 0.77623936 0.38811968
[60,] 0.5382082 0.92358357 0.46179178
[61,] 0.7048920 0.59021609 0.29510804
[62,] 0.6688041 0.66239181 0.33119591
[63,] 0.6441435 0.71171308 0.35585654
[64,] 0.5735977 0.85280450 0.42640225
[65,] 0.5084135 0.98317304 0.49158652
[66,] 0.4391541 0.87830830 0.56084585
[67,] 0.4574011 0.91480229 0.54259885
[68,] 0.3760103 0.75202059 0.62398970
[69,] 0.3684444 0.73688879 0.63155560
[70,] 0.2863141 0.57262813 0.71368594
[71,] 0.3065887 0.61317741 0.69341130
[72,] 0.2130227 0.42604534 0.78697733
[73,] 0.2428367 0.48567335 0.75716333
[74,] 0.1506802 0.30136035 0.84931983
[75,] 0.3917954 0.78359085 0.60820458
> postscript(file="/var/www/html/rcomp/tmp/1pkm51292787493.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/2pkm51292787493.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/30bm81292787493.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/40bm81292787493.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/50bm81292787493.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 = 86
Frequency = 1
1 2 3 4 5 6
-355.2418085 -356.0933521 148.3874969 -254.8209287 9.3277899 271.1792610
7 8 9 10 11 12
98.2486956 311.5549910 -21.1697532 5.3076741 11.7660621 -19.7299687
13 14 15 16 17 18
-23.3493504 96.8302817 58.6062624 -284.8914347 -40.6665080 177.1709076
19 20 21 22 23 24
146.6879609 -136.6096944 -2.5535858 -48.9570459 99.5990516 -347.3707815
25 26 27 28 29 30
-46.5871319 218.9265443 74.9038645 -284.2057225 -33.9498068 337.4434002
31 32 33 34 35 36
130.5853807 -149.8772874 -42.5152171 443.8099249 171.2038975 -324.6256202
37 38 39 40 41 42
2.4282324 308.6739084 150.3953222 -165.4222658 -45.4945595 -159.7099116
43 44 45 46 47 48
73.9944552 -414.9561567 -28.9331299 79.1930757 122.1482277 86.7183024
49 50 51 52 53 54
-36.8207746 -305.2355806 215.2285612 -334.1133129 1.4976714 -217.5030312
55 56 57 58 59 60
164.3860638 108.7935052 -35.7018668 -406.9547749 81.4616153 68.9937074
61 62 63 64 65 66
-2.6614555 -332.7871080 98.1812156 -0.8837738 -5.6350289 -343.1993478
67 68 69 70 71 72
203.1449318 -80.8131611 -16.5597295 117.1902643 122.0639136 302.7013206
73 74 75 76 77 78
-29.7672291 210.0595398 82.6278774 275.3688993 -3.7439035 -203.8240504
79 80 81 82 83 84
75.7740720 -325.0842714 -28.4685368 379.3835800 139.2867796 99.1226483
85 86
13.2557839 -96.1249629
> postscript(file="/var/www/html/rcomp/tmp/6t23t1292787493.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 = 86
Frequency = 1
lag(myerror, k = 1) myerror
0 -355.2418085 NA
1 -356.0933521 -355.2418085
2 148.3874969 -356.0933521
3 -254.8209287 148.3874969
4 9.3277899 -254.8209287
5 271.1792610 9.3277899
6 98.2486956 271.1792610
7 311.5549910 98.2486956
8 -21.1697532 311.5549910
9 5.3076741 -21.1697532
10 11.7660621 5.3076741
11 -19.7299687 11.7660621
12 -23.3493504 -19.7299687
13 96.8302817 -23.3493504
14 58.6062624 96.8302817
15 -284.8914347 58.6062624
16 -40.6665080 -284.8914347
17 177.1709076 -40.6665080
18 146.6879609 177.1709076
19 -136.6096944 146.6879609
20 -2.5535858 -136.6096944
21 -48.9570459 -2.5535858
22 99.5990516 -48.9570459
23 -347.3707815 99.5990516
24 -46.5871319 -347.3707815
25 218.9265443 -46.5871319
26 74.9038645 218.9265443
27 -284.2057225 74.9038645
28 -33.9498068 -284.2057225
29 337.4434002 -33.9498068
30 130.5853807 337.4434002
31 -149.8772874 130.5853807
32 -42.5152171 -149.8772874
33 443.8099249 -42.5152171
34 171.2038975 443.8099249
35 -324.6256202 171.2038975
36 2.4282324 -324.6256202
37 308.6739084 2.4282324
38 150.3953222 308.6739084
39 -165.4222658 150.3953222
40 -45.4945595 -165.4222658
41 -159.7099116 -45.4945595
42 73.9944552 -159.7099116
43 -414.9561567 73.9944552
44 -28.9331299 -414.9561567
45 79.1930757 -28.9331299
46 122.1482277 79.1930757
47 86.7183024 122.1482277
48 -36.8207746 86.7183024
49 -305.2355806 -36.8207746
50 215.2285612 -305.2355806
51 -334.1133129 215.2285612
52 1.4976714 -334.1133129
53 -217.5030312 1.4976714
54 164.3860638 -217.5030312
55 108.7935052 164.3860638
56 -35.7018668 108.7935052
57 -406.9547749 -35.7018668
58 81.4616153 -406.9547749
59 68.9937074 81.4616153
60 -2.6614555 68.9937074
61 -332.7871080 -2.6614555
62 98.1812156 -332.7871080
63 -0.8837738 98.1812156
64 -5.6350289 -0.8837738
65 -343.1993478 -5.6350289
66 203.1449318 -343.1993478
67 -80.8131611 203.1449318
68 -16.5597295 -80.8131611
69 117.1902643 -16.5597295
70 122.0639136 117.1902643
71 302.7013206 122.0639136
72 -29.7672291 302.7013206
73 210.0595398 -29.7672291
74 82.6278774 210.0595398
75 275.3688993 82.6278774
76 -3.7439035 275.3688993
77 -203.8240504 -3.7439035
78 75.7740720 -203.8240504
79 -325.0842714 75.7740720
80 -28.4685368 -325.0842714
81 379.3835800 -28.4685368
82 139.2867796 379.3835800
83 99.1226483 139.2867796
84 13.2557839 99.1226483
85 -96.1249629 13.2557839
86 NA -96.1249629
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -356.0933521 -355.2418085
[2,] 148.3874969 -356.0933521
[3,] -254.8209287 148.3874969
[4,] 9.3277899 -254.8209287
[5,] 271.1792610 9.3277899
[6,] 98.2486956 271.1792610
[7,] 311.5549910 98.2486956
[8,] -21.1697532 311.5549910
[9,] 5.3076741 -21.1697532
[10,] 11.7660621 5.3076741
[11,] -19.7299687 11.7660621
[12,] -23.3493504 -19.7299687
[13,] 96.8302817 -23.3493504
[14,] 58.6062624 96.8302817
[15,] -284.8914347 58.6062624
[16,] -40.6665080 -284.8914347
[17,] 177.1709076 -40.6665080
[18,] 146.6879609 177.1709076
[19,] -136.6096944 146.6879609
[20,] -2.5535858 -136.6096944
[21,] -48.9570459 -2.5535858
[22,] 99.5990516 -48.9570459
[23,] -347.3707815 99.5990516
[24,] -46.5871319 -347.3707815
[25,] 218.9265443 -46.5871319
[26,] 74.9038645 218.9265443
[27,] -284.2057225 74.9038645
[28,] -33.9498068 -284.2057225
[29,] 337.4434002 -33.9498068
[30,] 130.5853807 337.4434002
[31,] -149.8772874 130.5853807
[32,] -42.5152171 -149.8772874
[33,] 443.8099249 -42.5152171
[34,] 171.2038975 443.8099249
[35,] -324.6256202 171.2038975
[36,] 2.4282324 -324.6256202
[37,] 308.6739084 2.4282324
[38,] 150.3953222 308.6739084
[39,] -165.4222658 150.3953222
[40,] -45.4945595 -165.4222658
[41,] -159.7099116 -45.4945595
[42,] 73.9944552 -159.7099116
[43,] -414.9561567 73.9944552
[44,] -28.9331299 -414.9561567
[45,] 79.1930757 -28.9331299
[46,] 122.1482277 79.1930757
[47,] 86.7183024 122.1482277
[48,] -36.8207746 86.7183024
[49,] -305.2355806 -36.8207746
[50,] 215.2285612 -305.2355806
[51,] -334.1133129 215.2285612
[52,] 1.4976714 -334.1133129
[53,] -217.5030312 1.4976714
[54,] 164.3860638 -217.5030312
[55,] 108.7935052 164.3860638
[56,] -35.7018668 108.7935052
[57,] -406.9547749 -35.7018668
[58,] 81.4616153 -406.9547749
[59,] 68.9937074 81.4616153
[60,] -2.6614555 68.9937074
[61,] -332.7871080 -2.6614555
[62,] 98.1812156 -332.7871080
[63,] -0.8837738 98.1812156
[64,] -5.6350289 -0.8837738
[65,] -343.1993478 -5.6350289
[66,] 203.1449318 -343.1993478
[67,] -80.8131611 203.1449318
[68,] -16.5597295 -80.8131611
[69,] 117.1902643 -16.5597295
[70,] 122.0639136 117.1902643
[71,] 302.7013206 122.0639136
[72,] -29.7672291 302.7013206
[73,] 210.0595398 -29.7672291
[74,] 82.6278774 210.0595398
[75,] 275.3688993 82.6278774
[76,] -3.7439035 275.3688993
[77,] -203.8240504 -3.7439035
[78,] 75.7740720 -203.8240504
[79,] -325.0842714 75.7740720
[80,] -28.4685368 -325.0842714
[81,] 379.3835800 -28.4685368
[82,] 139.2867796 379.3835800
[83,] 99.1226483 139.2867796
[84,] 13.2557839 99.1226483
[85,] -96.1249629 13.2557839
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -356.0933521 -355.2418085
2 148.3874969 -356.0933521
3 -254.8209287 148.3874969
4 9.3277899 -254.8209287
5 271.1792610 9.3277899
6 98.2486956 271.1792610
7 311.5549910 98.2486956
8 -21.1697532 311.5549910
9 5.3076741 -21.1697532
10 11.7660621 5.3076741
11 -19.7299687 11.7660621
12 -23.3493504 -19.7299687
13 96.8302817 -23.3493504
14 58.6062624 96.8302817
15 -284.8914347 58.6062624
16 -40.6665080 -284.8914347
17 177.1709076 -40.6665080
18 146.6879609 177.1709076
19 -136.6096944 146.6879609
20 -2.5535858 -136.6096944
21 -48.9570459 -2.5535858
22 99.5990516 -48.9570459
23 -347.3707815 99.5990516
24 -46.5871319 -347.3707815
25 218.9265443 -46.5871319
26 74.9038645 218.9265443
27 -284.2057225 74.9038645
28 -33.9498068 -284.2057225
29 337.4434002 -33.9498068
30 130.5853807 337.4434002
31 -149.8772874 130.5853807
32 -42.5152171 -149.8772874
33 443.8099249 -42.5152171
34 171.2038975 443.8099249
35 -324.6256202 171.2038975
36 2.4282324 -324.6256202
37 308.6739084 2.4282324
38 150.3953222 308.6739084
39 -165.4222658 150.3953222
40 -45.4945595 -165.4222658
41 -159.7099116 -45.4945595
42 73.9944552 -159.7099116
43 -414.9561567 73.9944552
44 -28.9331299 -414.9561567
45 79.1930757 -28.9331299
46 122.1482277 79.1930757
47 86.7183024 122.1482277
48 -36.8207746 86.7183024
49 -305.2355806 -36.8207746
50 215.2285612 -305.2355806
51 -334.1133129 215.2285612
52 1.4976714 -334.1133129
53 -217.5030312 1.4976714
54 164.3860638 -217.5030312
55 108.7935052 164.3860638
56 -35.7018668 108.7935052
57 -406.9547749 -35.7018668
58 81.4616153 -406.9547749
59 68.9937074 81.4616153
60 -2.6614555 68.9937074
61 -332.7871080 -2.6614555
62 98.1812156 -332.7871080
63 -0.8837738 98.1812156
64 -5.6350289 -0.8837738
65 -343.1993478 -5.6350289
66 203.1449318 -343.1993478
67 -80.8131611 203.1449318
68 -16.5597295 -80.8131611
69 117.1902643 -16.5597295
70 122.0639136 117.1902643
71 302.7013206 122.0639136
72 -29.7672291 302.7013206
73 210.0595398 -29.7672291
74 82.6278774 210.0595398
75 275.3688993 82.6278774
76 -3.7439035 275.3688993
77 -203.8240504 -3.7439035
78 75.7740720 -203.8240504
79 -325.0842714 75.7740720
80 -28.4685368 -325.0842714
81 379.3835800 -28.4685368
82 139.2867796 379.3835800
83 99.1226483 139.2867796
84 13.2557839 99.1226483
85 -96.1249629 13.2557839
> 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/73bkw1292787493.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/83bkw1292787493.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/93bkw1292787493.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/10elkh1292787493.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/110lin1292787493.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/1234hb1292787493.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/13zvw11292787493.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/14r5d41292787493.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/15nff51292787494.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/161pce1292787494.tab")
+ }
>
> try(system("convert tmp/1pkm51292787493.ps tmp/1pkm51292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/2pkm51292787493.ps tmp/2pkm51292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/30bm81292787493.ps tmp/30bm81292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/40bm81292787493.ps tmp/40bm81292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/50bm81292787493.ps tmp/50bm81292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/6t23t1292787493.ps tmp/6t23t1292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/73bkw1292787493.ps tmp/73bkw1292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/83bkw1292787493.ps tmp/83bkw1292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/93bkw1292787493.ps tmp/93bkw1292787493.png",intern=TRUE))
character(0)
> try(system("convert tmp/10elkh1292787493.ps tmp/10elkh1292787493.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.825 1.668 6.859