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(1
+ ,216234.00
+ ,627
+ ,2
+ ,213586.00
+ ,696
+ ,3
+ ,209465.00
+ ,825
+ ,4
+ ,204045.00
+ ,677
+ ,5
+ ,200237.00
+ ,656
+ ,6
+ ,203666.00
+ ,785
+ ,7
+ ,241476.00
+ ,412
+ ,8
+ ,260307.00
+ ,352
+ ,9
+ ,243324.00
+ ,839
+ ,10
+ ,244460.00
+ ,729
+ ,11
+ ,233575.00
+ ,696
+ ,12
+ ,237217.00
+ ,641
+ ,1
+ ,235243.00
+ ,695
+ ,2
+ ,230354.00
+ ,638
+ ,3
+ ,227184.00
+ ,762
+ ,4
+ ,221678.00
+ ,635
+ ,5
+ ,217142.00
+ ,721
+ ,6
+ ,219452.00
+ ,854
+ ,7
+ ,256446.00
+ ,418
+ ,8
+ ,265845.00
+ ,367
+ ,9
+ ,248624.00
+ ,824
+ ,10
+ ,241114.00
+ ,687
+ ,11
+ ,229245.00
+ ,601
+ ,12
+ ,231805.00
+ ,676
+ ,1
+ ,219277.00
+ ,740
+ ,2
+ ,219313.00
+ ,691
+ ,3
+ ,212610.00
+ ,683
+ ,4
+ ,214771.00
+ ,594
+ ,5
+ ,211142.00
+ ,729
+ ,6
+ ,211457.00
+ ,731
+ ,7
+ ,240048.00
+ ,386
+ ,8
+ ,240636.00
+ ,331
+ ,9
+ ,230580.00
+ ,707
+ ,10
+ ,208795.00
+ ,715
+ ,11
+ ,197922.00
+ ,657
+ ,12
+ ,194596.00
+ ,653
+ ,1
+ ,194581.00
+ ,642
+ ,2
+ ,185686.00
+ ,643
+ ,3
+ ,178106.00
+ ,718
+ ,4
+ ,172608.00
+ ,654
+ ,5
+ ,167302.00
+ ,632
+ ,6
+ ,168053.00
+ ,731
+ ,7
+ ,202300.00
+ ,392
+ ,8
+ ,202388.00
+ ,344
+ ,9
+ ,182516.00
+ ,792
+ ,10
+ ,173476.00
+ ,852
+ ,11
+ ,166444.00
+ ,649
+ ,12
+ ,171297.00
+ ,629
+ ,1
+ ,169701.00
+ ,685
+ ,2
+ ,164182.00
+ ,617
+ ,3
+ ,161914.00
+ ,715
+ ,4
+ ,159612.00
+ ,715
+ ,5
+ ,151001.00
+ ,629
+ ,6
+ ,158114.00
+ ,916
+ ,7
+ ,186530.00
+ ,531
+ ,8
+ ,187069.00
+ ,357
+ ,9
+ ,174330.00
+ ,917
+ ,10
+ ,169362.00
+ ,828
+ ,11
+ ,166827.00
+ ,708
+ ,12
+ ,178037.00
+ ,858
+ ,1
+ ,186413.00
+ ,775
+ ,2
+ ,189226.00
+ ,785
+ ,3
+ ,191563.00
+ ,1006
+ ,4
+ ,188906.00
+ ,789
+ ,5
+ ,186005.00
+ ,734
+ ,6
+ ,195309.00
+ ,906
+ ,7
+ ,223532.00
+ ,532
+ ,8
+ ,226899.00
+ ,387
+ ,9
+ ,214126.00
+ ,991
+ ,10
+ ,206903.00
+ ,841
+ ,11
+ ,204442.00
+ ,892
+ ,12
+ ,220375.00
+ ,782)
+ ,dim=c(3
+ ,72)
+ ,dimnames=list(c('month'
+ ,'werklozen'
+ ,'faillissementen')
+ ,1:72))
> y <- array(NA,dim=c(3,72),dimnames=list(c('month','werklozen','faillissementen'),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 = 'Linear Trend'
> par2 = 'Do not include Seasonal Dummies'
> par1 = '2'
> #'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 month faillissementen t
1 216234 1 627 1
2 213586 2 696 2
3 209465 3 825 3
4 204045 4 677 4
5 200237 5 656 5
6 203666 6 785 6
7 241476 7 412 7
8 260307 8 352 8
9 243324 9 839 9
10 244460 10 729 10
11 233575 11 696 11
12 237217 12 641 12
13 235243 1 695 13
14 230354 2 638 14
15 227184 3 762 15
16 221678 4 635 16
17 217142 5 721 17
18 219452 6 854 18
19 256446 7 418 19
20 265845 8 367 20
21 248624 9 824 21
22 241114 10 687 22
23 229245 11 601 23
24 231805 12 676 24
25 219277 1 740 25
26 219313 2 691 26
27 212610 3 683 27
28 214771 4 594 28
29 211142 5 729 29
30 211457 6 731 30
31 240048 7 386 31
32 240636 8 331 32
33 230580 9 707 33
34 208795 10 715 34
35 197922 11 657 35
36 194596 12 653 36
37 194581 1 642 37
38 185686 2 643 38
39 178106 3 718 39
40 172608 4 654 40
41 167302 5 632 41
42 168053 6 731 42
43 202300 7 392 43
44 202388 8 344 44
45 182516 9 792 45
46 173476 10 852 46
47 166444 11 649 47
48 171297 12 629 48
49 169701 1 685 49
50 164182 2 617 50
51 161914 3 715 51
52 159612 4 715 52
53 151001 5 629 53
54 158114 6 916 54
55 186530 7 531 55
56 187069 8 357 56
57 174330 9 917 57
58 169362 10 828 58
59 166827 11 708 59
60 178037 12 858 60
61 186413 1 775 61
62 189226 2 785 62
63 191563 3 1006 63
64 188906 4 789 64
65 186005 5 734 65
66 195309 6 906 66
67 223532 7 532 67
68 226899 8 387 68
69 214126 9 991 69
70 206903 10 841 70
71 204442 11 892 71
72 220375 12 782 72
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) month faillissementen t
243874.21 1987.33 -34.91 -768.69
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-40504 -17387 1872 17207 40005
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 243874.21 12484.71 19.534 < 2e-16 ***
month 1987.33 746.82 2.661 0.0097 **
faillissementen -34.91 16.83 -2.074 0.0418 *
t -768.69 127.64 -6.023 7.77e-08 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 21530 on 68 degrees of freedom
Multiple R-squared: 0.4354, Adjusted R-squared: 0.4105
F-statistic: 17.48 on 3 and 68 DF, p-value: 1.615e-08
> 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.162896978 0.325793957 0.8371030217
[2,] 0.135308886 0.270617772 0.8646911140
[3,] 0.256452868 0.512905736 0.7435471320
[4,] 0.160303491 0.320606983 0.8396965087
[5,] 0.112893301 0.225786601 0.8871066994
[6,] 0.077440796 0.154881593 0.9225592036
[7,] 0.043777282 0.087554564 0.9562227179
[8,] 0.027029217 0.054058434 0.9729707832
[9,] 0.014660553 0.029321106 0.9853394468
[10,] 0.014334246 0.028668491 0.9856657544
[11,] 0.012162238 0.024324476 0.9878377620
[12,] 0.006658620 0.013317241 0.9933413796
[13,] 0.003989643 0.007979285 0.9960103573
[14,] 0.002949590 0.005899181 0.9970504095
[15,] 0.002826383 0.005652765 0.9971736175
[16,] 0.002126818 0.004253636 0.9978731822
[17,] 0.004927568 0.009855137 0.9950724315
[18,] 0.005238360 0.010476719 0.9947616403
[19,] 0.003949817 0.007899633 0.9960501833
[20,] 0.003366278 0.006732556 0.9966337219
[21,] 0.003677440 0.007354880 0.9963225601
[22,] 0.004628292 0.009256584 0.9953717080
[23,] 0.004814187 0.009628374 0.9951858129
[24,] 0.005326338 0.010652675 0.9946736623
[25,] 0.007545786 0.015091573 0.9924542136
[26,] 0.014691383 0.029382765 0.9853086174
[27,] 0.041719461 0.083438922 0.9582805391
[28,] 0.099825380 0.199650760 0.9001746201
[29,] 0.246652925 0.493305850 0.7533470748
[30,] 0.437998269 0.875996537 0.5620017313
[31,] 0.548613773 0.902772454 0.4513862268
[32,] 0.644024271 0.711951458 0.3559757292
[33,] 0.718915674 0.562168652 0.2810843260
[34,] 0.777207521 0.445584957 0.2227924786
[35,] 0.821201174 0.357597652 0.1787988262
[36,] 0.829164560 0.341670880 0.1708354401
[37,] 0.896568030 0.206863940 0.1034319698
[38,] 0.961439384 0.077121233 0.0385606163
[39,] 0.988093433 0.023813134 0.0119065671
[40,] 0.996691135 0.006617729 0.0033088647
[41,] 0.997790986 0.004418028 0.0022090140
[42,] 0.999248004 0.001503992 0.0007519960
[43,] 0.999456735 0.001086530 0.0005432652
[44,] 0.999071154 0.001857692 0.0009288461
[45,] 0.998379059 0.003241882 0.0016209408
[46,] 0.996760029 0.006479941 0.0032399706
[47,] 0.998169952 0.003660096 0.0018300482
[48,] 0.996115278 0.007769444 0.0038847218
[49,] 0.995630460 0.008739079 0.0043695397
[50,] 0.991877684 0.016244631 0.0081223157
[51,] 0.992251938 0.015496124 0.0077480622
[52,] 0.984285059 0.031429883 0.0157149413
[53,] 0.980477204 0.039045593 0.0195227963
[54,] 0.961392842 0.077214316 0.0386071579
[55,] 0.941772002 0.116455997 0.0582279984
[56,] 0.907082505 0.185834991 0.0929174955
[57,] 0.909114726 0.181770548 0.0908852738
[58,] 0.831710602 0.336578795 0.1682893977
[59,] 0.866746971 0.266506057 0.1332530286
> postscript(file="/var/www/html/rcomp/tmp/1w8jo1292614124.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/27z081292614124.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/37z081292614124.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/47z081292614124.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/5i9ic1292614124.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
-6968.3672 -8426.0037 -9261.8569 -21067.6260 -26827.4373 -20113.2905
7 8 9 10 11 12
3455.5032 18973.0828 17774.1029 13851.0298 595.2618 1098.4067
13 14 15 16 17 18
23639.0558 15541.3746 15481.9561 4323.3611 1571.2466 7306.0456
19 20 21 22 23 24
27859.3169 34259.1140 31774.7424 18263.0169 2172.8571 6132.6990
25 26 27 28 29 30
18468.4786 15575.1018 7374.1602 5209.2613 5074.8864 4241.0753
31 32 33 34 35 36
19568.4346 17017.5795 18870.2505 -3854.0822 -17970.6765 -22654.9659
37 38 39 40 41 42
-424.6653 -10503.3894 -16683.5475 -25634.6201 -32927.3445 -29938.5893
43 44 45 46 47 48
-8745.7517 -11552.2154 -17001.8045 -25165.6584 -40503.6455 -37567.5438
49 50 51 52 53 54
-14579.0686 -23690.7934 -23755.9513 -27276.5884 -40108.7482 -24194.3390
55 56 57 58 59 60
-10438.5018 -17193.0104 -11599.3374 -20893.2364 -28836.4400 -13608.1191
61 62 63 64 65 66
14499.4417 16442.9350 25277.0828 13825.3129 7785.4578 21875.8659
67 68 69 70 71 72
35822.7466 32908.7166 40004.5640 26325.9687 24426.8973 35300.8242
> postscript(file="/var/www/html/rcomp/tmp/6i9ic1292614124.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 -6968.3672 NA
1 -8426.0037 -6968.3672
2 -9261.8569 -8426.0037
3 -21067.6260 -9261.8569
4 -26827.4373 -21067.6260
5 -20113.2905 -26827.4373
6 3455.5032 -20113.2905
7 18973.0828 3455.5032
8 17774.1029 18973.0828
9 13851.0298 17774.1029
10 595.2618 13851.0298
11 1098.4067 595.2618
12 23639.0558 1098.4067
13 15541.3746 23639.0558
14 15481.9561 15541.3746
15 4323.3611 15481.9561
16 1571.2466 4323.3611
17 7306.0456 1571.2466
18 27859.3169 7306.0456
19 34259.1140 27859.3169
20 31774.7424 34259.1140
21 18263.0169 31774.7424
22 2172.8571 18263.0169
23 6132.6990 2172.8571
24 18468.4786 6132.6990
25 15575.1018 18468.4786
26 7374.1602 15575.1018
27 5209.2613 7374.1602
28 5074.8864 5209.2613
29 4241.0753 5074.8864
30 19568.4346 4241.0753
31 17017.5795 19568.4346
32 18870.2505 17017.5795
33 -3854.0822 18870.2505
34 -17970.6765 -3854.0822
35 -22654.9659 -17970.6765
36 -424.6653 -22654.9659
37 -10503.3894 -424.6653
38 -16683.5475 -10503.3894
39 -25634.6201 -16683.5475
40 -32927.3445 -25634.6201
41 -29938.5893 -32927.3445
42 -8745.7517 -29938.5893
43 -11552.2154 -8745.7517
44 -17001.8045 -11552.2154
45 -25165.6584 -17001.8045
46 -40503.6455 -25165.6584
47 -37567.5438 -40503.6455
48 -14579.0686 -37567.5438
49 -23690.7934 -14579.0686
50 -23755.9513 -23690.7934
51 -27276.5884 -23755.9513
52 -40108.7482 -27276.5884
53 -24194.3390 -40108.7482
54 -10438.5018 -24194.3390
55 -17193.0104 -10438.5018
56 -11599.3374 -17193.0104
57 -20893.2364 -11599.3374
58 -28836.4400 -20893.2364
59 -13608.1191 -28836.4400
60 14499.4417 -13608.1191
61 16442.9350 14499.4417
62 25277.0828 16442.9350
63 13825.3129 25277.0828
64 7785.4578 13825.3129
65 21875.8659 7785.4578
66 35822.7466 21875.8659
67 32908.7166 35822.7466
68 40004.5640 32908.7166
69 26325.9687 40004.5640
70 24426.8973 26325.9687
71 35300.8242 24426.8973
72 NA 35300.8242
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -8426.0037 -6968.3672
[2,] -9261.8569 -8426.0037
[3,] -21067.6260 -9261.8569
[4,] -26827.4373 -21067.6260
[5,] -20113.2905 -26827.4373
[6,] 3455.5032 -20113.2905
[7,] 18973.0828 3455.5032
[8,] 17774.1029 18973.0828
[9,] 13851.0298 17774.1029
[10,] 595.2618 13851.0298
[11,] 1098.4067 595.2618
[12,] 23639.0558 1098.4067
[13,] 15541.3746 23639.0558
[14,] 15481.9561 15541.3746
[15,] 4323.3611 15481.9561
[16,] 1571.2466 4323.3611
[17,] 7306.0456 1571.2466
[18,] 27859.3169 7306.0456
[19,] 34259.1140 27859.3169
[20,] 31774.7424 34259.1140
[21,] 18263.0169 31774.7424
[22,] 2172.8571 18263.0169
[23,] 6132.6990 2172.8571
[24,] 18468.4786 6132.6990
[25,] 15575.1018 18468.4786
[26,] 7374.1602 15575.1018
[27,] 5209.2613 7374.1602
[28,] 5074.8864 5209.2613
[29,] 4241.0753 5074.8864
[30,] 19568.4346 4241.0753
[31,] 17017.5795 19568.4346
[32,] 18870.2505 17017.5795
[33,] -3854.0822 18870.2505
[34,] -17970.6765 -3854.0822
[35,] -22654.9659 -17970.6765
[36,] -424.6653 -22654.9659
[37,] -10503.3894 -424.6653
[38,] -16683.5475 -10503.3894
[39,] -25634.6201 -16683.5475
[40,] -32927.3445 -25634.6201
[41,] -29938.5893 -32927.3445
[42,] -8745.7517 -29938.5893
[43,] -11552.2154 -8745.7517
[44,] -17001.8045 -11552.2154
[45,] -25165.6584 -17001.8045
[46,] -40503.6455 -25165.6584
[47,] -37567.5438 -40503.6455
[48,] -14579.0686 -37567.5438
[49,] -23690.7934 -14579.0686
[50,] -23755.9513 -23690.7934
[51,] -27276.5884 -23755.9513
[52,] -40108.7482 -27276.5884
[53,] -24194.3390 -40108.7482
[54,] -10438.5018 -24194.3390
[55,] -17193.0104 -10438.5018
[56,] -11599.3374 -17193.0104
[57,] -20893.2364 -11599.3374
[58,] -28836.4400 -20893.2364
[59,] -13608.1191 -28836.4400
[60,] 14499.4417 -13608.1191
[61,] 16442.9350 14499.4417
[62,] 25277.0828 16442.9350
[63,] 13825.3129 25277.0828
[64,] 7785.4578 13825.3129
[65,] 21875.8659 7785.4578
[66,] 35822.7466 21875.8659
[67,] 32908.7166 35822.7466
[68,] 40004.5640 32908.7166
[69,] 26325.9687 40004.5640
[70,] 24426.8973 26325.9687
[71,] 35300.8242 24426.8973
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -8426.0037 -6968.3672
2 -9261.8569 -8426.0037
3 -21067.6260 -9261.8569
4 -26827.4373 -21067.6260
5 -20113.2905 -26827.4373
6 3455.5032 -20113.2905
7 18973.0828 3455.5032
8 17774.1029 18973.0828
9 13851.0298 17774.1029
10 595.2618 13851.0298
11 1098.4067 595.2618
12 23639.0558 1098.4067
13 15541.3746 23639.0558
14 15481.9561 15541.3746
15 4323.3611 15481.9561
16 1571.2466 4323.3611
17 7306.0456 1571.2466
18 27859.3169 7306.0456
19 34259.1140 27859.3169
20 31774.7424 34259.1140
21 18263.0169 31774.7424
22 2172.8571 18263.0169
23 6132.6990 2172.8571
24 18468.4786 6132.6990
25 15575.1018 18468.4786
26 7374.1602 15575.1018
27 5209.2613 7374.1602
28 5074.8864 5209.2613
29 4241.0753 5074.8864
30 19568.4346 4241.0753
31 17017.5795 19568.4346
32 18870.2505 17017.5795
33 -3854.0822 18870.2505
34 -17970.6765 -3854.0822
35 -22654.9659 -17970.6765
36 -424.6653 -22654.9659
37 -10503.3894 -424.6653
38 -16683.5475 -10503.3894
39 -25634.6201 -16683.5475
40 -32927.3445 -25634.6201
41 -29938.5893 -32927.3445
42 -8745.7517 -29938.5893
43 -11552.2154 -8745.7517
44 -17001.8045 -11552.2154
45 -25165.6584 -17001.8045
46 -40503.6455 -25165.6584
47 -37567.5438 -40503.6455
48 -14579.0686 -37567.5438
49 -23690.7934 -14579.0686
50 -23755.9513 -23690.7934
51 -27276.5884 -23755.9513
52 -40108.7482 -27276.5884
53 -24194.3390 -40108.7482
54 -10438.5018 -24194.3390
55 -17193.0104 -10438.5018
56 -11599.3374 -17193.0104
57 -20893.2364 -11599.3374
58 -28836.4400 -20893.2364
59 -13608.1191 -28836.4400
60 14499.4417 -13608.1191
61 16442.9350 14499.4417
62 25277.0828 16442.9350
63 13825.3129 25277.0828
64 7785.4578 13825.3129
65 21875.8659 7785.4578
66 35822.7466 21875.8659
67 32908.7166 35822.7466
68 40004.5640 32908.7166
69 26325.9687 40004.5640
70 24426.8973 26325.9687
71 35300.8242 24426.8973
> 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/7bizf1292614124.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/8bizf1292614124.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/939yh1292614124.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/1039yh1292614124.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/116re51292614124.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/12asvb1292614124.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/13hts51292614124.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/14al9q1292614124.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/15v3qw1292614124.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/16rv6m1292614124.tab")
+ }
> try(system("convert tmp/1w8jo1292614124.ps tmp/1w8jo1292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/27z081292614124.ps tmp/27z081292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/37z081292614124.ps tmp/37z081292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/47z081292614124.ps tmp/47z081292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/5i9ic1292614124.ps tmp/5i9ic1292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/6i9ic1292614124.ps tmp/6i9ic1292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/7bizf1292614124.ps tmp/7bizf1292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/8bizf1292614124.ps tmp/8bizf1292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/939yh1292614124.ps tmp/939yh1292614124.png",intern=TRUE))
character(0)
> try(system("convert tmp/1039yh1292614124.ps tmp/1039yh1292614124.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.606 1.718 6.623