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.5,5.569,1.933,0.226,9.6,5.634,1.947,0.231,9.4,5.433,1.936,0.225,9.4,5.425,1.956,0.229,9.5,5.412,1.965,0.236,9.4,5.247,1.973,0.234,9.7,5.31,1.988,0.253,9.5,5.168,1.985,0.251,9.5,4.927,1.986,0.243,9.3,4.929,1.993,0.239,9.4,4.902,2.003,0.237,9.3,4.82,2,0.23,9.1,4.588,2.015,0.221,8.8,4.312,2.001,0.203,8.8,4.269,2.025,0.195,8.6,4.137,2.035,0.182,8.7,4.099,2.049,0.183,8.5,4.016,2.04,0.175,8.7,4.121,2.079,0.181,8.6,3.97,2.064,0.176,8.5,3.89,2.083,0.172,8.6,3.889,2.091,0.176,8.6,3.788,2.108,0.172,8.7,3.75,2.113,0.174,8.7,3.651,2.115,0.172,8.7,3.559,2.117,0.174,8.8,3.525,2.125,0.18,8.7,3.32,2.142,0.205,8.6,3.218,2.16,0.207,8.5,3.138,2.158,0.207,8.5,3.061,2.143,0.208,8.8,3.099,2.146,0.22,8.8,2.997,2.131,0.227,8.8,2.963,2.117,0.234,8.8,2.883,2.087,0.24,8.6,2.804,2.057,0.24,8.6,2.724,2.024,0.242,8.8,2.678,2.027,0.252,8.7,2.576,1.996,0.25,8.5,2.478,1.96,0.253),dim=c(4,40),dimnames=list(c('Rate','Heart_disease','Cancer','Diabetes
'),1:40))
> y <- array(NA,dim=c(4,40),dimnames=list(c('Rate','Heart_disease','Cancer','Diabetes
'),1:40))
> 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
Rate Heart_disease Cancer Diabetes\r
1 9.5 5.569 1.933 0.226
2 9.6 5.634 1.947 0.231
3 9.4 5.433 1.936 0.225
4 9.4 5.425 1.956 0.229
5 9.5 5.412 1.965 0.236
6 9.4 5.247 1.973 0.234
7 9.7 5.310 1.988 0.253
8 9.5 5.168 1.985 0.251
9 9.5 4.927 1.986 0.243
10 9.3 4.929 1.993 0.239
11 9.4 4.902 2.003 0.237
12 9.3 4.820 2.000 0.230
13 9.1 4.588 2.015 0.221
14 8.8 4.312 2.001 0.203
15 8.8 4.269 2.025 0.195
16 8.6 4.137 2.035 0.182
17 8.7 4.099 2.049 0.183
18 8.5 4.016 2.040 0.175
19 8.7 4.121 2.079 0.181
20 8.6 3.970 2.064 0.176
21 8.5 3.890 2.083 0.172
22 8.6 3.889 2.091 0.176
23 8.6 3.788 2.108 0.172
24 8.7 3.750 2.113 0.174
25 8.7 3.651 2.115 0.172
26 8.7 3.559 2.117 0.174
27 8.8 3.525 2.125 0.180
28 8.7 3.320 2.142 0.205
29 8.6 3.218 2.160 0.207
30 8.5 3.138 2.158 0.207
31 8.5 3.061 2.143 0.208
32 8.8 3.099 2.146 0.220
33 8.8 2.997 2.131 0.227
34 8.8 2.963 2.117 0.234
35 8.8 2.883 2.087 0.240
36 8.6 2.804 2.057 0.240
37 8.6 2.724 2.024 0.242
38 8.8 2.678 2.027 0.252
39 8.7 2.576 1.996 0.250
40 8.5 2.478 1.960 0.253
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Heart_disease Cancer `Diabetes\r`
3.1089 0.3709 1.2881 7.8710
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.18177 -0.04876 -0.01092 0.05011 0.22971
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.1089 0.9220 3.372 0.00180 **
Heart_disease 0.3709 0.0224 16.557 < 2e-16 ***
Cancer 1.2881 0.3704 3.478 0.00134 **
`Diabetes\r` 7.8710 0.6684 11.777 6.62e-14 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.09183 on 36 degrees of freedom
Multiple R-squared: 0.9468, Adjusted R-squared: 0.9424
F-statistic: 213.7 on 3 and 36 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.015277150 0.030554299 0.9847229
[2,] 0.027657851 0.055315702 0.9723421
[3,] 0.213974609 0.427949218 0.7860254
[4,] 0.140894687 0.281789375 0.8591053
[5,] 0.147055991 0.294111982 0.8529440
[6,] 0.110048236 0.220096471 0.8899518
[7,] 0.070207517 0.140415034 0.9297925
[8,] 0.039072702 0.078145404 0.9609273
[9,] 0.024982847 0.049965693 0.9750172
[10,] 0.013168659 0.026337318 0.9868313
[11,] 0.009358162 0.018716325 0.9906418
[12,] 0.006372174 0.012744347 0.9936278
[13,] 0.002859111 0.005718223 0.9971409
[14,] 0.001573041 0.003146083 0.9984270
[15,] 0.001856229 0.003712459 0.9981438
[16,] 0.003253268 0.006506536 0.9967467
[17,] 0.010483116 0.020966233 0.9895169
[18,] 0.071231197 0.142462394 0.9287688
[19,] 0.165957816 0.331915631 0.8340422
[20,] 0.222859290 0.445718580 0.7771407
[21,] 0.710312858 0.579374285 0.2896871
[22,] 0.756633919 0.486732162 0.2433661
[23,] 0.776315112 0.447369777 0.2236849
[24,] 0.756527307 0.486945387 0.2434727
[25,] 0.699006030 0.601987940 0.3009940
[26,] 0.788662687 0.422674625 0.2113373
[27,] 0.707068951 0.585862098 0.2929310
> postscript(file="/var/wessaorg/rcomp/tmp/1owxa1321943883.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/2l0sk1321943883.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/3b09i1321943883.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/412dr1321943883.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/5wp1m1321943883.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 = 40
Frequency = 1
1 2 3 4 5 6
0.056795070 0.075297271 0.011246780 -0.043031405 -0.004899199 -0.038260217
7 8 9 10 11 12
0.069501626 -0.058221695 0.092849675 -0.085424653 0.027451416 0.016828024
13 14 15 16 17 18
-0.045600914 -0.083516184 -0.035512416 -0.097108802 -0.008917955 -0.103570970
19 20 21 22 23 24
-0.039978395 -0.025293523 -0.088609427 -0.030027143 0.017022393 0.108934907
25 26 27 28 29 30
0.158821713 0.174628055 0.229708635 -0.012925547 -0.114019219 -0.181769583
31 32 33 34 35 36
-0.141758737 0.045830029 0.047887788 0.023435009 0.044524692 -0.087530492
37 38 39 40
-0.031092550 0.103395393 0.096901429 -0.043990879
> postscript(file="/var/wessaorg/rcomp/tmp/601951321943883.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 = 40
Frequency = 1
lag(myerror, k = 1) myerror
0 0.056795070 NA
1 0.075297271 0.056795070
2 0.011246780 0.075297271
3 -0.043031405 0.011246780
4 -0.004899199 -0.043031405
5 -0.038260217 -0.004899199
6 0.069501626 -0.038260217
7 -0.058221695 0.069501626
8 0.092849675 -0.058221695
9 -0.085424653 0.092849675
10 0.027451416 -0.085424653
11 0.016828024 0.027451416
12 -0.045600914 0.016828024
13 -0.083516184 -0.045600914
14 -0.035512416 -0.083516184
15 -0.097108802 -0.035512416
16 -0.008917955 -0.097108802
17 -0.103570970 -0.008917955
18 -0.039978395 -0.103570970
19 -0.025293523 -0.039978395
20 -0.088609427 -0.025293523
21 -0.030027143 -0.088609427
22 0.017022393 -0.030027143
23 0.108934907 0.017022393
24 0.158821713 0.108934907
25 0.174628055 0.158821713
26 0.229708635 0.174628055
27 -0.012925547 0.229708635
28 -0.114019219 -0.012925547
29 -0.181769583 -0.114019219
30 -0.141758737 -0.181769583
31 0.045830029 -0.141758737
32 0.047887788 0.045830029
33 0.023435009 0.047887788
34 0.044524692 0.023435009
35 -0.087530492 0.044524692
36 -0.031092550 -0.087530492
37 0.103395393 -0.031092550
38 0.096901429 0.103395393
39 -0.043990879 0.096901429
40 NA -0.043990879
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 0.075297271 0.056795070
[2,] 0.011246780 0.075297271
[3,] -0.043031405 0.011246780
[4,] -0.004899199 -0.043031405
[5,] -0.038260217 -0.004899199
[6,] 0.069501626 -0.038260217
[7,] -0.058221695 0.069501626
[8,] 0.092849675 -0.058221695
[9,] -0.085424653 0.092849675
[10,] 0.027451416 -0.085424653
[11,] 0.016828024 0.027451416
[12,] -0.045600914 0.016828024
[13,] -0.083516184 -0.045600914
[14,] -0.035512416 -0.083516184
[15,] -0.097108802 -0.035512416
[16,] -0.008917955 -0.097108802
[17,] -0.103570970 -0.008917955
[18,] -0.039978395 -0.103570970
[19,] -0.025293523 -0.039978395
[20,] -0.088609427 -0.025293523
[21,] -0.030027143 -0.088609427
[22,] 0.017022393 -0.030027143
[23,] 0.108934907 0.017022393
[24,] 0.158821713 0.108934907
[25,] 0.174628055 0.158821713
[26,] 0.229708635 0.174628055
[27,] -0.012925547 0.229708635
[28,] -0.114019219 -0.012925547
[29,] -0.181769583 -0.114019219
[30,] -0.141758737 -0.181769583
[31,] 0.045830029 -0.141758737
[32,] 0.047887788 0.045830029
[33,] 0.023435009 0.047887788
[34,] 0.044524692 0.023435009
[35,] -0.087530492 0.044524692
[36,] -0.031092550 -0.087530492
[37,] 0.103395393 -0.031092550
[38,] 0.096901429 0.103395393
[39,] -0.043990879 0.096901429
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 0.075297271 0.056795070
2 0.011246780 0.075297271
3 -0.043031405 0.011246780
4 -0.004899199 -0.043031405
5 -0.038260217 -0.004899199
6 0.069501626 -0.038260217
7 -0.058221695 0.069501626
8 0.092849675 -0.058221695
9 -0.085424653 0.092849675
10 0.027451416 -0.085424653
11 0.016828024 0.027451416
12 -0.045600914 0.016828024
13 -0.083516184 -0.045600914
14 -0.035512416 -0.083516184
15 -0.097108802 -0.035512416
16 -0.008917955 -0.097108802
17 -0.103570970 -0.008917955
18 -0.039978395 -0.103570970
19 -0.025293523 -0.039978395
20 -0.088609427 -0.025293523
21 -0.030027143 -0.088609427
22 0.017022393 -0.030027143
23 0.108934907 0.017022393
24 0.158821713 0.108934907
25 0.174628055 0.158821713
26 0.229708635 0.174628055
27 -0.012925547 0.229708635
28 -0.114019219 -0.012925547
29 -0.181769583 -0.114019219
30 -0.141758737 -0.181769583
31 0.045830029 -0.141758737
32 0.047887788 0.045830029
33 0.023435009 0.047887788
34 0.044524692 0.023435009
35 -0.087530492 0.044524692
36 -0.031092550 -0.087530492
37 0.103395393 -0.031092550
38 0.096901429 0.103395393
39 -0.043990879 0.096901429
> 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/7d3pb1321943883.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/8z72n1321943883.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/930xo1321943883.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/100q9c1321943883.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/11qjvs1321943883.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/12l9ji1321943883.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/13yl5z1321943883.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/14mf4z1321943883.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/15jtbi1321943884.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/161e641321943884.tab")
+ }
>
> try(system("convert tmp/1owxa1321943883.ps tmp/1owxa1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/2l0sk1321943883.ps tmp/2l0sk1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/3b09i1321943883.ps tmp/3b09i1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/412dr1321943883.ps tmp/412dr1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/5wp1m1321943883.ps tmp/5wp1m1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/601951321943883.ps tmp/601951321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/7d3pb1321943883.ps tmp/7d3pb1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/8z72n1321943883.ps tmp/8z72n1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/930xo1321943883.ps tmp/930xo1321943883.png",intern=TRUE))
character(0)
> try(system("convert tmp/100q9c1321943883.ps tmp/100q9c1321943883.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.802 0.486 3.375