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(9,1167,333,70,10,669,223,44,11,1053,371,35,10,1939,873,119,10,678,186,30,10,321,111,23,9,2667,1277,46,9,345,102,39,11,1367,580,58,11,1158,420,51,9,1385,521,65,9,1155,358,40,10,1120,435,41,10,1703,690,76,11,1189,393,31,11,3083,1149,82,9,1357,486,36,10,1892,767,62,10,883,338,28,10,1627,485,38,11,1412,465,70,9,1900,816,76,9,777,265,33,9,904,307,40,10,2115,850,126,10,1858,704,56,11,1781,693,63,11,1286,387,46,11,1035,406,35,9,1557,573,108,11,1527,595,34,11,1220,394,54,11,1368,521,35,10,564,172,23,10,1990,835,46,10,1557,669,49,10,2057,749,56,9,1111,368,38,11,686,216,19,9,2011,772,29,9,2232,1084,26,9,1032,445,52,10,1166,451,54,10,1020,300,45,10,1735,836,56,10,3623,1417,596,11,918,330,57,9,1579,477,55,10,2790,1028,99,9,1496,646,51,11,1108,342,21,11,496,218,20,9,1750,591,58,10,744,255,21,10,1101,434,66,11,1612,654,47,9,1805,478,55,10,2460,753,158,10,1653,689,46,11,1234,470,45),dim=c(4,60),dimnames=list(c('Month','TotalNrPV','TotalNrCC','TotalNrPRV'),1:60))
> y <- array(NA,dim=c(4,60),dimnames=list(c('Month','TotalNrPV','TotalNrCC','TotalNrPRV'),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 = '2'
> 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
TotalNrPV Month TotalNrCC TotalNrPRV
1 1167 9 333 70
2 669 10 223 44
3 1053 11 371 35
4 1939 10 873 119
5 678 10 186 30
6 321 10 111 23
7 2667 9 1277 46
8 345 9 102 39
9 1367 11 580 58
10 1158 11 420 51
11 1385 9 521 65
12 1155 9 358 40
13 1120 10 435 41
14 1703 10 690 76
15 1189 11 393 31
16 3083 11 1149 82
17 1357 9 486 36
18 1892 10 767 62
19 883 10 338 28
20 1627 10 485 38
21 1412 11 465 70
22 1900 9 816 76
23 777 9 265 33
24 904 9 307 40
25 2115 10 850 126
26 1858 10 704 56
27 1781 11 693 63
28 1286 11 387 46
29 1035 11 406 35
30 1557 9 573 108
31 1527 11 595 34
32 1220 11 394 54
33 1368 11 521 35
34 564 10 172 23
35 1990 10 835 46
36 1557 10 669 49
37 2057 10 749 56
38 1111 9 368 38
39 686 11 216 19
40 2011 9 772 29
41 2232 9 1084 26
42 1032 9 445 52
43 1166 10 451 54
44 1020 10 300 45
45 1735 10 836 56
46 3623 10 1417 596
47 918 11 330 57
48 1579 9 477 55
49 2790 10 1028 99
50 1496 9 646 51
51 1108 11 342 21
52 496 11 218 20
53 1750 9 591 58
54 744 10 255 21
55 1101 10 434 66
56 1612 11 654 47
57 1805 9 478 55
58 2460 10 753 158
59 1653 10 689 46
60 1234 11 470 45
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Month TotalNrCC TotalNrPRV
237.4614 2.5085 2.0661 0.9802
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-309.67 -103.72 -28.37 89.37 503.47
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 237.46144 310.49787 0.765 0.4476
Month 2.50849 30.14381 0.083 0.9340
TotalNrCC 2.06607 0.09968 20.727 <2e-16 ***
TotalNrPRV 0.98019 0.37110 2.641 0.0107 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 179.9 on 56 degrees of freedom
Multiple R-squared: 0.9266, Adjusted R-squared: 0.9226
F-statistic: 235.5 on 3 and 56 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.503365434 0.993269132 0.4966346
[2,] 0.472862677 0.945725354 0.5271373
[3,] 0.333191709 0.666383418 0.6668083
[4,] 0.236638007 0.473276014 0.7633620
[5,] 0.168592419 0.337184839 0.8314076
[6,] 0.179121202 0.358242403 0.8208788
[7,] 0.113941051 0.227882103 0.8860589
[8,] 0.073548518 0.147097036 0.9264515
[9,] 0.074661910 0.149323820 0.9253381
[10,] 0.454484262 0.908968524 0.5455157
[11,] 0.400759230 0.801518459 0.5992408
[12,] 0.314653814 0.629307629 0.6853462
[13,] 0.254696789 0.509393577 0.7453032
[14,] 0.453943791 0.907887581 0.5460562
[15,] 0.397419331 0.794838662 0.6025807
[16,] 0.339165615 0.678331230 0.6608344
[17,] 0.274506500 0.549012999 0.7254935
[18,] 0.215980017 0.431960034 0.7840200
[19,] 0.160711814 0.321423628 0.8392882
[20,] 0.126283865 0.252567729 0.8737161
[21,] 0.090151811 0.180303621 0.9098482
[22,] 0.083992371 0.167984742 0.9160076
[23,] 0.069775893 0.139551785 0.9302241
[24,] 0.048810725 0.097621449 0.9511893
[25,] 0.032324973 0.064649946 0.9676750
[26,] 0.022782092 0.045564184 0.9772179
[27,] 0.014392298 0.028784597 0.9856077
[28,] 0.009630086 0.019260171 0.9903699
[29,] 0.005575960 0.011151920 0.9944240
[30,] 0.004124909 0.008249818 0.9958751
[31,] 0.005025276 0.010050553 0.9949747
[32,] 0.003195557 0.006391113 0.9968044
[33,] 0.001776417 0.003552834 0.9982236
[34,] 0.001322443 0.002644886 0.9986776
[35,] 0.003679398 0.007358796 0.9963206
[36,] 0.006731493 0.013462986 0.9932685
[37,] 0.004474756 0.008949513 0.9955252
[38,] 0.002610228 0.005220457 0.9973898
[39,] 0.015446846 0.030893691 0.9845532
[40,] 0.095339704 0.190679408 0.9046603
[41,] 0.069487910 0.138975820 0.9305121
[42,] 0.077666690 0.155333380 0.9223333
[43,] 0.085966753 0.171933507 0.9140332
[44,] 0.153172062 0.306344123 0.8468279
[45,] 0.294588330 0.589176660 0.7054117
[46,] 0.210096665 0.420193329 0.7899033
[47,] 0.145477219 0.290954439 0.8545228
> postscript(file="/var/wessaorg/rcomp/tmp/15tlr1321899219.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/2d6xm1321899219.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/3qrgr1321899219.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/4xhja1321899219.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/5gmwx1321899219.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
150.3464676 -97.4090200 -12.8744447 -243.8705511 1.7583770 -193.4248452
7 8 9 10 11 12
-276.5010158 -164.0048190 -153.2280122 -24.7950975 -15.1741323 116.1005034
13 14 15 16 17 18
-81.4757421 -59.6309670 81.5927457 363.6522340 57.5640406 -15.9958007
19 20 21 22 23 24
-105.3242078 325.1612320 117.6079544 -120.4475727 -62.8934193 -29.5298145
25 26 27 28 29 30
-27.2122534 86.0479160 22.4048559 176.2862585 -103.1869717 7.2417417
31 32 33 34 35 36
-0.6944226 87.9821956 -7.7852746 -76.4552493 -42.8055952 -135.7781942
37 38 39 40 41 42
192.0746670 53.4001707 -43.9501388 127.5287545 -293.1451875 -198.4101142
43 44 45 46 47 48
-81.2754283 93.5232260 -309.6736143 -151.3666944 -84.7297677 279.5349913
49 50 51 52 53 54
306.4921518 -148.7104316 115.7643747 -239.0624779 212.0621766 -65.9788525
55 56 57 58 59 60
-122.9145372 -50.3352133 503.4689191 486.8305199 -78.1590541 -46.2175393
> postscript(file="/var/wessaorg/rcomp/tmp/6wkrc1321899219.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 150.3464676 NA
1 -97.4090200 150.3464676
2 -12.8744447 -97.4090200
3 -243.8705511 -12.8744447
4 1.7583770 -243.8705511
5 -193.4248452 1.7583770
6 -276.5010158 -193.4248452
7 -164.0048190 -276.5010158
8 -153.2280122 -164.0048190
9 -24.7950975 -153.2280122
10 -15.1741323 -24.7950975
11 116.1005034 -15.1741323
12 -81.4757421 116.1005034
13 -59.6309670 -81.4757421
14 81.5927457 -59.6309670
15 363.6522340 81.5927457
16 57.5640406 363.6522340
17 -15.9958007 57.5640406
18 -105.3242078 -15.9958007
19 325.1612320 -105.3242078
20 117.6079544 325.1612320
21 -120.4475727 117.6079544
22 -62.8934193 -120.4475727
23 -29.5298145 -62.8934193
24 -27.2122534 -29.5298145
25 86.0479160 -27.2122534
26 22.4048559 86.0479160
27 176.2862585 22.4048559
28 -103.1869717 176.2862585
29 7.2417417 -103.1869717
30 -0.6944226 7.2417417
31 87.9821956 -0.6944226
32 -7.7852746 87.9821956
33 -76.4552493 -7.7852746
34 -42.8055952 -76.4552493
35 -135.7781942 -42.8055952
36 192.0746670 -135.7781942
37 53.4001707 192.0746670
38 -43.9501388 53.4001707
39 127.5287545 -43.9501388
40 -293.1451875 127.5287545
41 -198.4101142 -293.1451875
42 -81.2754283 -198.4101142
43 93.5232260 -81.2754283
44 -309.6736143 93.5232260
45 -151.3666944 -309.6736143
46 -84.7297677 -151.3666944
47 279.5349913 -84.7297677
48 306.4921518 279.5349913
49 -148.7104316 306.4921518
50 115.7643747 -148.7104316
51 -239.0624779 115.7643747
52 212.0621766 -239.0624779
53 -65.9788525 212.0621766
54 -122.9145372 -65.9788525
55 -50.3352133 -122.9145372
56 503.4689191 -50.3352133
57 486.8305199 503.4689191
58 -78.1590541 486.8305199
59 -46.2175393 -78.1590541
60 NA -46.2175393
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -97.4090200 150.3464676
[2,] -12.8744447 -97.4090200
[3,] -243.8705511 -12.8744447
[4,] 1.7583770 -243.8705511
[5,] -193.4248452 1.7583770
[6,] -276.5010158 -193.4248452
[7,] -164.0048190 -276.5010158
[8,] -153.2280122 -164.0048190
[9,] -24.7950975 -153.2280122
[10,] -15.1741323 -24.7950975
[11,] 116.1005034 -15.1741323
[12,] -81.4757421 116.1005034
[13,] -59.6309670 -81.4757421
[14,] 81.5927457 -59.6309670
[15,] 363.6522340 81.5927457
[16,] 57.5640406 363.6522340
[17,] -15.9958007 57.5640406
[18,] -105.3242078 -15.9958007
[19,] 325.1612320 -105.3242078
[20,] 117.6079544 325.1612320
[21,] -120.4475727 117.6079544
[22,] -62.8934193 -120.4475727
[23,] -29.5298145 -62.8934193
[24,] -27.2122534 -29.5298145
[25,] 86.0479160 -27.2122534
[26,] 22.4048559 86.0479160
[27,] 176.2862585 22.4048559
[28,] -103.1869717 176.2862585
[29,] 7.2417417 -103.1869717
[30,] -0.6944226 7.2417417
[31,] 87.9821956 -0.6944226
[32,] -7.7852746 87.9821956
[33,] -76.4552493 -7.7852746
[34,] -42.8055952 -76.4552493
[35,] -135.7781942 -42.8055952
[36,] 192.0746670 -135.7781942
[37,] 53.4001707 192.0746670
[38,] -43.9501388 53.4001707
[39,] 127.5287545 -43.9501388
[40,] -293.1451875 127.5287545
[41,] -198.4101142 -293.1451875
[42,] -81.2754283 -198.4101142
[43,] 93.5232260 -81.2754283
[44,] -309.6736143 93.5232260
[45,] -151.3666944 -309.6736143
[46,] -84.7297677 -151.3666944
[47,] 279.5349913 -84.7297677
[48,] 306.4921518 279.5349913
[49,] -148.7104316 306.4921518
[50,] 115.7643747 -148.7104316
[51,] -239.0624779 115.7643747
[52,] 212.0621766 -239.0624779
[53,] -65.9788525 212.0621766
[54,] -122.9145372 -65.9788525
[55,] -50.3352133 -122.9145372
[56,] 503.4689191 -50.3352133
[57,] 486.8305199 503.4689191
[58,] -78.1590541 486.8305199
[59,] -46.2175393 -78.1590541
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -97.4090200 150.3464676
2 -12.8744447 -97.4090200
3 -243.8705511 -12.8744447
4 1.7583770 -243.8705511
5 -193.4248452 1.7583770
6 -276.5010158 -193.4248452
7 -164.0048190 -276.5010158
8 -153.2280122 -164.0048190
9 -24.7950975 -153.2280122
10 -15.1741323 -24.7950975
11 116.1005034 -15.1741323
12 -81.4757421 116.1005034
13 -59.6309670 -81.4757421
14 81.5927457 -59.6309670
15 363.6522340 81.5927457
16 57.5640406 363.6522340
17 -15.9958007 57.5640406
18 -105.3242078 -15.9958007
19 325.1612320 -105.3242078
20 117.6079544 325.1612320
21 -120.4475727 117.6079544
22 -62.8934193 -120.4475727
23 -29.5298145 -62.8934193
24 -27.2122534 -29.5298145
25 86.0479160 -27.2122534
26 22.4048559 86.0479160
27 176.2862585 22.4048559
28 -103.1869717 176.2862585
29 7.2417417 -103.1869717
30 -0.6944226 7.2417417
31 87.9821956 -0.6944226
32 -7.7852746 87.9821956
33 -76.4552493 -7.7852746
34 -42.8055952 -76.4552493
35 -135.7781942 -42.8055952
36 192.0746670 -135.7781942
37 53.4001707 192.0746670
38 -43.9501388 53.4001707
39 127.5287545 -43.9501388
40 -293.1451875 127.5287545
41 -198.4101142 -293.1451875
42 -81.2754283 -198.4101142
43 93.5232260 -81.2754283
44 -309.6736143 93.5232260
45 -151.3666944 -309.6736143
46 -84.7297677 -151.3666944
47 279.5349913 -84.7297677
48 306.4921518 279.5349913
49 -148.7104316 306.4921518
50 115.7643747 -148.7104316
51 -239.0624779 115.7643747
52 212.0621766 -239.0624779
53 -65.9788525 212.0621766
54 -122.9145372 -65.9788525
55 -50.3352133 -122.9145372
56 503.4689191 -50.3352133
57 486.8305199 503.4689191
58 -78.1590541 486.8305199
59 -46.2175393 -78.1590541
> 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/7o53c1321899219.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/8tohz1321899219.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/9h2qv1321899219.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/109oiz1321899219.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/11gzrd1321899219.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/12gaox1321899219.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/13xkuu1321899219.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/144l6t1321899220.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/15r95k1321899220.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/16ymjg1321899220.tab")
+ }
>
> try(system("convert tmp/15tlr1321899219.ps tmp/15tlr1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/2d6xm1321899219.ps tmp/2d6xm1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/3qrgr1321899219.ps tmp/3qrgr1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/4xhja1321899219.ps tmp/4xhja1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/5gmwx1321899219.ps tmp/5gmwx1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/6wkrc1321899219.ps tmp/6wkrc1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/7o53c1321899219.ps tmp/7o53c1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/8tohz1321899219.ps tmp/8tohz1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/9h2qv1321899219.ps tmp/9h2qv1321899219.png",intern=TRUE))
character(0)
> try(system("convert tmp/109oiz1321899219.ps tmp/109oiz1321899219.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.157 0.511 3.692