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(100
+ ,100
+ ,97.82226485
+ ,99.87129987
+ ,94.04971502
+ ,99.54459954
+ ,91.12460521
+ ,99.81189981
+ ,93.13202153
+ ,100.4851005
+ ,93.88342812
+ ,101.1385011
+ ,92.55349954
+ ,101.3662014
+ ,94.43494835
+ ,101.5147015
+ ,96.25017563
+ ,101.8216018
+ ,100.4355715
+ ,102.4354024
+ ,101.5036685
+ ,102.5344025
+ ,99.39789728
+ ,102.6532027
+ ,99.68990733
+ ,102.4651025
+ ,101.6895041
+ ,102.4354024
+ ,103.6652759
+ ,102.4156024
+ ,103.0532766
+ ,102.4453024
+ ,100.9500712
+ ,102.8908029
+ ,102.345366
+ ,102.8512029
+ ,101.6472299
+ ,103.3561034
+ ,99.56809393
+ ,103.7422037
+ ,95.67727392
+ ,103.7224037
+ ,96.58494865
+ ,104.0788041
+ ,96.32604937
+ ,104.2075042
+ ,95.37109101
+ ,103.9105039
+ ,96.00056203
+ ,103.7026037
+ ,96.88367859
+ ,103.960004
+ ,94.85280372
+ ,104.0986041
+ ,92.46943974
+ ,104.1481041
+ ,93.99180173
+ ,104.7124047
+ ,93.45262168
+ ,104.7223047
+ ,92.26698759
+ ,105.1975052
+ ,90.39653498
+ ,105.0688051
+ ,90.43001228
+ ,105.0589051
+ ,91.04995327
+ ,105.5044055
+ ,89.07845784
+ ,105.3757054
+ ,89.69314509
+ ,105.4747055
+ ,87.92459054
+ ,106.029106
+ ,85.8789319
+ ,107.019107
+ ,83.20612366
+ ,107.3161073
+ ,83.85722053
+ ,107.7517078
+ ,83.01393462
+ ,108.5239085
+ ,82.84508195
+ ,109.3159093
+ ,78.68864276
+ ,109.5634096
+ ,77.56959675
+ ,110.5435105
+ ,78.53689529
+ ,111.1573112
+ ,78.55717715
+ ,111.7414117
+ ,77.4761291
+ ,111.0583111
+ ,81.58931659
+ ,111.2365112
+ ,85.02428326
+ ,111.038511
+ ,91.71290159
+ ,110.3752104
+ ,95.96293061
+ ,110.1376101
+ ,90.84689022
+ ,110.2465102
+ ,92.28788036
+ ,110.6227106
+ ,95.56511274
+ ,109.98911
+ ,93.62452884
+ ,110.2168102
+ ,92.63071726
+ ,110.1376101
+ ,89.50914211
+ ,109.9297099
+ ,87.17171779
+ ,109.8604099
+ ,86.72624975
+ ,110.1970102
+ ,85.63212844
+ ,109.9099099)
+ ,dim=c(2
+ ,60)
+ ,dimnames=list(c('wisselkoers'
+ ,'consumptieprijzen')
+ ,1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('wisselkoers','consumptieprijzen'),1:60))
> 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
wisselkoers consumptieprijzen
1 100.00000 100.0000
2 97.82226 99.8713
3 94.04972 99.5446
4 91.12461 99.8119
5 93.13202 100.4851
6 93.88343 101.1385
7 92.55350 101.3662
8 94.43495 101.5147
9 96.25018 101.8216
10 100.43557 102.4354
11 101.50367 102.5344
12 99.39790 102.6532
13 99.68991 102.4651
14 101.68950 102.4354
15 103.66528 102.4156
16 103.05328 102.4453
17 100.95007 102.8908
18 102.34537 102.8512
19 101.64723 103.3561
20 99.56809 103.7422
21 95.67727 103.7224
22 96.58495 104.0788
23 96.32605 104.2075
24 95.37109 103.9105
25 96.00056 103.7026
26 96.88368 103.9600
27 94.85280 104.0986
28 92.46944 104.1481
29 93.99180 104.7124
30 93.45262 104.7223
31 92.26699 105.1975
32 90.39653 105.0688
33 90.43001 105.0589
34 91.04995 105.5044
35 89.07846 105.3757
36 89.69315 105.4747
37 87.92459 106.0291
38 85.87893 107.0191
39 83.20612 107.3161
40 83.85722 107.7517
41 83.01393 108.5239
42 82.84508 109.3159
43 78.68864 109.5634
44 77.56960 110.5435
45 78.53690 111.1573
46 78.55718 111.7414
47 77.47613 111.0583
48 81.58932 111.2365
49 85.02428 111.0385
50 91.71290 110.3752
51 95.96293 110.1376
52 90.84689 110.2465
53 92.28788 110.6227
54 95.56511 109.9891
55 93.62453 110.2168
56 92.63072 110.1376
57 89.50914 109.9297
58 87.17172 109.8604
59 86.72625 110.1970
60 85.63213 109.9099
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) consumptieprijzen
233.308 -1.336
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-8.8795 -3.9193 0.1305 4.1462 9.7493
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 233.3078 18.3956 12.683 < 2e-16 ***
consumptieprijzen -1.3355 0.1739 -7.679 2.11e-10 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 4.93 on 58 degrees of freedom
Multiple R-squared: 0.5042, Adjusted R-squared: 0.4956
F-statistic: 58.97 on 1 and 58 DF, p-value: 2.114e-10
> 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.47992908 0.95985816 0.5200709
[2,] 0.32394947 0.64789893 0.6760505
[3,] 0.21798300 0.43596599 0.7820170
[4,] 0.13919862 0.27839723 0.8608014
[5,] 0.10011620 0.20023241 0.8998838
[6,] 0.13222845 0.26445689 0.8677716
[7,] 0.12898604 0.25797209 0.8710140
[8,] 0.08078494 0.16156988 0.9192151
[9,] 0.04978178 0.09956356 0.9502182
[10,] 0.03946762 0.07893523 0.9605324
[11,] 0.04763526 0.09527052 0.9523647
[12,] 0.04315863 0.08631727 0.9568414
[13,] 0.02747488 0.05494977 0.9725251
[14,] 0.01987121 0.03974243 0.9801288
[15,] 0.01428420 0.02856839 0.9857158
[16,] 0.01379883 0.02759766 0.9862012
[17,] 0.02706331 0.05412662 0.9729367
[18,] 0.03316456 0.06632912 0.9668354
[19,] 0.03620064 0.07240129 0.9637994
[20,] 0.03696829 0.07393659 0.9630317
[21,] 0.03114733 0.06229465 0.9688527
[22,] 0.02534460 0.05068921 0.9746554
[23,] 0.02434558 0.04869116 0.9756544
[24,] 0.03111156 0.06222311 0.9688884
[25,] 0.02956003 0.05912006 0.9704400
[26,] 0.02780158 0.05560317 0.9721984
[27,] 0.02803717 0.05607433 0.9719628
[28,] 0.03191984 0.06383968 0.9680802
[29,] 0.03174574 0.06349149 0.9682543
[30,] 0.02724989 0.05449978 0.9727501
[31,] 0.02659530 0.05319060 0.9734047
[32,] 0.02250574 0.04501147 0.9774943
[33,] 0.02008159 0.04016318 0.9799184
[34,] 0.01837629 0.03675258 0.9816237
[35,] 0.02323203 0.04646406 0.9767680
[36,] 0.02547908 0.05095816 0.9745209
[37,] 0.03707921 0.07415841 0.9629208
[38,] 0.05294497 0.10588994 0.9470550
[39,] 0.28699533 0.57399066 0.7130047
[40,] 0.48910822 0.97821644 0.5108918
[41,] 0.46670405 0.93340810 0.5332959
[42,] 0.37699536 0.75399072 0.6230046
[43,] 0.51947892 0.96104215 0.4805211
[44,] 0.55173941 0.89652118 0.4482606
[45,] 0.70637581 0.58724837 0.2936242
[46,] 0.67489223 0.65021554 0.3251078
[47,] 0.80503089 0.38993822 0.1949691
[48,] 0.72456443 0.55087114 0.2754356
[49,] 0.64973652 0.70052697 0.3502635
[50,] 0.83231826 0.33536348 0.1676817
[51,] 0.78644901 0.42710197 0.2135510
> postscript(file="/var/www/html/rcomp/tmp/1ldt01258732662.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/2sdl31258732662.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/3ie931258732662.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/4c9o31258732662.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/51aob1258732662.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 = 60
Frequency = 1
1 2 3 4 5 6
0.247136208 -2.102484344 -6.311358654 -8.879475707 -5.972966505 -4.348911033
7 8 9 10 11 12
-5.374734579 -3.294956526 -1.069848686 3.935308303 5.135624842 3.188517177
13 14 15 16 17 18
3.229310088 5.189240903 7.138568821 6.566235343 5.058017938 6.400424976
19 20 21 22 23 24
6.376608515 4.813128629 0.895864738 2.279529864 2.192515945 0.840898966
25 26 27 28 29 30
1.192708967 2.419596383 0.573828815 -1.743425462 0.532587944 0.006629835
31 32 33 34 35 36
-0.544350438 -2.586688410 -2.566433050 -1.351504199 -3.494884990 -2.747978201
37 38 39 40 41 42
-3.776103409 -4.499566653 -6.775716274 -5.542853349 -5.354826957 -4.465923310
43 44 45 46 47 48
-8.291813584 -8.101886272 -6.314826480 -5.514449456 -7.507812209 -3.156629655
49 50 51 52 53 54
0.013897936 5.816645444 9.749347489 4.778748579 6.722172996 9.153200376
55 56 57 58 59 60
7.516721377 6.417134139 3.017897969 0.587920065 0.591998406 -0.885559582
> postscript(file="/var/www/html/rcomp/tmp/6xpjc1258732662.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 = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 0.247136208 NA
1 -2.102484344 0.247136208
2 -6.311358654 -2.102484344
3 -8.879475707 -6.311358654
4 -5.972966505 -8.879475707
5 -4.348911033 -5.972966505
6 -5.374734579 -4.348911033
7 -3.294956526 -5.374734579
8 -1.069848686 -3.294956526
9 3.935308303 -1.069848686
10 5.135624842 3.935308303
11 3.188517177 5.135624842
12 3.229310088 3.188517177
13 5.189240903 3.229310088
14 7.138568821 5.189240903
15 6.566235343 7.138568821
16 5.058017938 6.566235343
17 6.400424976 5.058017938
18 6.376608515 6.400424976
19 4.813128629 6.376608515
20 0.895864738 4.813128629
21 2.279529864 0.895864738
22 2.192515945 2.279529864
23 0.840898966 2.192515945
24 1.192708967 0.840898966
25 2.419596383 1.192708967
26 0.573828815 2.419596383
27 -1.743425462 0.573828815
28 0.532587944 -1.743425462
29 0.006629835 0.532587944
30 -0.544350438 0.006629835
31 -2.586688410 -0.544350438
32 -2.566433050 -2.586688410
33 -1.351504199 -2.566433050
34 -3.494884990 -1.351504199
35 -2.747978201 -3.494884990
36 -3.776103409 -2.747978201
37 -4.499566653 -3.776103409
38 -6.775716274 -4.499566653
39 -5.542853349 -6.775716274
40 -5.354826957 -5.542853349
41 -4.465923310 -5.354826957
42 -8.291813584 -4.465923310
43 -8.101886272 -8.291813584
44 -6.314826480 -8.101886272
45 -5.514449456 -6.314826480
46 -7.507812209 -5.514449456
47 -3.156629655 -7.507812209
48 0.013897936 -3.156629655
49 5.816645444 0.013897936
50 9.749347489 5.816645444
51 4.778748579 9.749347489
52 6.722172996 4.778748579
53 9.153200376 6.722172996
54 7.516721377 9.153200376
55 6.417134139 7.516721377
56 3.017897969 6.417134139
57 0.587920065 3.017897969
58 0.591998406 0.587920065
59 -0.885559582 0.591998406
60 NA -0.885559582
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -2.102484344 0.247136208
[2,] -6.311358654 -2.102484344
[3,] -8.879475707 -6.311358654
[4,] -5.972966505 -8.879475707
[5,] -4.348911033 -5.972966505
[6,] -5.374734579 -4.348911033
[7,] -3.294956526 -5.374734579
[8,] -1.069848686 -3.294956526
[9,] 3.935308303 -1.069848686
[10,] 5.135624842 3.935308303
[11,] 3.188517177 5.135624842
[12,] 3.229310088 3.188517177
[13,] 5.189240903 3.229310088
[14,] 7.138568821 5.189240903
[15,] 6.566235343 7.138568821
[16,] 5.058017938 6.566235343
[17,] 6.400424976 5.058017938
[18,] 6.376608515 6.400424976
[19,] 4.813128629 6.376608515
[20,] 0.895864738 4.813128629
[21,] 2.279529864 0.895864738
[22,] 2.192515945 2.279529864
[23,] 0.840898966 2.192515945
[24,] 1.192708967 0.840898966
[25,] 2.419596383 1.192708967
[26,] 0.573828815 2.419596383
[27,] -1.743425462 0.573828815
[28,] 0.532587944 -1.743425462
[29,] 0.006629835 0.532587944
[30,] -0.544350438 0.006629835
[31,] -2.586688410 -0.544350438
[32,] -2.566433050 -2.586688410
[33,] -1.351504199 -2.566433050
[34,] -3.494884990 -1.351504199
[35,] -2.747978201 -3.494884990
[36,] -3.776103409 -2.747978201
[37,] -4.499566653 -3.776103409
[38,] -6.775716274 -4.499566653
[39,] -5.542853349 -6.775716274
[40,] -5.354826957 -5.542853349
[41,] -4.465923310 -5.354826957
[42,] -8.291813584 -4.465923310
[43,] -8.101886272 -8.291813584
[44,] -6.314826480 -8.101886272
[45,] -5.514449456 -6.314826480
[46,] -7.507812209 -5.514449456
[47,] -3.156629655 -7.507812209
[48,] 0.013897936 -3.156629655
[49,] 5.816645444 0.013897936
[50,] 9.749347489 5.816645444
[51,] 4.778748579 9.749347489
[52,] 6.722172996 4.778748579
[53,] 9.153200376 6.722172996
[54,] 7.516721377 9.153200376
[55,] 6.417134139 7.516721377
[56,] 3.017897969 6.417134139
[57,] 0.587920065 3.017897969
[58,] 0.591998406 0.587920065
[59,] -0.885559582 0.591998406
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -2.102484344 0.247136208
2 -6.311358654 -2.102484344
3 -8.879475707 -6.311358654
4 -5.972966505 -8.879475707
5 -4.348911033 -5.972966505
6 -5.374734579 -4.348911033
7 -3.294956526 -5.374734579
8 -1.069848686 -3.294956526
9 3.935308303 -1.069848686
10 5.135624842 3.935308303
11 3.188517177 5.135624842
12 3.229310088 3.188517177
13 5.189240903 3.229310088
14 7.138568821 5.189240903
15 6.566235343 7.138568821
16 5.058017938 6.566235343
17 6.400424976 5.058017938
18 6.376608515 6.400424976
19 4.813128629 6.376608515
20 0.895864738 4.813128629
21 2.279529864 0.895864738
22 2.192515945 2.279529864
23 0.840898966 2.192515945
24 1.192708967 0.840898966
25 2.419596383 1.192708967
26 0.573828815 2.419596383
27 -1.743425462 0.573828815
28 0.532587944 -1.743425462
29 0.006629835 0.532587944
30 -0.544350438 0.006629835
31 -2.586688410 -0.544350438
32 -2.566433050 -2.586688410
33 -1.351504199 -2.566433050
34 -3.494884990 -1.351504199
35 -2.747978201 -3.494884990
36 -3.776103409 -2.747978201
37 -4.499566653 -3.776103409
38 -6.775716274 -4.499566653
39 -5.542853349 -6.775716274
40 -5.354826957 -5.542853349
41 -4.465923310 -5.354826957
42 -8.291813584 -4.465923310
43 -8.101886272 -8.291813584
44 -6.314826480 -8.101886272
45 -5.514449456 -6.314826480
46 -7.507812209 -5.514449456
47 -3.156629655 -7.507812209
48 0.013897936 -3.156629655
49 5.816645444 0.013897936
50 9.749347489 5.816645444
51 4.778748579 9.749347489
52 6.722172996 4.778748579
53 9.153200376 6.722172996
54 7.516721377 9.153200376
55 6.417134139 7.516721377
56 3.017897969 6.417134139
57 0.587920065 3.017897969
58 0.591998406 0.587920065
59 -0.885559582 0.591998406
> 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/7ccmm1258732662.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/89amh1258732662.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/9s4nc1258732662.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/10bb2u1258732662.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/111r4s1258732662.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/12fv8q1258732662.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/135bsh1258732662.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/14gjao1258732662.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/15e2ew1258732662.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/16ol4d1258732662.tab")
+ }
>
> system("convert tmp/1ldt01258732662.ps tmp/1ldt01258732662.png")
> system("convert tmp/2sdl31258732662.ps tmp/2sdl31258732662.png")
> system("convert tmp/3ie931258732662.ps tmp/3ie931258732662.png")
> system("convert tmp/4c9o31258732662.ps tmp/4c9o31258732662.png")
> system("convert tmp/51aob1258732662.ps tmp/51aob1258732662.png")
> system("convert tmp/6xpjc1258732662.ps tmp/6xpjc1258732662.png")
> system("convert tmp/7ccmm1258732662.ps tmp/7ccmm1258732662.png")
> system("convert tmp/89amh1258732662.ps tmp/89amh1258732662.png")
> system("convert tmp/9s4nc1258732662.ps tmp/9s4nc1258732662.png")
> system("convert tmp/10bb2u1258732662.ps tmp/10bb2u1258732662.png")
>
>
> proc.time()
user system elapsed
2.482 1.576 2.888