R version 2.9.0 (2009-04-17)
Copyright (C) 2009 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> x <- array(list(269645,0,267037,0,258113,0,262813,0,267413,0,267366,0,264777,0,258863,0,254844,0,254868,0,277267,0,285351,0,286602,0,283042,0,276687,0,277915,0,277128,0,277103,0,275037,0,270150,0,267140,0,264993,0,287259,0,291186,0,292300,0,288186,0,281477,0,282656,0,280190,0,280408,0,276836,0,275216,0,274352,0,271311,0,289802,0,290726,0,292300,0,278506,0,269826,0,265861,0,269034,0,264176,0,255198,0,253353,0,246057,0,235372,0,258556,0,260993,0,254663,0,250643,0,243422,0,247105,0,248541,0,245039,0,237080,0,237085,0,225554,0,226839,0,247934,0,248333,1,246969,1,245098,1,246263,1,255765,1,264319,1,268347,1,273046,1,273963,1,267430,1,271993,1,292710,1,295881,1),dim=c(2,72),dimnames=list(c('Y','x'),1:72))
> y <- array(NA,dim=c(2,72),dimnames=list(c('Y','x'),1:72))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = 'No Linear Trend'
> par2 = 'Do not include Seasonal Dummies'
> par1 = '1'
> #'GNU S' R Code compiled by R2WASP v. 1.0.44 ()
> #Author: Prof. Dr. P. Wessa
> #To cite this work: AUTHOR(S), (YEAR), YOUR SOFTWARE TITLE (vNUMBER) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_YOURPAGE.wasp/
> #Source of accompanying publication: Office for Research, Development, and Education
> #Technical description: Write here your technical program description (don't use hard returns!)
> library(lattice)
> library(lmtest)
Loading required package: zoo
Attaching package: 'zoo'
The following object(s) are masked from package:base :
as.Date.numeric
> n25 <- 25 #minimum number of obs. for Goldfeld-Quandt test
> par1 <- as.numeric(par1)
> x <- t(y)
> k <- length(x[1,])
> n <- length(x[,1])
> x1 <- cbind(x[,par1], x[,1:k!=par1])
> mycolnames <- c(colnames(x)[par1], colnames(x)[1:k!=par1])
> colnames(x1) <- mycolnames #colnames(x)[par1]
> x <- x1
> if (par3 == 'First Differences'){
+ x2 <- array(0, dim=c(n-1,k), dimnames=list(1:(n-1), paste('(1-B)',colnames(x),sep='')))
+ for (i in 1:n-1) {
+ for (j in 1:k) {
+ x2[i,j] <- x[i+1,j] - x[i,j]
+ }
+ }
+ x <- x2
+ }
> if (par2 == 'Include Monthly Dummies'){
+ x2 <- array(0, dim=c(n,11), dimnames=list(1:n, paste('M', seq(1:11), sep ='')))
+ for (i in 1:11){
+ x2[seq(i,n,12),i] <- 1
+ }
+ x <- cbind(x, x2)
+ }
> if (par2 == 'Include Quarterly Dummies'){
+ x2 <- array(0, dim=c(n,3), dimnames=list(1:n, paste('Q', seq(1:3), sep ='')))
+ for (i in 1:3){
+ x2[seq(i,n,4),i] <- 1
+ }
+ x <- cbind(x, x2)
+ }
> k <- length(x[1,])
> if (par3 == 'Linear Trend'){
+ x <- cbind(x, c(1:n))
+ colnames(x)[k+1] <- 't'
+ }
> x
Y x
1 269645 0
2 267037 0
3 258113 0
4 262813 0
5 267413 0
6 267366 0
7 264777 0
8 258863 0
9 254844 0
10 254868 0
11 277267 0
12 285351 0
13 286602 0
14 283042 0
15 276687 0
16 277915 0
17 277128 0
18 277103 0
19 275037 0
20 270150 0
21 267140 0
22 264993 0
23 287259 0
24 291186 0
25 292300 0
26 288186 0
27 281477 0
28 282656 0
29 280190 0
30 280408 0
31 276836 0
32 275216 0
33 274352 0
34 271311 0
35 289802 0
36 290726 0
37 292300 0
38 278506 0
39 269826 0
40 265861 0
41 269034 0
42 264176 0
43 255198 0
44 253353 0
45 246057 0
46 235372 0
47 258556 0
48 260993 0
49 254663 0
50 250643 0
51 243422 0
52 247105 0
53 248541 0
54 245039 0
55 237080 0
56 237085 0
57 225554 0
58 226839 0
59 247934 0
60 248333 1
61 246969 1
62 245098 1
63 246263 1
64 255765 1
65 264319 1
66 268347 1
67 273046 1
68 273963 1
69 267430 1
70 271993 1
71 292710 1
72 295881 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) x
266427 -1033
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-40873 -11628 1511 11002 30487
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 266427 2198 121.2 <2e-16 ***
x -1033 5173 -0.2 0.842
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 16880 on 70 degrees of freedom
Multiple R-squared: 0.0005698, Adjusted R-squared: -0.01371
F-statistic: 0.03991 on 1 and 70 DF, p-value: 0.8422
> 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.0354627184 0.070925437 0.9645373
[2,] 0.0094682924 0.018936585 0.9905317
[3,] 0.0021296381 0.004259276 0.9978704
[4,] 0.0011737165 0.002347433 0.9988263
[5,] 0.0014330231 0.002866046 0.9985670
[6,] 0.0010942107 0.002188421 0.9989058
[7,] 0.0028218183 0.005643637 0.9971782
[8,] 0.0148938128 0.029787626 0.9851062
[9,] 0.0338157795 0.067631559 0.9661842
[10,] 0.0374604193 0.074920839 0.9625396
[11,] 0.0257332076 0.051466415 0.9742668
[12,] 0.0182505602 0.036501120 0.9817494
[13,] 0.0121117193 0.024223439 0.9878883
[14,] 0.0078335053 0.015667011 0.9921665
[15,] 0.0045388036 0.009077607 0.9954612
[16,] 0.0023388007 0.004677601 0.9976612
[17,] 0.0012052501 0.002410500 0.9987947
[18,] 0.0006439250 0.001287850 0.9993561
[19,] 0.0010850521 0.002170104 0.9989149
[20,] 0.0026272337 0.005254467 0.9973728
[21,] 0.0059766977 0.011953395 0.9940233
[22,] 0.0079809438 0.015961888 0.9920191
[23,] 0.0065991453 0.013198291 0.9934009
[24,] 0.0059344295 0.011868859 0.9940656
[25,] 0.0047506167 0.009501233 0.9952494
[26,] 0.0039342150 0.007868430 0.9960658
[27,] 0.0028478176 0.005695635 0.9971522
[28,] 0.0019894292 0.003978858 0.9980106
[29,] 0.0013787244 0.002757449 0.9986213
[30,] 0.0009201421 0.001840284 0.9990799
[31,] 0.0023875806 0.004775161 0.9976124
[32,] 0.0076271222 0.015254244 0.9923729
[33,] 0.0328435667 0.065687133 0.9671564
[34,] 0.0459838503 0.091967701 0.9540161
[35,] 0.0501041289 0.100208258 0.9498959
[36,] 0.0541302360 0.108260472 0.9458698
[37,] 0.0657443098 0.131488620 0.9342557
[38,] 0.0774209777 0.154841955 0.9225790
[39,] 0.0962355169 0.192471034 0.9037645
[40,] 0.1178789624 0.235757925 0.8821210
[41,] 0.1653036742 0.330607348 0.8346963
[42,] 0.3064802249 0.612960450 0.6935198
[43,] 0.3100907086 0.620181417 0.6899093
[44,] 0.3330386942 0.666077388 0.6669613
[45,] 0.3434460685 0.686892137 0.6565539
[46,] 0.3519864176 0.703972835 0.6480136
[47,] 0.3692739767 0.738547953 0.6307260
[48,] 0.3695241716 0.739048343 0.6304758
[49,] 0.3704495929 0.740899186 0.6295504
[50,] 0.3704983805 0.740996761 0.6295016
[51,] 0.3781334284 0.756266857 0.6218666
[52,] 0.3733655113 0.746731023 0.6266345
[53,] 0.4298681806 0.859736361 0.5701318
[54,] 0.4900948050 0.980189610 0.5099052
[55,] 0.4155514968 0.831102994 0.5844485
[56,] 0.4021088335 0.804217667 0.5978912
[57,] 0.4230148033 0.846029607 0.5769852
[58,] 0.5146644476 0.970671105 0.4853356
[59,] 0.6711632982 0.657673404 0.3288367
[60,] 0.7269999758 0.546000048 0.2730000
[61,] 0.6923335891 0.615332822 0.3076664
[62,] 0.6154759199 0.769048160 0.3845241
[63,] 0.4802551621 0.960510324 0.5197448
> postscript(file="/var/www/html/rcomp/tmp/16tkk1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index')
> points(x[,1]-mysum$resid)
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/rcomp/tmp/24yor1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/rcomp/tmp/353uc1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/rcomp/tmp/4p06l1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals')
> dev.off()
null device
1
> postscript(file="/var/www/html/rcomp/tmp/5fp1j1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> qqnorm(mysum$resid, main='Residual Normal Q-Q Plot')
> qqline(mysum$resid)
> grid()
> dev.off()
null device
1
> (myerror <- as.ts(mysum$resid))
Time Series:
Start = 1
End = 72
Frequency = 1
1 2 3 4 5 6
3217.9492 609.9492 -8314.0508 -3614.0508 985.9492 938.9492
7 8 9 10 11 12
-1650.0508 -7564.0508 -11583.0508 -11559.0508 10839.9492 18923.9492
13 14 15 16 17 18
20174.9492 16614.9492 10259.9492 11487.9492 10700.9492 10675.9492
19 20 21 22 23 24
8609.9492 3722.9492 712.9492 -1434.0508 20831.9492 24758.9492
25 26 27 28 29 30
25872.9492 21758.9492 15049.9492 16228.9492 13762.9492 13980.9492
31 32 33 34 35 36
10408.9492 8788.9492 7924.9492 4883.9492 23374.9492 24298.9492
37 38 39 40 41 42
25872.9492 12078.9492 3398.9492 -566.0508 2606.9492 -2251.0508
43 44 45 46 47 48
-11229.0508 -13074.0508 -20370.0508 -31055.0508 -7871.0508 -5434.0508
49 50 51 52 53 54
-11764.0508 -15784.0508 -23005.0508 -19322.0508 -17886.0508 -21388.0508
55 56 57 58 59 60
-29347.0508 -29342.0508 -40873.0508 -39588.0508 -18493.0508 -17060.6154
61 62 63 64 65 66
-18424.6154 -20295.6154 -19130.6154 -9628.6154 -1074.6154 2953.3846
67 68 69 70 71 72
7652.3846 8569.3846 2036.3846 6599.3846 27316.3846 30487.3846
> postscript(file="/var/www/html/rcomp/tmp/633m91258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> dum <- cbind(lag(myerror,k=1),myerror)
> dum
Time Series:
Start = 0
End = 72
Frequency = 1
lag(myerror, k = 1) myerror
0 3217.9492 NA
1 609.9492 3217.9492
2 -8314.0508 609.9492
3 -3614.0508 -8314.0508
4 985.9492 -3614.0508
5 938.9492 985.9492
6 -1650.0508 938.9492
7 -7564.0508 -1650.0508
8 -11583.0508 -7564.0508
9 -11559.0508 -11583.0508
10 10839.9492 -11559.0508
11 18923.9492 10839.9492
12 20174.9492 18923.9492
13 16614.9492 20174.9492
14 10259.9492 16614.9492
15 11487.9492 10259.9492
16 10700.9492 11487.9492
17 10675.9492 10700.9492
18 8609.9492 10675.9492
19 3722.9492 8609.9492
20 712.9492 3722.9492
21 -1434.0508 712.9492
22 20831.9492 -1434.0508
23 24758.9492 20831.9492
24 25872.9492 24758.9492
25 21758.9492 25872.9492
26 15049.9492 21758.9492
27 16228.9492 15049.9492
28 13762.9492 16228.9492
29 13980.9492 13762.9492
30 10408.9492 13980.9492
31 8788.9492 10408.9492
32 7924.9492 8788.9492
33 4883.9492 7924.9492
34 23374.9492 4883.9492
35 24298.9492 23374.9492
36 25872.9492 24298.9492
37 12078.9492 25872.9492
38 3398.9492 12078.9492
39 -566.0508 3398.9492
40 2606.9492 -566.0508
41 -2251.0508 2606.9492
42 -11229.0508 -2251.0508
43 -13074.0508 -11229.0508
44 -20370.0508 -13074.0508
45 -31055.0508 -20370.0508
46 -7871.0508 -31055.0508
47 -5434.0508 -7871.0508
48 -11764.0508 -5434.0508
49 -15784.0508 -11764.0508
50 -23005.0508 -15784.0508
51 -19322.0508 -23005.0508
52 -17886.0508 -19322.0508
53 -21388.0508 -17886.0508
54 -29347.0508 -21388.0508
55 -29342.0508 -29347.0508
56 -40873.0508 -29342.0508
57 -39588.0508 -40873.0508
58 -18493.0508 -39588.0508
59 -17060.6154 -18493.0508
60 -18424.6154 -17060.6154
61 -20295.6154 -18424.6154
62 -19130.6154 -20295.6154
63 -9628.6154 -19130.6154
64 -1074.6154 -9628.6154
65 2953.3846 -1074.6154
66 7652.3846 2953.3846
67 8569.3846 7652.3846
68 2036.3846 8569.3846
69 6599.3846 2036.3846
70 27316.3846 6599.3846
71 30487.3846 27316.3846
72 NA 30487.3846
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 609.9492 3217.9492
[2,] -8314.0508 609.9492
[3,] -3614.0508 -8314.0508
[4,] 985.9492 -3614.0508
[5,] 938.9492 985.9492
[6,] -1650.0508 938.9492
[7,] -7564.0508 -1650.0508
[8,] -11583.0508 -7564.0508
[9,] -11559.0508 -11583.0508
[10,] 10839.9492 -11559.0508
[11,] 18923.9492 10839.9492
[12,] 20174.9492 18923.9492
[13,] 16614.9492 20174.9492
[14,] 10259.9492 16614.9492
[15,] 11487.9492 10259.9492
[16,] 10700.9492 11487.9492
[17,] 10675.9492 10700.9492
[18,] 8609.9492 10675.9492
[19,] 3722.9492 8609.9492
[20,] 712.9492 3722.9492
[21,] -1434.0508 712.9492
[22,] 20831.9492 -1434.0508
[23,] 24758.9492 20831.9492
[24,] 25872.9492 24758.9492
[25,] 21758.9492 25872.9492
[26,] 15049.9492 21758.9492
[27,] 16228.9492 15049.9492
[28,] 13762.9492 16228.9492
[29,] 13980.9492 13762.9492
[30,] 10408.9492 13980.9492
[31,] 8788.9492 10408.9492
[32,] 7924.9492 8788.9492
[33,] 4883.9492 7924.9492
[34,] 23374.9492 4883.9492
[35,] 24298.9492 23374.9492
[36,] 25872.9492 24298.9492
[37,] 12078.9492 25872.9492
[38,] 3398.9492 12078.9492
[39,] -566.0508 3398.9492
[40,] 2606.9492 -566.0508
[41,] -2251.0508 2606.9492
[42,] -11229.0508 -2251.0508
[43,] -13074.0508 -11229.0508
[44,] -20370.0508 -13074.0508
[45,] -31055.0508 -20370.0508
[46,] -7871.0508 -31055.0508
[47,] -5434.0508 -7871.0508
[48,] -11764.0508 -5434.0508
[49,] -15784.0508 -11764.0508
[50,] -23005.0508 -15784.0508
[51,] -19322.0508 -23005.0508
[52,] -17886.0508 -19322.0508
[53,] -21388.0508 -17886.0508
[54,] -29347.0508 -21388.0508
[55,] -29342.0508 -29347.0508
[56,] -40873.0508 -29342.0508
[57,] -39588.0508 -40873.0508
[58,] -18493.0508 -39588.0508
[59,] -17060.6154 -18493.0508
[60,] -18424.6154 -17060.6154
[61,] -20295.6154 -18424.6154
[62,] -19130.6154 -20295.6154
[63,] -9628.6154 -19130.6154
[64,] -1074.6154 -9628.6154
[65,] 2953.3846 -1074.6154
[66,] 7652.3846 2953.3846
[67,] 8569.3846 7652.3846
[68,] 2036.3846 8569.3846
[69,] 6599.3846 2036.3846
[70,] 27316.3846 6599.3846
[71,] 30487.3846 27316.3846
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 609.9492 3217.9492
2 -8314.0508 609.9492
3 -3614.0508 -8314.0508
4 985.9492 -3614.0508
5 938.9492 985.9492
6 -1650.0508 938.9492
7 -7564.0508 -1650.0508
8 -11583.0508 -7564.0508
9 -11559.0508 -11583.0508
10 10839.9492 -11559.0508
11 18923.9492 10839.9492
12 20174.9492 18923.9492
13 16614.9492 20174.9492
14 10259.9492 16614.9492
15 11487.9492 10259.9492
16 10700.9492 11487.9492
17 10675.9492 10700.9492
18 8609.9492 10675.9492
19 3722.9492 8609.9492
20 712.9492 3722.9492
21 -1434.0508 712.9492
22 20831.9492 -1434.0508
23 24758.9492 20831.9492
24 25872.9492 24758.9492
25 21758.9492 25872.9492
26 15049.9492 21758.9492
27 16228.9492 15049.9492
28 13762.9492 16228.9492
29 13980.9492 13762.9492
30 10408.9492 13980.9492
31 8788.9492 10408.9492
32 7924.9492 8788.9492
33 4883.9492 7924.9492
34 23374.9492 4883.9492
35 24298.9492 23374.9492
36 25872.9492 24298.9492
37 12078.9492 25872.9492
38 3398.9492 12078.9492
39 -566.0508 3398.9492
40 2606.9492 -566.0508
41 -2251.0508 2606.9492
42 -11229.0508 -2251.0508
43 -13074.0508 -11229.0508
44 -20370.0508 -13074.0508
45 -31055.0508 -20370.0508
46 -7871.0508 -31055.0508
47 -5434.0508 -7871.0508
48 -11764.0508 -5434.0508
49 -15784.0508 -11764.0508
50 -23005.0508 -15784.0508
51 -19322.0508 -23005.0508
52 -17886.0508 -19322.0508
53 -21388.0508 -17886.0508
54 -29347.0508 -21388.0508
55 -29342.0508 -29347.0508
56 -40873.0508 -29342.0508
57 -39588.0508 -40873.0508
58 -18493.0508 -39588.0508
59 -17060.6154 -18493.0508
60 -18424.6154 -17060.6154
61 -20295.6154 -18424.6154
62 -19130.6154 -20295.6154
63 -9628.6154 -19130.6154
64 -1074.6154 -9628.6154
65 2953.3846 -1074.6154
66 7652.3846 2953.3846
67 8569.3846 7652.3846
68 2036.3846 8569.3846
69 6599.3846 2036.3846
70 27316.3846 6599.3846
71 30487.3846 27316.3846
> plot(z,main=paste('Residual Lag plot, lowess, and regression line'), ylab='values of Residuals', xlab='lagged values of Residuals')
> lines(lowess(z))
> abline(lm(z))
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/rcomp/tmp/79wdk1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/rcomp/tmp/81my11258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/rcomp/tmp/9kp8w1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
> plot(mylm, las = 1, sub='Residual Diagnostics')
> par(opar)
> dev.off()
null device
1
> if (n > n25) {
+ postscript(file="/var/www/html/rcomp/tmp/10d0cv1258739090.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
+ plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint')
+ grid()
+ dev.off()
+ }
null device
1
>
> #Note: the /var/www/html/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/html/rcomp/createtable")
>
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Estimated Regression Equation', 1, TRUE)
> a<-table.row.end(a)
> myeq <- colnames(x)[1]
> myeq <- paste(myeq, '[t] = ', sep='')
> for (i in 1:k){
+ if (mysum$coefficients[i,1] > 0) myeq <- paste(myeq, '+', '')
+ myeq <- paste(myeq, mysum$coefficients[i,1], sep=' ')
+ if (rownames(mysum$coefficients)[i] != '(Intercept)') {
+ myeq <- paste(myeq, rownames(mysum$coefficients)[i], sep='')
+ if (rownames(mysum$coefficients)[i] != 't') myeq <- paste(myeq, '[t]', sep='')
+ }
+ }
> myeq <- paste(myeq, ' + e[t]')
> a<-table.row.start(a)
> a<-table.element(a, myeq)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/www/html/rcomp/tmp/119w7z1258739090.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,hyperlink('http://www.xycoon.com/ols1.htm','Multiple Linear Regression - Ordinary Least Squares',''), 6, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Variable',header=TRUE)
> a<-table.element(a,'Parameter',header=TRUE)
> a<-table.element(a,'S.D.',header=TRUE)
> a<-table.element(a,'T-STAT
H0: parameter = 0',header=TRUE)
> a<-table.element(a,'2-tail p-value',header=TRUE)
> a<-table.element(a,'1-tail p-value',header=TRUE)
> a<-table.row.end(a)
> for (i in 1:k){
+ a<-table.row.start(a)
+ a<-table.element(a,rownames(mysum$coefficients)[i],header=TRUE)
+ a<-table.element(a,mysum$coefficients[i,1])
+ a<-table.element(a, round(mysum$coefficients[i,2],6))
+ a<-table.element(a, round(mysum$coefficients[i,3],4))
+ a<-table.element(a, round(mysum$coefficients[i,4],6))
+ a<-table.element(a, round(mysum$coefficients[i,4]/2,6))
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/www/html/rcomp/tmp/12ujzt1258739090.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Regression Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple R',1,TRUE)
> a<-table.element(a, sqrt(mysum$r.squared))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'R-squared',1,TRUE)
> a<-table.element(a, mysum$r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Adjusted R-squared',1,TRUE)
> a<-table.element(a, mysum$adj.r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (value)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[3])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'p-value',1,TRUE)
> a<-table.element(a, 1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3]))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Residual Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Residual Standard Deviation',1,TRUE)
> a<-table.element(a, mysum$sigma)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Sum Squared Residuals',1,TRUE)
> a<-table.element(a, sum(myerror*myerror))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/www/html/rcomp/tmp/13vwn51258739090.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Actuals, Interpolation, and Residuals', 4, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Time or Index', 1, TRUE)
> a<-table.element(a, 'Actuals', 1, TRUE)
> a<-table.element(a, 'Interpolation
Forecast', 1, TRUE)
> a<-table.element(a, 'Residuals
Prediction Error', 1, TRUE)
> a<-table.row.end(a)
> for (i in 1:n) {
+ a<-table.row.start(a)
+ a<-table.element(a,i, 1, TRUE)
+ a<-table.element(a,x[i])
+ a<-table.element(a,x[i]-mysum$resid[i])
+ a<-table.element(a,mysum$resid[i])
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/www/html/rcomp/tmp/14jnxp1258739090.tab")
> if (n > n25) {
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'p-values',header=TRUE)
+ a<-table.element(a,'Alternative Hypothesis',3,header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'breakpoint index',header=TRUE)
+ a<-table.element(a,'greater',header=TRUE)
+ a<-table.element(a,'2-sided',header=TRUE)
+ a<-table.element(a,'less',header=TRUE)
+ a<-table.row.end(a)
+ for (mypoint in kp3:nmkm3) {
+ a<-table.row.start(a)
+ a<-table.element(a,mypoint,header=TRUE)
+ a<-table.element(a,gqarr[mypoint-kp3+1,1])
+ a<-table.element(a,gqarr[mypoint-kp3+1,2])
+ a<-table.element(a,gqarr[mypoint-kp3+1,3])
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/www/html/rcomp/tmp/15mhwu1258739090.tab")
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Meta Analysis of Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'Description',header=TRUE)
+ a<-table.element(a,'# significant tests',header=TRUE)
+ a<-table.element(a,'% significant tests',header=TRUE)
+ a<-table.element(a,'OK/NOK',header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'1% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant1)
+ a<-table.element(a,numsignificant1/numgqtests)
+ if (numsignificant1/numgqtests < 0.01) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'5% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant5)
+ a<-table.element(a,numsignificant5/numgqtests)
+ if (numsignificant5/numgqtests < 0.05) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'10% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant10)
+ a<-table.element(a,numsignificant10/numgqtests)
+ if (numsignificant10/numgqtests < 0.1) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.end(a)
+ table.save(a,file="/var/www/html/rcomp/tmp/16nd671258739090.tab")
+ }
>
> system("convert tmp/16tkk1258739090.ps tmp/16tkk1258739090.png")
> system("convert tmp/24yor1258739090.ps tmp/24yor1258739090.png")
> system("convert tmp/353uc1258739090.ps tmp/353uc1258739090.png")
> system("convert tmp/4p06l1258739090.ps tmp/4p06l1258739090.png")
> system("convert tmp/5fp1j1258739090.ps tmp/5fp1j1258739090.png")
> system("convert tmp/633m91258739090.ps tmp/633m91258739090.png")
> system("convert tmp/79wdk1258739090.ps tmp/79wdk1258739090.png")
> system("convert tmp/81my11258739090.ps tmp/81my11258739090.png")
> system("convert tmp/9kp8w1258739090.ps tmp/9kp8w1258739090.png")
> system("convert tmp/10d0cv1258739090.ps tmp/10d0cv1258739090.png")
>
>
> proc.time()
user system elapsed
2.621 1.587 5.752