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(593530
+ ,3922
+ ,18004
+ ,707169
+ ,610763
+ ,3759
+ ,17537
+ ,703434
+ ,612613
+ ,4138
+ ,20366
+ ,701017
+ ,611324
+ ,4634
+ ,22782
+ ,696968
+ ,594167
+ ,3996
+ ,19169
+ ,688558
+ ,595454
+ ,4308
+ ,13807
+ ,679237
+ ,590865
+ ,4143
+ ,29743
+ ,677362
+ ,589379
+ ,4429
+ ,25591
+ ,676693
+ ,584428
+ ,5219
+ ,29096
+ ,670009
+ ,573100
+ ,4929
+ ,26482
+ ,667209
+ ,567456
+ ,5761
+ ,22405
+ ,662976
+ ,569028
+ ,5592
+ ,27044
+ ,660194
+ ,620735
+ ,4163
+ ,17970
+ ,652270
+ ,628884
+ ,4962
+ ,18730
+ ,648024
+ ,628232
+ ,5208
+ ,19684
+ ,629295
+ ,612117
+ ,4755
+ ,19785
+ ,624961
+ ,595404
+ ,4491
+ ,18479
+ ,617306
+ ,597141
+ ,5732
+ ,10698
+ ,607691
+ ,593408
+ ,5731
+ ,31956
+ ,596219
+ ,590072
+ ,5040
+ ,29506
+ ,591130
+ ,579799
+ ,6102
+ ,34506
+ ,584528
+ ,574205
+ ,4904
+ ,27165
+ ,576798
+ ,572775
+ ,5369
+ ,26736
+ ,575683
+ ,572942
+ ,5578
+ ,23691
+ ,574369
+ ,619567
+ ,4619
+ ,18157
+ ,566815
+ ,625809
+ ,4731
+ ,17328
+ ,573074
+ ,619916
+ ,5011
+ ,18205
+ ,567739
+ ,587625
+ ,5299
+ ,20995
+ ,571942
+ ,565742
+ ,4146
+ ,17382
+ ,570274
+ ,557274
+ ,4625
+ ,9367
+ ,568800
+ ,560576
+ ,4736
+ ,31124
+ ,558115
+ ,548854
+ ,4219
+ ,26551
+ ,550591
+ ,531673
+ ,5116
+ ,30651
+ ,548872
+ ,525919
+ ,4205
+ ,25859
+ ,547009
+ ,511038
+ ,4121
+ ,25100
+ ,545946
+ ,498662
+ ,5103
+ ,25778
+ ,539702
+ ,555362
+ ,4300
+ ,20418
+ ,542427
+ ,564591
+ ,4578
+ ,18688
+ ,542968
+ ,541657
+ ,3809
+ ,20424
+ ,536640
+ ,527070
+ ,5526
+ ,24776
+ ,533653
+ ,509846
+ ,4248
+ ,19814
+ ,540996
+ ,514258
+ ,3830
+ ,12738
+ ,538316
+ ,516922
+ ,4428
+ ,31566
+ ,532646
+ ,507561
+ ,4834
+ ,30111
+ ,533390
+ ,492622
+ ,4406
+ ,30019
+ ,528715
+ ,490243
+ ,4565
+ ,31934
+ ,530664
+ ,469357
+ ,4104
+ ,25826
+ ,528564
+ ,477580
+ ,4798
+ ,26835
+ ,519107
+ ,528379
+ ,3935
+ ,20205
+ ,518703
+ ,533590
+ ,3792
+ ,17789
+ ,519059
+ ,517945
+ ,4387
+ ,20520
+ ,518498
+ ,506174
+ ,4006
+ ,22518
+ ,524575
+ ,501866
+ ,4078
+ ,15572
+ ,536046
+ ,516141
+ ,4724
+ ,11509
+ ,552006
+ ,528222
+ ,3157
+ ,25447
+ ,560687
+ ,532638
+ ,3558
+ ,24090
+ ,578884
+ ,536322
+ ,3899
+ ,27786
+ ,591491
+ ,536535
+ ,4118
+ ,26195
+ ,599228
+ ,523597
+ ,3790
+ ,20516
+ ,633019
+ ,536214
+ ,4278
+ ,22759
+ ,649918
+ ,586570
+ ,4035
+ ,19028
+ ,655509)
+ ,dim=c(4
+ ,61)
+ ,dimnames=list(c('Werkzoekend'
+ ,'Bouw'
+ ,'Auto'
+ ,'Krediet')
+ ,1:61))
> y <- array(NA,dim=c(4,61),dimnames=list(c('Werkzoekend','Bouw','Auto','Krediet'),1:61))
> 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
Werkzoekend Bouw Auto Krediet t
1 593530 3922 18004 707169 1
2 610763 3759 17537 703434 2
3 612613 4138 20366 701017 3
4 611324 4634 22782 696968 4
5 594167 3996 19169 688558 5
6 595454 4308 13807 679237 6
7 590865 4143 29743 677362 7
8 589379 4429 25591 676693 8
9 584428 5219 29096 670009 9
10 573100 4929 26482 667209 10
11 567456 5761 22405 662976 11
12 569028 5592 27044 660194 12
13 620735 4163 17970 652270 13
14 628884 4962 18730 648024 14
15 628232 5208 19684 629295 15
16 612117 4755 19785 624961 16
17 595404 4491 18479 617306 17
18 597141 5732 10698 607691 18
19 593408 5731 31956 596219 19
20 590072 5040 29506 591130 20
21 579799 6102 34506 584528 21
22 574205 4904 27165 576798 22
23 572775 5369 26736 575683 23
24 572942 5578 23691 574369 24
25 619567 4619 18157 566815 25
26 625809 4731 17328 573074 26
27 619916 5011 18205 567739 27
28 587625 5299 20995 571942 28
29 565742 4146 17382 570274 29
30 557274 4625 9367 568800 30
31 560576 4736 31124 558115 31
32 548854 4219 26551 550591 32
33 531673 5116 30651 548872 33
34 525919 4205 25859 547009 34
35 511038 4121 25100 545946 35
36 498662 5103 25778 539702 36
37 555362 4300 20418 542427 37
38 564591 4578 18688 542968 38
39 541657 3809 20424 536640 39
40 527070 5526 24776 533653 40
41 509846 4248 19814 540996 41
42 514258 3830 12738 538316 42
43 516922 4428 31566 532646 43
44 507561 4834 30111 533390 44
45 492622 4406 30019 528715 45
46 490243 4565 31934 530664 46
47 469357 4104 25826 528564 47
48 477580 4798 26835 519107 48
49 528379 3935 20205 518703 49
50 533590 3792 17789 519059 50
51 517945 4387 20520 518498 51
52 506174 4006 22518 524575 52
53 501866 4078 15572 536046 53
54 516141 4724 11509 552006 54
55 528222 3157 25447 560687 55
56 532638 3558 24090 578884 56
57 536322 3899 27786 591491 57
58 536535 4118 26195 599228 58
59 523597 3790 20516 633019 59
60 536214 4278 22759 649918 60
61 586570 4035 19028 655509 61
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Bouw Auto Krediet t
4.345e+05 1.025e+01 -1.706e+00 2.522e-01 -1.066e+03
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-46362 -20136 -2925 16588 55554
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 4.345e+05 7.019e+04 6.191 7.42e-08 ***
Bouw 1.025e+01 6.153e+00 1.665 0.101409
Auto -1.706e+00 5.904e-01 -2.890 0.005471 **
Krediet 2.522e-01 8.158e-02 3.091 0.003103 **
t -1.066e+03 2.906e+02 -3.668 0.000546 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 25000 on 56 degrees of freedom
Multiple R-squared: 0.6765, Adjusted R-squared: 0.6534
F-statistic: 29.28 on 4 and 56 DF, p-value: 3.772e-13
> 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.07343641 0.14687281 0.9265636
[2,] 0.02847113 0.05694227 0.9715289
[3,] 0.02787386 0.05574772 0.9721261
[4,] 0.02063738 0.04127475 0.9793626
[5,] 0.01328415 0.02656829 0.9867159
[6,] 0.04049335 0.08098669 0.9595067
[7,] 0.10919911 0.21839822 0.8908009
[8,] 0.18656970 0.37313940 0.8134303
[9,] 0.14098365 0.28196730 0.8590164
[10,] 0.20318043 0.40636085 0.7968196
[11,] 0.18104022 0.36208043 0.8189598
[12,] 0.14451419 0.28902839 0.8554858
[13,] 0.10533722 0.21067444 0.8946628
[14,] 0.07016552 0.14033104 0.9298345
[15,] 0.08720596 0.17441191 0.9127940
[16,] 0.07824803 0.15649606 0.9217520
[17,] 0.07015735 0.14031469 0.9298427
[18,] 0.06894381 0.13788762 0.9310562
[19,] 0.08311960 0.16623921 0.9168804
[20,] 0.13451915 0.26903830 0.8654809
[21,] 0.17863296 0.35726592 0.8213670
[22,] 0.35471830 0.70943659 0.6452817
[23,] 0.48423974 0.96847948 0.5157603
[24,] 0.47286199 0.94572399 0.5271380
[25,] 0.45756913 0.91513826 0.5424309
[26,] 0.43793569 0.87587139 0.5620643
[27,] 0.42181893 0.84363787 0.5781811
[28,] 0.46878355 0.93756709 0.5312165
[29,] 0.55858852 0.88282296 0.4414115
[30,] 0.54750026 0.90499947 0.4524997
[31,] 0.67545774 0.64908452 0.3245423
[32,] 0.66148236 0.67703527 0.3385176
[33,] 0.70803686 0.58392628 0.2919631
[34,] 0.63176119 0.73647762 0.3682388
[35,] 0.56661752 0.86676495 0.4333825
[36,] 0.58877266 0.82245468 0.4112273
[37,] 0.61327068 0.77345864 0.3867293
[38,] 0.53724740 0.92550519 0.4627526
[39,] 0.49387307 0.98774614 0.5061269
[40,] 0.54343071 0.91313859 0.4565693
[41,] 0.58485319 0.83029362 0.4151468
[42,] 0.51847847 0.96304306 0.4815215
[43,] 0.59425015 0.81149970 0.4057498
[44,] 0.71650724 0.56698552 0.2834928
[45,] 0.75830864 0.48338272 0.2416914
[46,] 0.61015589 0.77968821 0.3898441
> postscript(file="/var/www/html/rcomp/tmp/1okoj1260641383.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/2kg1t1260641383.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/3h1nb1260641383.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/4gd6s1260641383.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/5pcf91260641383.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 = 61
Frequency = 1
1 2 3 4 5 6 7
-27748.963 -7634.596 -3165.614 -3327.654 -16924.885 -24567.414 1264.677
8 9 10 11 12 13 14
-9001.833 -13315.825 -24360.292 -43352.912 -30366.077 23565.516 26960.667
15 16 17 18 19 20 21
31205.071 22063.338 8823.699 -11941.762 24567.040 26480.774 16587.797
22 23 24 25 26 27 28
13759.388 8179.656 2406.708 52387.126 55554.394 50700.070 20224.440
29 30 31 32 33 34 35
5478.146 -20136.443 22912.644 11649.026 -6228.326 -9287.960 -23269.225
36 37 38 39 40 41 42
-41910.317 14251.104 18609.052 9179.074 -13757.127 -27137.881 -28774.460
43 44 45 46 47 48 49
2383.784 -12741.862 -21207.088 -21373.393 -46361.931 -40077.763 9419.620
50 51 52 53 54 55 56
12949.738 -2924.950 -7849.273 -26573.970 -28810.379 21986.769 16454.961
57 58 59 60 61
20837.663 15206.547 -11516.501 -3268.801 42866.995
> postscript(file="/var/www/html/rcomp/tmp/661xp1260641383.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 = 61
Frequency = 1
lag(myerror, k = 1) myerror
0 -27748.963 NA
1 -7634.596 -27748.963
2 -3165.614 -7634.596
3 -3327.654 -3165.614
4 -16924.885 -3327.654
5 -24567.414 -16924.885
6 1264.677 -24567.414
7 -9001.833 1264.677
8 -13315.825 -9001.833
9 -24360.292 -13315.825
10 -43352.912 -24360.292
11 -30366.077 -43352.912
12 23565.516 -30366.077
13 26960.667 23565.516
14 31205.071 26960.667
15 22063.338 31205.071
16 8823.699 22063.338
17 -11941.762 8823.699
18 24567.040 -11941.762
19 26480.774 24567.040
20 16587.797 26480.774
21 13759.388 16587.797
22 8179.656 13759.388
23 2406.708 8179.656
24 52387.126 2406.708
25 55554.394 52387.126
26 50700.070 55554.394
27 20224.440 50700.070
28 5478.146 20224.440
29 -20136.443 5478.146
30 22912.644 -20136.443
31 11649.026 22912.644
32 -6228.326 11649.026
33 -9287.960 -6228.326
34 -23269.225 -9287.960
35 -41910.317 -23269.225
36 14251.104 -41910.317
37 18609.052 14251.104
38 9179.074 18609.052
39 -13757.127 9179.074
40 -27137.881 -13757.127
41 -28774.460 -27137.881
42 2383.784 -28774.460
43 -12741.862 2383.784
44 -21207.088 -12741.862
45 -21373.393 -21207.088
46 -46361.931 -21373.393
47 -40077.763 -46361.931
48 9419.620 -40077.763
49 12949.738 9419.620
50 -2924.950 12949.738
51 -7849.273 -2924.950
52 -26573.970 -7849.273
53 -28810.379 -26573.970
54 21986.769 -28810.379
55 16454.961 21986.769
56 20837.663 16454.961
57 15206.547 20837.663
58 -11516.501 15206.547
59 -3268.801 -11516.501
60 42866.995 -3268.801
61 NA 42866.995
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -7634.596 -27748.963
[2,] -3165.614 -7634.596
[3,] -3327.654 -3165.614
[4,] -16924.885 -3327.654
[5,] -24567.414 -16924.885
[6,] 1264.677 -24567.414
[7,] -9001.833 1264.677
[8,] -13315.825 -9001.833
[9,] -24360.292 -13315.825
[10,] -43352.912 -24360.292
[11,] -30366.077 -43352.912
[12,] 23565.516 -30366.077
[13,] 26960.667 23565.516
[14,] 31205.071 26960.667
[15,] 22063.338 31205.071
[16,] 8823.699 22063.338
[17,] -11941.762 8823.699
[18,] 24567.040 -11941.762
[19,] 26480.774 24567.040
[20,] 16587.797 26480.774
[21,] 13759.388 16587.797
[22,] 8179.656 13759.388
[23,] 2406.708 8179.656
[24,] 52387.126 2406.708
[25,] 55554.394 52387.126
[26,] 50700.070 55554.394
[27,] 20224.440 50700.070
[28,] 5478.146 20224.440
[29,] -20136.443 5478.146
[30,] 22912.644 -20136.443
[31,] 11649.026 22912.644
[32,] -6228.326 11649.026
[33,] -9287.960 -6228.326
[34,] -23269.225 -9287.960
[35,] -41910.317 -23269.225
[36,] 14251.104 -41910.317
[37,] 18609.052 14251.104
[38,] 9179.074 18609.052
[39,] -13757.127 9179.074
[40,] -27137.881 -13757.127
[41,] -28774.460 -27137.881
[42,] 2383.784 -28774.460
[43,] -12741.862 2383.784
[44,] -21207.088 -12741.862
[45,] -21373.393 -21207.088
[46,] -46361.931 -21373.393
[47,] -40077.763 -46361.931
[48,] 9419.620 -40077.763
[49,] 12949.738 9419.620
[50,] -2924.950 12949.738
[51,] -7849.273 -2924.950
[52,] -26573.970 -7849.273
[53,] -28810.379 -26573.970
[54,] 21986.769 -28810.379
[55,] 16454.961 21986.769
[56,] 20837.663 16454.961
[57,] 15206.547 20837.663
[58,] -11516.501 15206.547
[59,] -3268.801 -11516.501
[60,] 42866.995 -3268.801
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -7634.596 -27748.963
2 -3165.614 -7634.596
3 -3327.654 -3165.614
4 -16924.885 -3327.654
5 -24567.414 -16924.885
6 1264.677 -24567.414
7 -9001.833 1264.677
8 -13315.825 -9001.833
9 -24360.292 -13315.825
10 -43352.912 -24360.292
11 -30366.077 -43352.912
12 23565.516 -30366.077
13 26960.667 23565.516
14 31205.071 26960.667
15 22063.338 31205.071
16 8823.699 22063.338
17 -11941.762 8823.699
18 24567.040 -11941.762
19 26480.774 24567.040
20 16587.797 26480.774
21 13759.388 16587.797
22 8179.656 13759.388
23 2406.708 8179.656
24 52387.126 2406.708
25 55554.394 52387.126
26 50700.070 55554.394
27 20224.440 50700.070
28 5478.146 20224.440
29 -20136.443 5478.146
30 22912.644 -20136.443
31 11649.026 22912.644
32 -6228.326 11649.026
33 -9287.960 -6228.326
34 -23269.225 -9287.960
35 -41910.317 -23269.225
36 14251.104 -41910.317
37 18609.052 14251.104
38 9179.074 18609.052
39 -13757.127 9179.074
40 -27137.881 -13757.127
41 -28774.460 -27137.881
42 2383.784 -28774.460
43 -12741.862 2383.784
44 -21207.088 -12741.862
45 -21373.393 -21207.088
46 -46361.931 -21373.393
47 -40077.763 -46361.931
48 9419.620 -40077.763
49 12949.738 9419.620
50 -2924.950 12949.738
51 -7849.273 -2924.950
52 -26573.970 -7849.273
53 -28810.379 -26573.970
54 21986.769 -28810.379
55 16454.961 21986.769
56 20837.663 16454.961
57 15206.547 20837.663
58 -11516.501 15206.547
59 -3268.801 -11516.501
60 42866.995 -3268.801
> 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/7nc2x1260641383.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/8v4wm1260641383.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/9rz3k1260641383.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/10r73f1260641383.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/11x8pa1260641383.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/125f0z1260641383.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/133uub1260641383.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/14bjns1260641383.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/15kk4e1260641383.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/16zo0w1260641383.tab")
+ }
>
> try(system("convert tmp/1okoj1260641383.ps tmp/1okoj1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/2kg1t1260641383.ps tmp/2kg1t1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/3h1nb1260641383.ps tmp/3h1nb1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/4gd6s1260641383.ps tmp/4gd6s1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/5pcf91260641383.ps tmp/5pcf91260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/661xp1260641383.ps tmp/661xp1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/7nc2x1260641383.ps tmp/7nc2x1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/8v4wm1260641383.ps tmp/8v4wm1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/9rz3k1260641383.ps tmp/9rz3k1260641383.png",intern=TRUE))
character(0)
> try(system("convert tmp/10r73f1260641383.ps tmp/10r73f1260641383.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.460 1.547 3.877