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 = '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
Werkzoekend Bouw Auto Krediet
1 593530 3922 18004 707169
2 610763 3759 17537 703434
3 612613 4138 20366 701017
4 611324 4634 22782 696968
5 594167 3996 19169 688558
6 595454 4308 13807 679237
7 590865 4143 29743 677362
8 589379 4429 25591 676693
9 584428 5219 29096 670009
10 573100 4929 26482 667209
11 567456 5761 22405 662976
12 569028 5592 27044 660194
13 620735 4163 17970 652270
14 628884 4962 18730 648024
15 628232 5208 19684 629295
16 612117 4755 19785 624961
17 595404 4491 18479 617306
18 597141 5732 10698 607691
19 593408 5731 31956 596219
20 590072 5040 29506 591130
21 579799 6102 34506 584528
22 574205 4904 27165 576798
23 572775 5369 26736 575683
24 572942 5578 23691 574369
25 619567 4619 18157 566815
26 625809 4731 17328 573074
27 619916 5011 18205 567739
28 587625 5299 20995 571942
29 565742 4146 17382 570274
30 557274 4625 9367 568800
31 560576 4736 31124 558115
32 548854 4219 26551 550591
33 531673 5116 30651 548872
34 525919 4205 25859 547009
35 511038 4121 25100 545946
36 498662 5103 25778 539702
37 555362 4300 20418 542427
38 564591 4578 18688 542968
39 541657 3809 20424 536640
40 527070 5526 24776 533653
41 509846 4248 19814 540996
42 514258 3830 12738 538316
43 516922 4428 31566 532646
44 507561 4834 30111 533390
45 492622 4406 30019 528715
46 490243 4565 31934 530664
47 469357 4104 25826 528564
48 477580 4798 26835 519107
49 528379 3935 20205 518703
50 533590 3792 17789 519059
51 517945 4387 20520 518498
52 506174 4006 22518 524575
53 501866 4078 15572 536046
54 516141 4724 11509 552006
55 528222 3157 25447 560687
56 532638 3558 24090 578884
57 536322 3899 27786 591491
58 536535 4118 26195 599228
59 523597 3790 20516 633019
60 536214 4278 22759 649918
61 586570 4035 19028 655509
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Bouw Auto Krediet
2.265e+05 2.146e+01 -1.903e+00 4.696e-01
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-51332 -17721 1547 19131 62357
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.265e+05 4.563e+04 4.963 6.61e-06 ***
Bouw 2.146e+01 5.894e+00 3.640 0.000588 ***
Auto -1.903e+00 6.490e-01 -2.932 0.004847 **
Krediet 4.696e-01 6.188e-02 7.589 3.32e-10 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 27600 on 57 degrees of freedom
Multiple R-squared: 0.5988, Adjusted R-squared: 0.5777
F-statistic: 28.36 on 3 and 57 DF, p-value: 2.374e-11
> 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.0499639627 0.099927925 0.95003604
[2,] 0.0157654347 0.031530869 0.98423457
[3,] 0.0063654113 0.012730823 0.99363459
[4,] 0.0035856808 0.007171362 0.99641432
[5,] 0.0019834350 0.003966870 0.99801656
[6,] 0.0007789305 0.001557861 0.99922107
[7,] 0.0070541444 0.014108289 0.99294586
[8,] 0.0229533431 0.045906686 0.97704666
[9,] 0.0201257554 0.040251511 0.97987424
[10,] 0.0123788161 0.024757632 0.98762118
[11,] 0.0176833512 0.035366702 0.98231665
[12,] 0.0128242890 0.025648578 0.98717571
[13,] 0.0071460942 0.014292188 0.99285391
[14,] 0.0049094238 0.009818848 0.99509058
[15,] 0.0024519644 0.004903929 0.99754804
[16,] 0.0035638702 0.007127740 0.99643613
[17,] 0.0028091453 0.005618291 0.99719085
[18,] 0.0020219043 0.004043809 0.99797810
[19,] 0.0033896804 0.006779361 0.99661032
[20,] 0.0086967353 0.017393471 0.99130326
[21,] 0.0255640651 0.051128130 0.97443593
[22,] 0.0376553090 0.075310618 0.96234469
[23,] 0.1415942287 0.283188457 0.85840577
[24,] 0.2990754063 0.598150813 0.70092459
[25,] 0.3963401298 0.792680260 0.60365987
[26,] 0.4877317943 0.975463589 0.51226821
[27,] 0.5685241372 0.862951726 0.43147586
[28,] 0.6157873416 0.768425317 0.38421266
[29,] 0.6737230714 0.652553857 0.32627693
[30,] 0.7810802333 0.437839533 0.21891977
[31,] 0.8101435876 0.379712825 0.18985641
[32,] 0.9001812241 0.199637552 0.09981878
[33,] 0.9050888583 0.189822283 0.09491114
[34,] 0.9365847056 0.126830589 0.06341529
[35,] 0.9317517197 0.136496561 0.06824828
[36,] 0.9251967674 0.149606465 0.07480323
[37,] 0.9175288342 0.164942332 0.08247117
[38,] 0.9156998317 0.168600337 0.08430017
[39,] 0.8839539159 0.232092168 0.11604608
[40,] 0.8463391838 0.307321632 0.15366082
[41,] 0.9248779814 0.150244037 0.07512202
[42,] 0.9251837043 0.149632591 0.07481630
[43,] 0.8941103587 0.211779283 0.10588964
[44,] 0.8902493495 0.219501301 0.10975065
[45,] 0.8449707628 0.310058474 0.15502924
[46,] 0.7436805414 0.512638917 0.25631946
[47,] 0.6450495395 0.709900921 0.35495046
[48,] 0.5283504822 0.943299036 0.47164952
> postscript(file="/var/www/html/rcomp/tmp/1nbex1260641519.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/2xbxr1260641519.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/36t9t1260641519.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/4ijov1260641519.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/5slb11260641519.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
-14923.681 6672.298 6907.426 1473.578 -4918.450 -16151.219 14001.400
8 9 10 11 12 13 14
-1207.145 -13302.202 -22066.150 -51332.280 -36001.144 32825.271 27269.465
15 16 17 18 19 20 21
31948.926 27781.623 17843.294 -17338.118 24783.982 34003.490 13555.896
22 23 24 25 26 27 28
23330.777 11630.307 2136.121 62357.054 61679.342 53952.068 18815.965
29 30 31 32 33 34 35
15582.681 -17721.114 29612.581 23816.575 -4003.886 1547.415 -12476.077
36 37 38 39 40 41 42
-41701.419 20751.274 20469.406 20310.924 -21435.935 -24126.145 -22949.482
43 44 45 46 47 48 49
5368.473 -15822.078 -19556.904 -22619.331 -44248.584 -44556.547 12335.586
50 51 52 53 54 55 56
15851.061 -7101.689 -9749.506 -34204.913 -49016.679 19131.148 3815.616
57 58 59 60 61
1294.633 -9851.923 -42424.816 -43947.103 1898.861
> postscript(file="/var/www/html/rcomp/tmp/6bv721260641519.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 -14923.681 NA
1 6672.298 -14923.681
2 6907.426 6672.298
3 1473.578 6907.426
4 -4918.450 1473.578
5 -16151.219 -4918.450
6 14001.400 -16151.219
7 -1207.145 14001.400
8 -13302.202 -1207.145
9 -22066.150 -13302.202
10 -51332.280 -22066.150
11 -36001.144 -51332.280
12 32825.271 -36001.144
13 27269.465 32825.271
14 31948.926 27269.465
15 27781.623 31948.926
16 17843.294 27781.623
17 -17338.118 17843.294
18 24783.982 -17338.118
19 34003.490 24783.982
20 13555.896 34003.490
21 23330.777 13555.896
22 11630.307 23330.777
23 2136.121 11630.307
24 62357.054 2136.121
25 61679.342 62357.054
26 53952.068 61679.342
27 18815.965 53952.068
28 15582.681 18815.965
29 -17721.114 15582.681
30 29612.581 -17721.114
31 23816.575 29612.581
32 -4003.886 23816.575
33 1547.415 -4003.886
34 -12476.077 1547.415
35 -41701.419 -12476.077
36 20751.274 -41701.419
37 20469.406 20751.274
38 20310.924 20469.406
39 -21435.935 20310.924
40 -24126.145 -21435.935
41 -22949.482 -24126.145
42 5368.473 -22949.482
43 -15822.078 5368.473
44 -19556.904 -15822.078
45 -22619.331 -19556.904
46 -44248.584 -22619.331
47 -44556.547 -44248.584
48 12335.586 -44556.547
49 15851.061 12335.586
50 -7101.689 15851.061
51 -9749.506 -7101.689
52 -34204.913 -9749.506
53 -49016.679 -34204.913
54 19131.148 -49016.679
55 3815.616 19131.148
56 1294.633 3815.616
57 -9851.923 1294.633
58 -42424.816 -9851.923
59 -43947.103 -42424.816
60 1898.861 -43947.103
61 NA 1898.861
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 6672.298 -14923.681
[2,] 6907.426 6672.298
[3,] 1473.578 6907.426
[4,] -4918.450 1473.578
[5,] -16151.219 -4918.450
[6,] 14001.400 -16151.219
[7,] -1207.145 14001.400
[8,] -13302.202 -1207.145
[9,] -22066.150 -13302.202
[10,] -51332.280 -22066.150
[11,] -36001.144 -51332.280
[12,] 32825.271 -36001.144
[13,] 27269.465 32825.271
[14,] 31948.926 27269.465
[15,] 27781.623 31948.926
[16,] 17843.294 27781.623
[17,] -17338.118 17843.294
[18,] 24783.982 -17338.118
[19,] 34003.490 24783.982
[20,] 13555.896 34003.490
[21,] 23330.777 13555.896
[22,] 11630.307 23330.777
[23,] 2136.121 11630.307
[24,] 62357.054 2136.121
[25,] 61679.342 62357.054
[26,] 53952.068 61679.342
[27,] 18815.965 53952.068
[28,] 15582.681 18815.965
[29,] -17721.114 15582.681
[30,] 29612.581 -17721.114
[31,] 23816.575 29612.581
[32,] -4003.886 23816.575
[33,] 1547.415 -4003.886
[34,] -12476.077 1547.415
[35,] -41701.419 -12476.077
[36,] 20751.274 -41701.419
[37,] 20469.406 20751.274
[38,] 20310.924 20469.406
[39,] -21435.935 20310.924
[40,] -24126.145 -21435.935
[41,] -22949.482 -24126.145
[42,] 5368.473 -22949.482
[43,] -15822.078 5368.473
[44,] -19556.904 -15822.078
[45,] -22619.331 -19556.904
[46,] -44248.584 -22619.331
[47,] -44556.547 -44248.584
[48,] 12335.586 -44556.547
[49,] 15851.061 12335.586
[50,] -7101.689 15851.061
[51,] -9749.506 -7101.689
[52,] -34204.913 -9749.506
[53,] -49016.679 -34204.913
[54,] 19131.148 -49016.679
[55,] 3815.616 19131.148
[56,] 1294.633 3815.616
[57,] -9851.923 1294.633
[58,] -42424.816 -9851.923
[59,] -43947.103 -42424.816
[60,] 1898.861 -43947.103
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 6672.298 -14923.681
2 6907.426 6672.298
3 1473.578 6907.426
4 -4918.450 1473.578
5 -16151.219 -4918.450
6 14001.400 -16151.219
7 -1207.145 14001.400
8 -13302.202 -1207.145
9 -22066.150 -13302.202
10 -51332.280 -22066.150
11 -36001.144 -51332.280
12 32825.271 -36001.144
13 27269.465 32825.271
14 31948.926 27269.465
15 27781.623 31948.926
16 17843.294 27781.623
17 -17338.118 17843.294
18 24783.982 -17338.118
19 34003.490 24783.982
20 13555.896 34003.490
21 23330.777 13555.896
22 11630.307 23330.777
23 2136.121 11630.307
24 62357.054 2136.121
25 61679.342 62357.054
26 53952.068 61679.342
27 18815.965 53952.068
28 15582.681 18815.965
29 -17721.114 15582.681
30 29612.581 -17721.114
31 23816.575 29612.581
32 -4003.886 23816.575
33 1547.415 -4003.886
34 -12476.077 1547.415
35 -41701.419 -12476.077
36 20751.274 -41701.419
37 20469.406 20751.274
38 20310.924 20469.406
39 -21435.935 20310.924
40 -24126.145 -21435.935
41 -22949.482 -24126.145
42 5368.473 -22949.482
43 -15822.078 5368.473
44 -19556.904 -15822.078
45 -22619.331 -19556.904
46 -44248.584 -22619.331
47 -44556.547 -44248.584
48 12335.586 -44556.547
49 15851.061 12335.586
50 -7101.689 15851.061
51 -9749.506 -7101.689
52 -34204.913 -9749.506
53 -49016.679 -34204.913
54 19131.148 -49016.679
55 3815.616 19131.148
56 1294.633 3815.616
57 -9851.923 1294.633
58 -42424.816 -9851.923
59 -43947.103 -42424.816
60 1898.861 -43947.103
> 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/7qr991260641519.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/8fmnv1260641519.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/926b31260641520.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/108txm1260641520.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/117uyz1260641520.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/127b0q1260641520.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/13kcs31260641520.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/14m3h11260641520.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/150xdo1260641520.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/16ttri1260641520.tab")
+ }
>
> try(system("convert tmp/1nbex1260641519.ps tmp/1nbex1260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/2xbxr1260641519.ps tmp/2xbxr1260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/36t9t1260641519.ps tmp/36t9t1260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/4ijov1260641519.ps tmp/4ijov1260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/5slb11260641519.ps tmp/5slb11260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/6bv721260641519.ps tmp/6bv721260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/7qr991260641519.ps tmp/7qr991260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/8fmnv1260641519.ps tmp/8fmnv1260641519.png",intern=TRUE))
character(0)
> try(system("convert tmp/926b31260641520.ps tmp/926b31260641520.png",intern=TRUE))
character(0)
> try(system("convert tmp/108txm1260641520.ps tmp/108txm1260641520.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.487 1.547 2.955