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(2350.44,0,2440.25,0,2408.64,0,2472.81,0,2407.6,0,2454.62,0,2448.05,0,2497.84,0,2645.64,0,2756.76,0,2849.27,0,2921.44,0,2981.85,0,3080.58,0,3106.22,0,3119.31,0,3061.26,0,3097.31,0,3161.69,0,3257.16,0,3277.01,0,3295.32,0,3363.99,0,3494.17,0,3667.03,0,3813.06,0,3917.96,0,3895.51,0,3801.06,0,3570.12,0,3701.61,0,3862.27,0,3970.1,0,4138.52,0,4199.75,0,4290.89,0,4443.91,0,4502.64,0,4356.98,0,4591.27,0,4696.96,0,4621.4,0,4562.84,0,4202.52,0,4296.49,0,4435.23,0,4105.18,0,4116.68,0,3844.49,0,3720.98,0,3674.4,0,3857.62,0,3801.06,0,3504.37,0,3032.6,0,3047.03,0,2962.34,1,2197.82,1,2014.45,1,1862.83,1,1905.41,1,1810.99,1,1670.07,1,1864.44,1,2052.02,1,2029.6,1,2070.83,1,2293.41,1,2443.27,1,2513.17,1,2466.92,1),dim=c(2,71),dimnames=list(c('Y','X'),1:71))
> y <- array(NA,dim=c(2,71),dimnames=list(c('Y','X'),1:71))
> 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 2350.44 0
2 2440.25 0
3 2408.64 0
4 2472.81 0
5 2407.60 0
6 2454.62 0
7 2448.05 0
8 2497.84 0
9 2645.64 0
10 2756.76 0
11 2849.27 0
12 2921.44 0
13 2981.85 0
14 3080.58 0
15 3106.22 0
16 3119.31 0
17 3061.26 0
18 3097.31 0
19 3161.69 0
20 3257.16 0
21 3277.01 0
22 3295.32 0
23 3363.99 0
24 3494.17 0
25 3667.03 0
26 3813.06 0
27 3917.96 0
28 3895.51 0
29 3801.06 0
30 3570.12 0
31 3701.61 0
32 3862.27 0
33 3970.10 0
34 4138.52 0
35 4199.75 0
36 4290.89 0
37 4443.91 0
38 4502.64 0
39 4356.98 0
40 4591.27 0
41 4696.96 0
42 4621.40 0
43 4562.84 0
44 4202.52 0
45 4296.49 0
46 4435.23 0
47 4105.18 0
48 4116.68 0
49 3844.49 0
50 3720.98 0
51 3674.40 0
52 3857.62 0
53 3801.06 0
54 3504.37 0
55 3032.60 0
56 3047.03 0
57 2962.34 1
58 2197.82 1
59 2014.45 1
60 1862.83 1
61 1905.41 1
62 1810.99 1
63 1670.07 1
64 1864.44 1
65 2052.02 1
66 2029.60 1
67 2070.83 1
68 2293.41 1
69 2443.27 1
70 2513.17 1
71 2466.92 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X
3521 -1377
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1170.84 -432.34 -16.91 385.45 1175.68
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3521.28 84.83 41.512 < 2e-16 ***
X -1377.44 184.55 -7.464 1.91e-10 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 634.8 on 69 degrees of freedom
Multiple R-squared: 0.4467, Adjusted R-squared: 0.4387
F-statistic: 55.71 on 1 and 69 DF, p-value: 1.907e-10
> 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,] 1.055483e-03 2.110965e-03 0.9989445173
[2,] 1.176344e-04 2.352688e-04 0.9998823656
[3,] 1.208118e-05 2.416237e-05 0.9999879188
[4,] 3.327283e-06 6.654567e-06 0.9999966727
[5,] 3.237881e-05 6.475762e-05 0.9999676212
[6,] 2.097001e-04 4.194002e-04 0.9997902999
[7,] 8.602189e-04 1.720438e-03 0.9991397811
[8,] 2.494453e-03 4.988906e-03 0.9975055472
[9,] 5.794197e-03 1.158839e-02 0.9942058034
[10,] 1.412982e-02 2.825964e-02 0.9858701824
[11,] 2.582169e-02 5.164337e-02 0.9741783140
[12,] 3.950561e-02 7.901121e-02 0.9604943932
[13,] 4.864649e-02 9.729298e-02 0.9513535075
[14,] 6.195352e-02 1.239070e-01 0.9380464774
[15,] 8.396589e-02 1.679318e-01 0.9160341087
[16,] 1.225830e-01 2.451660e-01 0.8774170094
[17,] 1.686878e-01 3.373756e-01 0.8313121758
[18,] 2.223418e-01 4.446836e-01 0.7776581961
[19,] 2.925671e-01 5.851341e-01 0.7074329390
[20,] 3.919415e-01 7.838830e-01 0.6080585130
[21,] 5.280651e-01 9.438698e-01 0.4719348790
[22,] 6.715542e-01 6.568916e-01 0.3284458250
[23,] 7.868676e-01 4.262649e-01 0.2131324376
[24,] 8.473281e-01 3.053437e-01 0.1526718699
[25,] 8.723483e-01 2.553035e-01 0.1276517411
[26,] 8.797799e-01 2.404401e-01 0.1202200668
[27,] 8.881638e-01 2.236724e-01 0.1118361944
[28,] 9.003494e-01 1.993011e-01 0.0996505585
[29,] 9.135834e-01 1.728332e-01 0.0864166093
[30,] 9.323603e-01 1.352794e-01 0.0676397026
[31,] 9.467490e-01 1.065019e-01 0.0532509551
[32,] 9.599105e-01 8.017908e-02 0.0400895387
[33,] 9.747879e-01 5.042411e-02 0.0252120575
[34,] 9.849130e-01 3.017403e-02 0.0150870148
[35,] 9.872880e-01 2.542393e-02 0.0127119643
[36,] 9.932844e-01 1.343117e-02 0.0067155843
[37,] 9.975454e-01 4.909161e-03 0.0024545804
[38,] 9.990009e-01 1.998108e-03 0.0009990542
[39,] 9.995838e-01 8.323943e-04 0.0004161972
[40,] 9.995246e-01 9.507356e-04 0.0004753678
[41,] 9.996051e-01 7.898517e-04 0.0003949258
[42,] 9.998449e-01 3.101391e-04 0.0001550695
[43,] 9.998349e-01 3.301246e-04 0.0001650623
[44,] 9.998588e-01 2.824773e-04 0.0001412387
[45,] 9.997687e-01 4.626919e-04 0.0002313459
[46,] 9.995587e-01 8.826496e-04 0.0004413248
[47,] 9.991566e-01 1.686727e-03 0.0008433633
[48,] 9.990530e-01 1.894041e-03 0.0009470207
[49,] 9.991634e-01 1.673181e-03 0.0008365905
[50,] 9.988586e-01 2.282892e-03 0.0011414458
[51,] 9.975886e-01 4.822724e-03 0.0024113618
[52,] 9.950376e-01 9.924880e-03 0.0049624400
[53,] 9.993664e-01 1.267233e-03 0.0006336164
[54,] 9.985294e-01 2.941283e-03 0.0014706416
[55,] 9.964737e-01 7.052699e-03 0.0035263497
[56,] 9.933886e-01 1.322289e-02 0.0066114436
[57,] 9.869795e-01 2.604091e-02 0.0130204542
[58,] 9.805706e-01 3.885875e-02 0.0194293764
[59,] 9.874113e-01 2.517750e-02 0.0125887488
[60,] 9.864959e-01 2.700823e-02 0.0135041148
[61,] 9.723433e-01 5.531331e-02 0.0276566573
[62,] 9.589526e-01 8.209487e-02 0.0410474338
> postscript(file="/var/www/html/rcomp/tmp/1pt4j1260884326.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/2ea901260884326.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/3s89c1260884326.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/45irx1260884326.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/502kj1260884326.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 = 71
Frequency = 1
1 2 3 4 5 6
-1170.84143 -1081.03143 -1112.64143 -1048.47143 -1113.68143 -1066.66143
7 8 9 10 11 12
-1073.23143 -1023.44143 -875.64143 -764.52143 -672.01143 -599.84143
13 14 15 16 17 18
-539.43143 -440.70143 -415.06143 -401.97143 -460.02143 -423.97143
19 20 21 22 23 24
-359.59143 -264.12143 -244.27143 -225.96143 -157.29143 -27.11143
25 26 27 28 29 30
145.74857 291.77857 396.67857 374.22857 279.77857 48.83857
31 32 33 34 35 36
180.32857 340.98857 448.81857 617.23857 678.46857 769.60857
37 38 39 40 41 42
922.62857 981.35857 835.69857 1069.98857 1175.67857 1100.11857
43 44 45 46 47 48
1041.55857 681.23857 775.20857 913.94857 583.89857 595.39857
49 50 51 52 53 54
323.20857 199.69857 153.11857 336.33857 279.77857 -16.91143
55 56 57 58 59 60
-488.68143 -474.25143 818.50200 53.98200 -129.38800 -281.00800
61 62 63 64 65 66
-238.42800 -332.84800 -473.76800 -279.39800 -91.81800 -114.23800
67 68 69 70 71
-73.00800 149.57200 299.43200 369.33200 323.08200
> postscript(file="/var/www/html/rcomp/tmp/6pwno1260884326.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 = 71
Frequency = 1
lag(myerror, k = 1) myerror
0 -1170.84143 NA
1 -1081.03143 -1170.84143
2 -1112.64143 -1081.03143
3 -1048.47143 -1112.64143
4 -1113.68143 -1048.47143
5 -1066.66143 -1113.68143
6 -1073.23143 -1066.66143
7 -1023.44143 -1073.23143
8 -875.64143 -1023.44143
9 -764.52143 -875.64143
10 -672.01143 -764.52143
11 -599.84143 -672.01143
12 -539.43143 -599.84143
13 -440.70143 -539.43143
14 -415.06143 -440.70143
15 -401.97143 -415.06143
16 -460.02143 -401.97143
17 -423.97143 -460.02143
18 -359.59143 -423.97143
19 -264.12143 -359.59143
20 -244.27143 -264.12143
21 -225.96143 -244.27143
22 -157.29143 -225.96143
23 -27.11143 -157.29143
24 145.74857 -27.11143
25 291.77857 145.74857
26 396.67857 291.77857
27 374.22857 396.67857
28 279.77857 374.22857
29 48.83857 279.77857
30 180.32857 48.83857
31 340.98857 180.32857
32 448.81857 340.98857
33 617.23857 448.81857
34 678.46857 617.23857
35 769.60857 678.46857
36 922.62857 769.60857
37 981.35857 922.62857
38 835.69857 981.35857
39 1069.98857 835.69857
40 1175.67857 1069.98857
41 1100.11857 1175.67857
42 1041.55857 1100.11857
43 681.23857 1041.55857
44 775.20857 681.23857
45 913.94857 775.20857
46 583.89857 913.94857
47 595.39857 583.89857
48 323.20857 595.39857
49 199.69857 323.20857
50 153.11857 199.69857
51 336.33857 153.11857
52 279.77857 336.33857
53 -16.91143 279.77857
54 -488.68143 -16.91143
55 -474.25143 -488.68143
56 818.50200 -474.25143
57 53.98200 818.50200
58 -129.38800 53.98200
59 -281.00800 -129.38800
60 -238.42800 -281.00800
61 -332.84800 -238.42800
62 -473.76800 -332.84800
63 -279.39800 -473.76800
64 -91.81800 -279.39800
65 -114.23800 -91.81800
66 -73.00800 -114.23800
67 149.57200 -73.00800
68 299.43200 149.57200
69 369.33200 299.43200
70 323.08200 369.33200
71 NA 323.08200
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -1081.03143 -1170.84143
[2,] -1112.64143 -1081.03143
[3,] -1048.47143 -1112.64143
[4,] -1113.68143 -1048.47143
[5,] -1066.66143 -1113.68143
[6,] -1073.23143 -1066.66143
[7,] -1023.44143 -1073.23143
[8,] -875.64143 -1023.44143
[9,] -764.52143 -875.64143
[10,] -672.01143 -764.52143
[11,] -599.84143 -672.01143
[12,] -539.43143 -599.84143
[13,] -440.70143 -539.43143
[14,] -415.06143 -440.70143
[15,] -401.97143 -415.06143
[16,] -460.02143 -401.97143
[17,] -423.97143 -460.02143
[18,] -359.59143 -423.97143
[19,] -264.12143 -359.59143
[20,] -244.27143 -264.12143
[21,] -225.96143 -244.27143
[22,] -157.29143 -225.96143
[23,] -27.11143 -157.29143
[24,] 145.74857 -27.11143
[25,] 291.77857 145.74857
[26,] 396.67857 291.77857
[27,] 374.22857 396.67857
[28,] 279.77857 374.22857
[29,] 48.83857 279.77857
[30,] 180.32857 48.83857
[31,] 340.98857 180.32857
[32,] 448.81857 340.98857
[33,] 617.23857 448.81857
[34,] 678.46857 617.23857
[35,] 769.60857 678.46857
[36,] 922.62857 769.60857
[37,] 981.35857 922.62857
[38,] 835.69857 981.35857
[39,] 1069.98857 835.69857
[40,] 1175.67857 1069.98857
[41,] 1100.11857 1175.67857
[42,] 1041.55857 1100.11857
[43,] 681.23857 1041.55857
[44,] 775.20857 681.23857
[45,] 913.94857 775.20857
[46,] 583.89857 913.94857
[47,] 595.39857 583.89857
[48,] 323.20857 595.39857
[49,] 199.69857 323.20857
[50,] 153.11857 199.69857
[51,] 336.33857 153.11857
[52,] 279.77857 336.33857
[53,] -16.91143 279.77857
[54,] -488.68143 -16.91143
[55,] -474.25143 -488.68143
[56,] 818.50200 -474.25143
[57,] 53.98200 818.50200
[58,] -129.38800 53.98200
[59,] -281.00800 -129.38800
[60,] -238.42800 -281.00800
[61,] -332.84800 -238.42800
[62,] -473.76800 -332.84800
[63,] -279.39800 -473.76800
[64,] -91.81800 -279.39800
[65,] -114.23800 -91.81800
[66,] -73.00800 -114.23800
[67,] 149.57200 -73.00800
[68,] 299.43200 149.57200
[69,] 369.33200 299.43200
[70,] 323.08200 369.33200
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -1081.03143 -1170.84143
2 -1112.64143 -1081.03143
3 -1048.47143 -1112.64143
4 -1113.68143 -1048.47143
5 -1066.66143 -1113.68143
6 -1073.23143 -1066.66143
7 -1023.44143 -1073.23143
8 -875.64143 -1023.44143
9 -764.52143 -875.64143
10 -672.01143 -764.52143
11 -599.84143 -672.01143
12 -539.43143 -599.84143
13 -440.70143 -539.43143
14 -415.06143 -440.70143
15 -401.97143 -415.06143
16 -460.02143 -401.97143
17 -423.97143 -460.02143
18 -359.59143 -423.97143
19 -264.12143 -359.59143
20 -244.27143 -264.12143
21 -225.96143 -244.27143
22 -157.29143 -225.96143
23 -27.11143 -157.29143
24 145.74857 -27.11143
25 291.77857 145.74857
26 396.67857 291.77857
27 374.22857 396.67857
28 279.77857 374.22857
29 48.83857 279.77857
30 180.32857 48.83857
31 340.98857 180.32857
32 448.81857 340.98857
33 617.23857 448.81857
34 678.46857 617.23857
35 769.60857 678.46857
36 922.62857 769.60857
37 981.35857 922.62857
38 835.69857 981.35857
39 1069.98857 835.69857
40 1175.67857 1069.98857
41 1100.11857 1175.67857
42 1041.55857 1100.11857
43 681.23857 1041.55857
44 775.20857 681.23857
45 913.94857 775.20857
46 583.89857 913.94857
47 595.39857 583.89857
48 323.20857 595.39857
49 199.69857 323.20857
50 153.11857 199.69857
51 336.33857 153.11857
52 279.77857 336.33857
53 -16.91143 279.77857
54 -488.68143 -16.91143
55 -474.25143 -488.68143
56 818.50200 -474.25143
57 53.98200 818.50200
58 -129.38800 53.98200
59 -281.00800 -129.38800
60 -238.42800 -281.00800
61 -332.84800 -238.42800
62 -473.76800 -332.84800
63 -279.39800 -473.76800
64 -91.81800 -279.39800
65 -114.23800 -91.81800
66 -73.00800 -114.23800
67 149.57200 -73.00800
68 299.43200 149.57200
69 369.33200 299.43200
70 323.08200 369.33200
> 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/7jn5k1260884326.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/8x0r81260884326.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/99p1n1260884326.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/103a9l1260884326.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/11kpuf1260884326.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/12kew81260884326.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/13z5ay1260884326.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/142ktm1260884326.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/1569lh1260884326.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/166mvi1260884326.tab")
+ }
>
> try(system("convert tmp/1pt4j1260884326.ps tmp/1pt4j1260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/2ea901260884326.ps tmp/2ea901260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/3s89c1260884326.ps tmp/3s89c1260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/45irx1260884326.ps tmp/45irx1260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/502kj1260884326.ps tmp/502kj1260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/6pwno1260884326.ps tmp/6pwno1260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/7jn5k1260884326.ps tmp/7jn5k1260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/8x0r81260884326.ps tmp/8x0r81260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/99p1n1260884326.ps tmp/99p1n1260884326.png",intern=TRUE))
character(0)
> try(system("convert tmp/103a9l1260884326.ps tmp/103a9l1260884326.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.593 1.611 6.894