R version 2.13.0 (2011-04-13)
Copyright (C) 2011 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: i486-pc-linux-gnu (32-bit)
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(129.99,30,94,59.99,12,85.5,49.99,15,86,84.99,40,94,179.99,512,109,329.99,1500,118,25.99,16,72,499.99,8000,140,89.99,7,102.8,119.99,20,99.8,79.99,128,80,199.99,256,106,449.99,256,122,549.99,4000,161,529.99,8000,135,639.99,16000,140,749.99,32000,140,399.99,130,135,169.99,256,109,189.99,8000,135,199.99,8000,135,69.99,20,90,69.99,20,90,109.99,5,81,159.99,128,104,159.99,128,104,199.99,1000,135,75,30,81,349.99,512,126,439.99,8000,140,309.99,512,120,379.99,512,120,349.99,512,110,169.99,256,108,239.99,192,120,229.99,512,118,69.99,64,85,99.99,20,94,29.99,8,72.6,39.99,12,78,21.99,8,65,499.99,60,130,29.99,1,70,29.99,4,78.5,49.99,32,93.5,49.99,10,80,55.99,10,78.8,59.99,9,90.3,79.99,30,87.7,139.99,51,107,159.99,16000,90,169.99,46,103,229.99,32000,126,249.99,16000,98,309.99,256,128,499.99,16000,132,65.99,7,94,89.99,48,111,89.99,100,95,449.99,16000,155),dim=c(3,60),dimnames=list(c('Prijs','Geheugen','Gewicht'),1:60))
> y <- array(NA,dim=c(3,60),dimnames=list(c('Prijs','Geheugen','Gewicht'),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'
> library(lattice)
> library(lmtest)
Loading required package: zoo
> 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
Prijs Geheugen Gewicht
1 129.99 30 94.0
2 59.99 12 85.5
3 49.99 15 86.0
4 84.99 40 94.0
5 179.99 512 109.0
6 329.99 1500 118.0
7 25.99 16 72.0
8 499.99 8000 140.0
9 89.99 7 102.8
10 119.99 20 99.8
11 79.99 128 80.0
12 199.99 256 106.0
13 449.99 256 122.0
14 549.99 4000 161.0
15 529.99 8000 135.0
16 639.99 16000 140.0
17 749.99 32000 140.0
18 399.99 130 135.0
19 169.99 256 109.0
20 189.99 8000 135.0
21 199.99 8000 135.0
22 69.99 20 90.0
23 69.99 20 90.0
24 109.99 5 81.0
25 159.99 128 104.0
26 159.99 128 104.0
27 199.99 1000 135.0
28 75.00 30 81.0
29 349.99 512 126.0
30 439.99 8000 140.0
31 309.99 512 120.0
32 379.99 512 120.0
33 349.99 512 110.0
34 169.99 256 108.0
35 239.99 192 120.0
36 229.99 512 118.0
37 69.99 64 85.0
38 99.99 20 94.0
39 29.99 8 72.6
40 39.99 12 78.0
41 21.99 8 65.0
42 499.99 60 130.0
43 29.99 1 70.0
44 29.99 4 78.5
45 49.99 32 93.5
46 49.99 10 80.0
47 55.99 10 78.8
48 59.99 9 90.3
49 79.99 30 87.7
50 139.99 51 107.0
51 159.99 16000 90.0
52 169.99 46 103.0
53 229.99 32000 126.0
54 249.99 16000 98.0
55 309.99 256 128.0
56 499.99 16000 132.0
57 65.99 7 94.0
58 89.99 48 111.0
59 89.99 100 95.0
60 449.99 16000 155.0
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Geheugen Gewicht
-4.250e+02 4.931e-03 5.815e+00
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-235.541 -34.619 -2.879 39.165 203.045
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -4.250e+02 5.776e+01 -7.358 8.05e-10 ***
Geheugen 4.931e-03 1.797e-03 2.744 0.0081 **
Gewicht 5.815e+00 5.553e-01 10.472 6.74e-15 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 87.59 on 57 degrees of freedom
Multiple R-squared: 0.7627, Adjusted R-squared: 0.7544
F-statistic: 91.63 on 2 and 57 DF, p-value: < 2.2e-16
> if (n > n25) {
+ kp3 <- k + 3
+ nmkm3 <- n - k - 3
+ gqarr <- array(NA, dim=c(nmkm3-kp3+1,3))
+ numgqtests <- 0
+ numsignificant1 <- 0
+ numsignificant5 <- 0
+ numsignificant10 <- 0
+ for (mypoint in kp3:nmkm3) {
+ j <- 0
+ numgqtests <- numgqtests + 1
+ for (myalt in c('greater', 'two.sided', 'less')) {
+ j <- j + 1
+ gqarr[mypoint-kp3+1,j] <- gqtest(mylm, point=mypoint, alternative=myalt)$p.value
+ }
+ if (gqarr[mypoint-kp3+1,2] < 0.01) numsignificant1 <- numsignificant1 + 1
+ if (gqarr[mypoint-kp3+1,2] < 0.05) numsignificant5 <- numsignificant5 + 1
+ if (gqarr[mypoint-kp3+1,2] < 0.10) numsignificant10 <- numsignificant10 + 1
+ }
+ gqarr
+ }
[,1] [,2] [,3]
[1,] 0.025825704 0.051651408 0.97417430
[2,] 0.005885898 0.011771796 0.99411410
[3,] 0.022031385 0.044062770 0.97796862
[4,] 0.029253792 0.058507584 0.97074621
[5,] 0.012328346 0.024656693 0.98767165
[6,] 0.006753958 0.013507917 0.99324604
[7,] 0.003264222 0.006528443 0.99673578
[8,] 0.084627875 0.169255750 0.91537212
[9,] 0.063714329 0.127428658 0.93628567
[10,] 0.056617565 0.113235130 0.94338244
[11,] 0.052098690 0.104197380 0.94790131
[12,] 0.130952693 0.261905386 0.86904731
[13,] 0.094501210 0.189002421 0.90549879
[14,] 0.079508670 0.159017339 0.92049133
[15,] 0.766897332 0.466205336 0.23310267
[16,] 0.963963856 0.072072289 0.03603614
[17,] 0.947615305 0.104769390 0.05238470
[18,] 0.926106728 0.147786543 0.07389327
[19,] 0.911873014 0.176253972 0.08812699
[20,] 0.877988355 0.244023291 0.12201165
[21,] 0.836138783 0.327722434 0.16386122
[22,] 0.931509588 0.136980824 0.06849041
[23,] 0.904274164 0.191451673 0.09572584
[24,] 0.880339282 0.239321436 0.11966072
[25,] 0.841305399 0.317389201 0.15869460
[26,] 0.802634420 0.394731160 0.19736558
[27,] 0.843955514 0.312088971 0.15604449
[28,] 0.918072808 0.163854385 0.08192719
[29,] 0.888223393 0.223553214 0.11177661
[30,] 0.849211997 0.301576006 0.15078800
[31,] 0.802619508 0.394760984 0.19738049
[32,] 0.741237867 0.517524267 0.25876213
[33,] 0.676768616 0.646462768 0.32323138
[34,] 0.605475523 0.789048953 0.39452448
[35,] 0.523086351 0.953827297 0.47691365
[36,] 0.486516477 0.973032954 0.51348352
[37,] 0.836069302 0.327861395 0.16393070
[38,] 0.793461085 0.413077831 0.20653892
[39,] 0.721541360 0.556917280 0.27845864
[40,] 0.680639202 0.638721595 0.31936080
[41,] 0.591376390 0.817247219 0.40862361
[42,] 0.505224070 0.989551859 0.49477593
[43,] 0.413334228 0.826668455 0.58666577
[44,] 0.314276851 0.628553703 0.68572315
[45,] 0.237632260 0.475264519 0.76236774
[46,] 0.186070252 0.372140504 0.81392975
[47,] 0.117169818 0.234339635 0.88283018
[48,] 0.593418364 0.813163272 0.40658164
[49,] 0.429609682 0.859219364 0.57039032
> postscript(file="/var/wessaorg/rcomp/tmp/1xa4a1323898689.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/wessaorg/rcomp/tmp/2jhox1323898689.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/wessaorg/rcomp/tmp/3sris1323898689.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/wessaorg/rcomp/tmp/429p51323898689.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/wessaorg/rcomp/tmp/538f31323898689.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 = 60
Frequency = 1
1 2 3 4 5 6
8.1966032 -12.2847983 -25.2072292 -36.8527082 -31.4093335 61.3812216
7 8 9 10 11 12
32.2016904 71.3927402 -82.8644010 -35.4826807 39.1272019 7.2988642
13 14 15 16 17 18
164.2544634 18.9965365 130.4691155 171.9435954 203.0453057 39.2772117
19 20 21 22 23 24
-40.1469609 -209.5308845 -199.5308845 -28.4929851 -28.4929851 63.9184575
25 26 27 28 29 30
-20.4393993 -20.4393993 -165.0128828 28.8051789 39.7309906 11.3927402
31 32 33 34 35 36
34.6226409 104.6226409 132.7753914 -34.3316859 -33.7993933 -33.7468090
37 38 39 40 41 42
0.3664198 -21.7540853 32.7519745 11.3297646 68.9480649 168.6987670
43 44 45 46 47 48
47.9062076 -1.5384237 -68.9056215 9.7090768 22.6874069 -40.1833251
49 50 51 52 53 54
-5.1671639 -57.5055265 -17.2926520 -4.2197705 -235.5408436 26.1851476
55 56 57 58 59 60
-10.6371869 78.4657958 -55.6899805 -130.7518333 -37.9638518 -105.2855304
> postscript(file="/var/wessaorg/rcomp/tmp/6mz4x1323898689.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 = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 8.1966032 NA
1 -12.2847983 8.1966032
2 -25.2072292 -12.2847983
3 -36.8527082 -25.2072292
4 -31.4093335 -36.8527082
5 61.3812216 -31.4093335
6 32.2016904 61.3812216
7 71.3927402 32.2016904
8 -82.8644010 71.3927402
9 -35.4826807 -82.8644010
10 39.1272019 -35.4826807
11 7.2988642 39.1272019
12 164.2544634 7.2988642
13 18.9965365 164.2544634
14 130.4691155 18.9965365
15 171.9435954 130.4691155
16 203.0453057 171.9435954
17 39.2772117 203.0453057
18 -40.1469609 39.2772117
19 -209.5308845 -40.1469609
20 -199.5308845 -209.5308845
21 -28.4929851 -199.5308845
22 -28.4929851 -28.4929851
23 63.9184575 -28.4929851
24 -20.4393993 63.9184575
25 -20.4393993 -20.4393993
26 -165.0128828 -20.4393993
27 28.8051789 -165.0128828
28 39.7309906 28.8051789
29 11.3927402 39.7309906
30 34.6226409 11.3927402
31 104.6226409 34.6226409
32 132.7753914 104.6226409
33 -34.3316859 132.7753914
34 -33.7993933 -34.3316859
35 -33.7468090 -33.7993933
36 0.3664198 -33.7468090
37 -21.7540853 0.3664198
38 32.7519745 -21.7540853
39 11.3297646 32.7519745
40 68.9480649 11.3297646
41 168.6987670 68.9480649
42 47.9062076 168.6987670
43 -1.5384237 47.9062076
44 -68.9056215 -1.5384237
45 9.7090768 -68.9056215
46 22.6874069 9.7090768
47 -40.1833251 22.6874069
48 -5.1671639 -40.1833251
49 -57.5055265 -5.1671639
50 -17.2926520 -57.5055265
51 -4.2197705 -17.2926520
52 -235.5408436 -4.2197705
53 26.1851476 -235.5408436
54 -10.6371869 26.1851476
55 78.4657958 -10.6371869
56 -55.6899805 78.4657958
57 -130.7518333 -55.6899805
58 -37.9638518 -130.7518333
59 -105.2855304 -37.9638518
60 NA -105.2855304
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -12.2847983 8.1966032
[2,] -25.2072292 -12.2847983
[3,] -36.8527082 -25.2072292
[4,] -31.4093335 -36.8527082
[5,] 61.3812216 -31.4093335
[6,] 32.2016904 61.3812216
[7,] 71.3927402 32.2016904
[8,] -82.8644010 71.3927402
[9,] -35.4826807 -82.8644010
[10,] 39.1272019 -35.4826807
[11,] 7.2988642 39.1272019
[12,] 164.2544634 7.2988642
[13,] 18.9965365 164.2544634
[14,] 130.4691155 18.9965365
[15,] 171.9435954 130.4691155
[16,] 203.0453057 171.9435954
[17,] 39.2772117 203.0453057
[18,] -40.1469609 39.2772117
[19,] -209.5308845 -40.1469609
[20,] -199.5308845 -209.5308845
[21,] -28.4929851 -199.5308845
[22,] -28.4929851 -28.4929851
[23,] 63.9184575 -28.4929851
[24,] -20.4393993 63.9184575
[25,] -20.4393993 -20.4393993
[26,] -165.0128828 -20.4393993
[27,] 28.8051789 -165.0128828
[28,] 39.7309906 28.8051789
[29,] 11.3927402 39.7309906
[30,] 34.6226409 11.3927402
[31,] 104.6226409 34.6226409
[32,] 132.7753914 104.6226409
[33,] -34.3316859 132.7753914
[34,] -33.7993933 -34.3316859
[35,] -33.7468090 -33.7993933
[36,] 0.3664198 -33.7468090
[37,] -21.7540853 0.3664198
[38,] 32.7519745 -21.7540853
[39,] 11.3297646 32.7519745
[40,] 68.9480649 11.3297646
[41,] 168.6987670 68.9480649
[42,] 47.9062076 168.6987670
[43,] -1.5384237 47.9062076
[44,] -68.9056215 -1.5384237
[45,] 9.7090768 -68.9056215
[46,] 22.6874069 9.7090768
[47,] -40.1833251 22.6874069
[48,] -5.1671639 -40.1833251
[49,] -57.5055265 -5.1671639
[50,] -17.2926520 -57.5055265
[51,] -4.2197705 -17.2926520
[52,] -235.5408436 -4.2197705
[53,] 26.1851476 -235.5408436
[54,] -10.6371869 26.1851476
[55,] 78.4657958 -10.6371869
[56,] -55.6899805 78.4657958
[57,] -130.7518333 -55.6899805
[58,] -37.9638518 -130.7518333
[59,] -105.2855304 -37.9638518
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -12.2847983 8.1966032
2 -25.2072292 -12.2847983
3 -36.8527082 -25.2072292
4 -31.4093335 -36.8527082
5 61.3812216 -31.4093335
6 32.2016904 61.3812216
7 71.3927402 32.2016904
8 -82.8644010 71.3927402
9 -35.4826807 -82.8644010
10 39.1272019 -35.4826807
11 7.2988642 39.1272019
12 164.2544634 7.2988642
13 18.9965365 164.2544634
14 130.4691155 18.9965365
15 171.9435954 130.4691155
16 203.0453057 171.9435954
17 39.2772117 203.0453057
18 -40.1469609 39.2772117
19 -209.5308845 -40.1469609
20 -199.5308845 -209.5308845
21 -28.4929851 -199.5308845
22 -28.4929851 -28.4929851
23 63.9184575 -28.4929851
24 -20.4393993 63.9184575
25 -20.4393993 -20.4393993
26 -165.0128828 -20.4393993
27 28.8051789 -165.0128828
28 39.7309906 28.8051789
29 11.3927402 39.7309906
30 34.6226409 11.3927402
31 104.6226409 34.6226409
32 132.7753914 104.6226409
33 -34.3316859 132.7753914
34 -33.7993933 -34.3316859
35 -33.7468090 -33.7993933
36 0.3664198 -33.7468090
37 -21.7540853 0.3664198
38 32.7519745 -21.7540853
39 11.3297646 32.7519745
40 68.9480649 11.3297646
41 168.6987670 68.9480649
42 47.9062076 168.6987670
43 -1.5384237 47.9062076
44 -68.9056215 -1.5384237
45 9.7090768 -68.9056215
46 22.6874069 9.7090768
47 -40.1833251 22.6874069
48 -5.1671639 -40.1833251
49 -57.5055265 -5.1671639
50 -17.2926520 -57.5055265
51 -4.2197705 -17.2926520
52 -235.5408436 -4.2197705
53 26.1851476 -235.5408436
54 -10.6371869 26.1851476
55 78.4657958 -10.6371869
56 -55.6899805 78.4657958
57 -130.7518333 -55.6899805
58 -37.9638518 -130.7518333
59 -105.2855304 -37.9638518
> 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/wessaorg/rcomp/tmp/7n9rk1323898689.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/wessaorg/rcomp/tmp/8mzmk1323898689.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/wessaorg/rcomp/tmp/99rut1323898689.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/wessaorg/rcomp/tmp/109jr11323898689.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/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/wessaorg/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/wessaorg/rcomp/tmp/11a9ud1323898689.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/wessaorg/rcomp/tmp/12jjjm1323898689.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/wessaorg/rcomp/tmp/13syqf1323898689.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/wessaorg/rcomp/tmp/14yucl1323898689.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/wessaorg/rcomp/tmp/15446q1323898689.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/wessaorg/rcomp/tmp/163sbd1323898689.tab")
+ }
>
> try(system("convert tmp/1xa4a1323898689.ps tmp/1xa4a1323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/2jhox1323898689.ps tmp/2jhox1323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/3sris1323898689.ps tmp/3sris1323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/429p51323898689.ps tmp/429p51323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/538f31323898689.ps tmp/538f31323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/6mz4x1323898689.ps tmp/6mz4x1323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/7n9rk1323898689.ps tmp/7n9rk1323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/8mzmk1323898689.ps tmp/8mzmk1323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/99rut1323898689.ps tmp/99rut1323898689.png",intern=TRUE))
character(0)
> try(system("convert tmp/109jr11323898689.ps tmp/109jr11323898689.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.163 0.642 3.955