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(210907,0,2,149061,0,0,237213,1,0,133131,1,4,324799,1,0,230964,0,-1,236785,1,0,344297,1,1,174724,1,0,174415,1,3,223632,1,-1,294424,0,4,325107,1,3,106408,0,1,96560,0,0,265769,1,-2,149112,0,-4,152871,0,2,362301,1,2,183167,0,-4,218946,1,2,244052,1,2,341570,1,0,196553,1,-3,143246,0,2,167488,0,0,143756,0,4,152299,1,2,193339,1,2,130585,0,-4,112611,1,3,148446,1,3,182079,0,2,243060,1,-1,162765,1,-3,85574,1,0,225060,0,1,133328,1,-3,100750,1,3,101523,1,0,243511,1,0,152474,1,0,132487,1,3,317394,0,-3,244749,1,0,128423,0,2,97839,0,-1,229242,1,2,324598,0,2,195838,0,-2,254488,0,0,92499,1,-2,224330,0,0,181633,1,6,271856,1,-3,95227,1,3,98146,0,0,118612,0,-2,65475,1,1,108446,0,0,121848,0,2,76302,1,2,98104,0,-3,30989,1,-2,31774,0,1,150580,1,-4,59382,0,1,84105,0,0),dim=c(3,68),dimnames=list(c('RFCseconds','Gender','Testscore'),1:68))
> y <- array(NA,dim=c(3,68),dimnames=list(c('RFCseconds','Gender','Testscore'),1:68))
> 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
> 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
RFCseconds Gender Testscore
1 210907 0 2
2 149061 0 0
3 237213 1 0
4 133131 1 4
5 324799 1 0
6 230964 0 -1
7 236785 1 0
8 344297 1 1
9 174724 1 0
10 174415 1 3
11 223632 1 -1
12 294424 0 4
13 325107 1 3
14 106408 0 1
15 96560 0 0
16 265769 1 -2
17 149112 0 -4
18 152871 0 2
19 362301 1 2
20 183167 0 -4
21 218946 1 2
22 244052 1 2
23 341570 1 0
24 196553 1 -3
25 143246 0 2
26 167488 0 0
27 143756 0 4
28 152299 1 2
29 193339 1 2
30 130585 0 -4
31 112611 1 3
32 148446 1 3
33 182079 0 2
34 243060 1 -1
35 162765 1 -3
36 85574 1 0
37 225060 0 1
38 133328 1 -3
39 100750 1 3
40 101523 1 0
41 243511 1 0
42 152474 1 0
43 132487 1 3
44 317394 0 -3
45 244749 1 0
46 128423 0 2
47 97839 0 -1
48 229242 1 2
49 324598 0 2
50 195838 0 -2
51 254488 0 0
52 92499 1 -2
53 224330 0 0
54 181633 1 6
55 271856 1 -3
56 95227 1 3
57 98146 0 0
58 118612 0 -2
59 65475 1 1
60 108446 0 0
61 121848 0 2
62 76302 1 2
63 98104 0 -3
64 30989 1 -2
65 31774 0 1
66 150580 1 -4
67 59382 0 1
68 84105 0 0
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Gender Testscore
160952.0 28489.8 227.1
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-157999 -57324 -13304 53901 172405
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 160952.0 14738.2 10.921 2.39e-16 ***
Gender 28489.8 19853.3 1.435 0.156
Testscore 227.1 4396.8 0.052 0.959
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 80710 on 65 degrees of freedom
Multiple R-squared: 0.03144, Adjusted R-squared: 0.001634
F-statistic: 1.055 on 2 and 65 DF, p-value: 0.3541
> 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.3719035 0.74380705 0.628096475
[2,] 0.2324774 0.46495483 0.767522587
[3,] 0.3719094 0.74381879 0.628090606
[4,] 0.4225187 0.84503736 0.577481319
[5,] 0.3114723 0.62294464 0.688527680
[6,] 0.2499454 0.49989088 0.750054558
[7,] 0.4019487 0.80389734 0.598051329
[8,] 0.4802090 0.96041794 0.519791028
[9,] 0.5264548 0.94709037 0.473545187
[10,] 0.5404069 0.91918624 0.459593122
[11,] 0.4818710 0.96374192 0.518129042
[12,] 0.3986982 0.79739647 0.601301766
[13,] 0.3239050 0.64781009 0.676094956
[14,] 0.4958353 0.99167058 0.504164710
[15,] 0.4201205 0.84024090 0.579879548
[16,] 0.3682967 0.73659348 0.631703260
[17,] 0.3227379 0.64547575 0.677262126
[18,] 0.4706803 0.94136067 0.529319665
[19,] 0.4270658 0.85413161 0.572934196
[20,] 0.3653381 0.73067627 0.634661867
[21,] 0.2968833 0.59376653 0.703116737
[22,] 0.2436989 0.48739788 0.756301060
[23,] 0.2602381 0.52047620 0.739761899
[24,] 0.2314668 0.46293350 0.768533248
[25,] 0.1903757 0.38075131 0.809624345
[26,] 0.2422070 0.48441399 0.757793004
[27,] 0.2287317 0.45746331 0.771268345
[28,] 0.1823613 0.36472254 0.817638731
[29,] 0.1667826 0.33356527 0.833217364
[30,] 0.1466958 0.29339169 0.853304153
[31,] 0.2001552 0.40031033 0.799844833
[32,] 0.1858514 0.37170281 0.814148595
[33,] 0.1677227 0.33544540 0.832277301
[34,] 0.1815252 0.36305047 0.818474767
[35,] 0.1866663 0.37333268 0.813333658
[36,] 0.1775249 0.35504971 0.822475145
[37,] 0.1429525 0.28590491 0.857047546
[38,] 0.1198308 0.23966153 0.880169236
[39,] 0.2594659 0.51893184 0.740534079
[40,] 0.2709649 0.54192985 0.729035076
[41,] 0.2192503 0.43850055 0.780749726
[42,] 0.1926596 0.38531919 0.807340406
[43,] 0.1949123 0.38982454 0.805087729
[44,] 0.5120502 0.97589965 0.487949826
[45,] 0.4769328 0.95386566 0.523067171
[46,] 0.6538336 0.69233285 0.346166423
[47,] 0.6249386 0.75012277 0.375061387
[48,] 0.7631388 0.47372236 0.236861182
[49,] 0.8453916 0.30921670 0.154608352
[50,] 0.9912686 0.01746276 0.008731379
[51,] 0.9878920 0.02421592 0.012107959
[52,] 0.9760507 0.04789863 0.023949314
[53,] 0.9550802 0.08983951 0.044919757
[54,] 0.9210871 0.15782587 0.078912933
[55,] 0.8679063 0.26418748 0.132093740
[56,] 0.8844378 0.23112435 0.115562175
[57,] 0.8552302 0.28953962 0.144769812
> postscript(file="/var/wessaorg/rcomp/tmp/1si8z1323789687.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/2mv9l1323789687.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/3y4qy1323789687.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/4vium1323789687.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/5fdna1323789687.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 = 68
Frequency = 1
1 2 3 4 5 6
49500.718 -11891.025 47771.183 -57219.332 135357.183 70239.104
7 8 9 10 11 12
47343.183 154628.054 -14717.817 -15708.203 34417.312 132563.461
13 14 15 16 17 18
134983.797 -54771.153 -64392.025 76781.440 -10931.510 -8535.282
19 20 21 22 23 24
172404.926 23123.490 29049.926 54155.926 152128.183 7792.569
25 26 27 28 29 30
-18160.282 6535.975 -18104.539 -37597.074 3442.926 -29458.510
31 32 33 34 35 36
-77512.203 -41677.203 20672.718 53845.312 -25995.431 -103867.817
37 38 39 40 41 42
63880.847 -55432.431 -89373.203 -87918.817 54069.183 -36967.817
43 44 45 46 47 48
-57636.203 157123.361 55307.183 -32983.282 -62885.896 39345.926
49 50 51 52 53 54
163191.718 35340.233 93535.975 -96488.560 63377.975 -9171.589
55 56 57 58 59 60
83095.569 -94896.203 -62806.025 -41885.767 -124193.946 -52506.025
61 62 63 64 65 66
-39558.282 -113594.074 -62166.639 -157998.560 -129405.153 -37953.302
67 68
-101797.153 -76847.025
> postscript(file="/var/wessaorg/rcomp/tmp/6gynk1323789688.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 = 68
Frequency = 1
lag(myerror, k = 1) myerror
0 49500.718 NA
1 -11891.025 49500.718
2 47771.183 -11891.025
3 -57219.332 47771.183
4 135357.183 -57219.332
5 70239.104 135357.183
6 47343.183 70239.104
7 154628.054 47343.183
8 -14717.817 154628.054
9 -15708.203 -14717.817
10 34417.312 -15708.203
11 132563.461 34417.312
12 134983.797 132563.461
13 -54771.153 134983.797
14 -64392.025 -54771.153
15 76781.440 -64392.025
16 -10931.510 76781.440
17 -8535.282 -10931.510
18 172404.926 -8535.282
19 23123.490 172404.926
20 29049.926 23123.490
21 54155.926 29049.926
22 152128.183 54155.926
23 7792.569 152128.183
24 -18160.282 7792.569
25 6535.975 -18160.282
26 -18104.539 6535.975
27 -37597.074 -18104.539
28 3442.926 -37597.074
29 -29458.510 3442.926
30 -77512.203 -29458.510
31 -41677.203 -77512.203
32 20672.718 -41677.203
33 53845.312 20672.718
34 -25995.431 53845.312
35 -103867.817 -25995.431
36 63880.847 -103867.817
37 -55432.431 63880.847
38 -89373.203 -55432.431
39 -87918.817 -89373.203
40 54069.183 -87918.817
41 -36967.817 54069.183
42 -57636.203 -36967.817
43 157123.361 -57636.203
44 55307.183 157123.361
45 -32983.282 55307.183
46 -62885.896 -32983.282
47 39345.926 -62885.896
48 163191.718 39345.926
49 35340.233 163191.718
50 93535.975 35340.233
51 -96488.560 93535.975
52 63377.975 -96488.560
53 -9171.589 63377.975
54 83095.569 -9171.589
55 -94896.203 83095.569
56 -62806.025 -94896.203
57 -41885.767 -62806.025
58 -124193.946 -41885.767
59 -52506.025 -124193.946
60 -39558.282 -52506.025
61 -113594.074 -39558.282
62 -62166.639 -113594.074
63 -157998.560 -62166.639
64 -129405.153 -157998.560
65 -37953.302 -129405.153
66 -101797.153 -37953.302
67 -76847.025 -101797.153
68 NA -76847.025
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -11891.025 49500.718
[2,] 47771.183 -11891.025
[3,] -57219.332 47771.183
[4,] 135357.183 -57219.332
[5,] 70239.104 135357.183
[6,] 47343.183 70239.104
[7,] 154628.054 47343.183
[8,] -14717.817 154628.054
[9,] -15708.203 -14717.817
[10,] 34417.312 -15708.203
[11,] 132563.461 34417.312
[12,] 134983.797 132563.461
[13,] -54771.153 134983.797
[14,] -64392.025 -54771.153
[15,] 76781.440 -64392.025
[16,] -10931.510 76781.440
[17,] -8535.282 -10931.510
[18,] 172404.926 -8535.282
[19,] 23123.490 172404.926
[20,] 29049.926 23123.490
[21,] 54155.926 29049.926
[22,] 152128.183 54155.926
[23,] 7792.569 152128.183
[24,] -18160.282 7792.569
[25,] 6535.975 -18160.282
[26,] -18104.539 6535.975
[27,] -37597.074 -18104.539
[28,] 3442.926 -37597.074
[29,] -29458.510 3442.926
[30,] -77512.203 -29458.510
[31,] -41677.203 -77512.203
[32,] 20672.718 -41677.203
[33,] 53845.312 20672.718
[34,] -25995.431 53845.312
[35,] -103867.817 -25995.431
[36,] 63880.847 -103867.817
[37,] -55432.431 63880.847
[38,] -89373.203 -55432.431
[39,] -87918.817 -89373.203
[40,] 54069.183 -87918.817
[41,] -36967.817 54069.183
[42,] -57636.203 -36967.817
[43,] 157123.361 -57636.203
[44,] 55307.183 157123.361
[45,] -32983.282 55307.183
[46,] -62885.896 -32983.282
[47,] 39345.926 -62885.896
[48,] 163191.718 39345.926
[49,] 35340.233 163191.718
[50,] 93535.975 35340.233
[51,] -96488.560 93535.975
[52,] 63377.975 -96488.560
[53,] -9171.589 63377.975
[54,] 83095.569 -9171.589
[55,] -94896.203 83095.569
[56,] -62806.025 -94896.203
[57,] -41885.767 -62806.025
[58,] -124193.946 -41885.767
[59,] -52506.025 -124193.946
[60,] -39558.282 -52506.025
[61,] -113594.074 -39558.282
[62,] -62166.639 -113594.074
[63,] -157998.560 -62166.639
[64,] -129405.153 -157998.560
[65,] -37953.302 -129405.153
[66,] -101797.153 -37953.302
[67,] -76847.025 -101797.153
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -11891.025 49500.718
2 47771.183 -11891.025
3 -57219.332 47771.183
4 135357.183 -57219.332
5 70239.104 135357.183
6 47343.183 70239.104
7 154628.054 47343.183
8 -14717.817 154628.054
9 -15708.203 -14717.817
10 34417.312 -15708.203
11 132563.461 34417.312
12 134983.797 132563.461
13 -54771.153 134983.797
14 -64392.025 -54771.153
15 76781.440 -64392.025
16 -10931.510 76781.440
17 -8535.282 -10931.510
18 172404.926 -8535.282
19 23123.490 172404.926
20 29049.926 23123.490
21 54155.926 29049.926
22 152128.183 54155.926
23 7792.569 152128.183
24 -18160.282 7792.569
25 6535.975 -18160.282
26 -18104.539 6535.975
27 -37597.074 -18104.539
28 3442.926 -37597.074
29 -29458.510 3442.926
30 -77512.203 -29458.510
31 -41677.203 -77512.203
32 20672.718 -41677.203
33 53845.312 20672.718
34 -25995.431 53845.312
35 -103867.817 -25995.431
36 63880.847 -103867.817
37 -55432.431 63880.847
38 -89373.203 -55432.431
39 -87918.817 -89373.203
40 54069.183 -87918.817
41 -36967.817 54069.183
42 -57636.203 -36967.817
43 157123.361 -57636.203
44 55307.183 157123.361
45 -32983.282 55307.183
46 -62885.896 -32983.282
47 39345.926 -62885.896
48 163191.718 39345.926
49 35340.233 163191.718
50 93535.975 35340.233
51 -96488.560 93535.975
52 63377.975 -96488.560
53 -9171.589 63377.975
54 83095.569 -9171.589
55 -94896.203 83095.569
56 -62806.025 -94896.203
57 -41885.767 -62806.025
58 -124193.946 -41885.767
59 -52506.025 -124193.946
60 -39558.282 -52506.025
61 -113594.074 -39558.282
62 -62166.639 -113594.074
63 -157998.560 -62166.639
64 -129405.153 -157998.560
65 -37953.302 -129405.153
66 -101797.153 -37953.302
67 -76847.025 -101797.153
> 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/7d7rl1323789688.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/89qi11323789688.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/99r1v1323789688.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/10mbhl1323789688.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/11f8ka1323789688.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/12qew61323789688.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/13skh91323789688.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/14j4tw1323789688.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/15ckim1323789688.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/16a1oj1323789688.tab")
+ }
>
> try(system("convert tmp/1si8z1323789687.ps tmp/1si8z1323789687.png",intern=TRUE))
character(0)
> try(system("convert tmp/2mv9l1323789687.ps tmp/2mv9l1323789687.png",intern=TRUE))
character(0)
> try(system("convert tmp/3y4qy1323789687.ps tmp/3y4qy1323789687.png",intern=TRUE))
character(0)
> try(system("convert tmp/4vium1323789687.ps tmp/4vium1323789687.png",intern=TRUE))
character(0)
> try(system("convert tmp/5fdna1323789687.ps tmp/5fdna1323789687.png",intern=TRUE))
character(0)
> try(system("convert tmp/6gynk1323789688.ps tmp/6gynk1323789688.png",intern=TRUE))
character(0)
> try(system("convert tmp/7d7rl1323789688.ps tmp/7d7rl1323789688.png",intern=TRUE))
character(0)
> try(system("convert tmp/89qi11323789688.ps tmp/89qi11323789688.png",intern=TRUE))
character(0)
> try(system("convert tmp/99r1v1323789688.ps tmp/99r1v1323789688.png",intern=TRUE))
character(0)
> try(system("convert tmp/10mbhl1323789688.ps tmp/10mbhl1323789688.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.665 0.645 4.909