R version 2.7.0 (2008-04-22)
Copyright (C) 2008 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(13698.3,0,12477.6,0,13139.7,0,14532.2,0,15167,0,16071.1,0,14827.5,0,15082,0,14772.7,0,16083,0,14272.5,0,15223.3,0,14897.3,0,13062.6,0,12603.8,0,13629.8,0,14421.1,0,13978.3,0,12927.9,0,13429.9,0,13470.1,0,14785.8,0,14292,0,14308.8,0,14013,0,13240.9,0,12153.4,0,14289.7,0,15669.2,0,14169.5,0,14569.8,0,14469.1,0,14264.9,0,15320.9,0,14433.5,0,13691.5,0,14194.1,0,13519.2,0,11857.9,0,14616,0,15643.4,0,14077.2,0,14887.5,0,14159.9,0,14643,0,17192.5,1,15386.1,1,14287.1,1,17526.6,1,14497,1,14398.3,1,16629.6,1,16670.7,1,16614.8,1,16869.2,1,15663.9,1,16359.9,1,18447.7,1,16889,1,16505,1,18320.9,1,15052.1,1,15699.8,1,18135.3,1,16768.7,1,18883,1,19021,1,18101.9,1,17776.1,1,21489.9,1,17065.3,1,18690,1,18953.1,1,16398.9,1,16895.7,1,18553,1,19270,1,19422.1,1,17579.4,1,18637.3,1,18076.7,1,20438.6,1,18075.2,1,19563,1,19899.2,1,19227.5,1,17789.6,1,19220.8,1,22058.6,1,21230.8,1,19504.4,1,23913.1,1,23165.7,1,23574.3,1,25002,1,22603.9,1,23408.6,1),dim=c(2,97),dimnames=list(c('y','x'),1:97))
> y <- array(NA,dim=c(2,97),dimnames=list(c('y','x'),1:97))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = 'Linear Trend'
> par2 = 'Include Monthly 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)
> 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 M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 t
1 13698.3 0 1 0 0 0 0 0 0 0 0 0 0 1
2 12477.6 0 0 1 0 0 0 0 0 0 0 0 0 2
3 13139.7 0 0 0 1 0 0 0 0 0 0 0 0 3
4 14532.2 0 0 0 0 1 0 0 0 0 0 0 0 4
5 15167.0 0 0 0 0 0 1 0 0 0 0 0 0 5
6 16071.1 0 0 0 0 0 0 1 0 0 0 0 0 6
7 14827.5 0 0 0 0 0 0 0 1 0 0 0 0 7
8 15082.0 0 0 0 0 0 0 0 0 1 0 0 0 8
9 14772.7 0 0 0 0 0 0 0 0 0 1 0 0 9
10 16083.0 0 0 0 0 0 0 0 0 0 0 1 0 10
11 14272.5 0 0 0 0 0 0 0 0 0 0 0 1 11
12 15223.3 0 0 0 0 0 0 0 0 0 0 0 0 12
13 14897.3 0 1 0 0 0 0 0 0 0 0 0 0 13
14 13062.6 0 0 1 0 0 0 0 0 0 0 0 0 14
15 12603.8 0 0 0 1 0 0 0 0 0 0 0 0 15
16 13629.8 0 0 0 0 1 0 0 0 0 0 0 0 16
17 14421.1 0 0 0 0 0 1 0 0 0 0 0 0 17
18 13978.3 0 0 0 0 0 0 1 0 0 0 0 0 18
19 12927.9 0 0 0 0 0 0 0 1 0 0 0 0 19
20 13429.9 0 0 0 0 0 0 0 0 1 0 0 0 20
21 13470.1 0 0 0 0 0 0 0 0 0 1 0 0 21
22 14785.8 0 0 0 0 0 0 0 0 0 0 1 0 22
23 14292.0 0 0 0 0 0 0 0 0 0 0 0 1 23
24 14308.8 0 0 0 0 0 0 0 0 0 0 0 0 24
25 14013.0 0 1 0 0 0 0 0 0 0 0 0 0 25
26 13240.9 0 0 1 0 0 0 0 0 0 0 0 0 26
27 12153.4 0 0 0 1 0 0 0 0 0 0 0 0 27
28 14289.7 0 0 0 0 1 0 0 0 0 0 0 0 28
29 15669.2 0 0 0 0 0 1 0 0 0 0 0 0 29
30 14169.5 0 0 0 0 0 0 1 0 0 0 0 0 30
31 14569.8 0 0 0 0 0 0 0 1 0 0 0 0 31
32 14469.1 0 0 0 0 0 0 0 0 1 0 0 0 32
33 14264.9 0 0 0 0 0 0 0 0 0 1 0 0 33
34 15320.9 0 0 0 0 0 0 0 0 0 0 1 0 34
35 14433.5 0 0 0 0 0 0 0 0 0 0 0 1 35
36 13691.5 0 0 0 0 0 0 0 0 0 0 0 0 36
37 14194.1 0 1 0 0 0 0 0 0 0 0 0 0 37
38 13519.2 0 0 1 0 0 0 0 0 0 0 0 0 38
39 11857.9 0 0 0 1 0 0 0 0 0 0 0 0 39
40 14616.0 0 0 0 0 1 0 0 0 0 0 0 0 40
41 15643.4 0 0 0 0 0 1 0 0 0 0 0 0 41
42 14077.2 0 0 0 0 0 0 1 0 0 0 0 0 42
43 14887.5 0 0 0 0 0 0 0 1 0 0 0 0 43
44 14159.9 0 0 0 0 0 0 0 0 1 0 0 0 44
45 14643.0 0 0 0 0 0 0 0 0 0 1 0 0 45
46 17192.5 1 0 0 0 0 0 0 0 0 0 1 0 46
47 15386.1 1 0 0 0 0 0 0 0 0 0 0 1 47
48 14287.1 1 0 0 0 0 0 0 0 0 0 0 0 48
49 17526.6 1 1 0 0 0 0 0 0 0 0 0 0 49
50 14497.0 1 0 1 0 0 0 0 0 0 0 0 0 50
51 14398.3 1 0 0 1 0 0 0 0 0 0 0 0 51
52 16629.6 1 0 0 0 1 0 0 0 0 0 0 0 52
53 16670.7 1 0 0 0 0 1 0 0 0 0 0 0 53
54 16614.8 1 0 0 0 0 0 1 0 0 0 0 0 54
55 16869.2 1 0 0 0 0 0 0 1 0 0 0 0 55
56 15663.9 1 0 0 0 0 0 0 0 1 0 0 0 56
57 16359.9 1 0 0 0 0 0 0 0 0 1 0 0 57
58 18447.7 1 0 0 0 0 0 0 0 0 0 1 0 58
59 16889.0 1 0 0 0 0 0 0 0 0 0 0 1 59
60 16505.0 1 0 0 0 0 0 0 0 0 0 0 0 60
61 18320.9 1 1 0 0 0 0 0 0 0 0 0 0 61
62 15052.1 1 0 1 0 0 0 0 0 0 0 0 0 62
63 15699.8 1 0 0 1 0 0 0 0 0 0 0 0 63
64 18135.3 1 0 0 0 1 0 0 0 0 0 0 0 64
65 16768.7 1 0 0 0 0 1 0 0 0 0 0 0 65
66 18883.0 1 0 0 0 0 0 1 0 0 0 0 0 66
67 19021.0 1 0 0 0 0 0 0 1 0 0 0 0 67
68 18101.9 1 0 0 0 0 0 0 0 1 0 0 0 68
69 17776.1 1 0 0 0 0 0 0 0 0 1 0 0 69
70 21489.9 1 0 0 0 0 0 0 0 0 0 1 0 70
71 17065.3 1 0 0 0 0 0 0 0 0 0 0 1 71
72 18690.0 1 0 0 0 0 0 0 0 0 0 0 0 72
73 18953.1 1 1 0 0 0 0 0 0 0 0 0 0 73
74 16398.9 1 0 1 0 0 0 0 0 0 0 0 0 74
75 16895.7 1 0 0 1 0 0 0 0 0 0 0 0 75
76 18553.0 1 0 0 0 1 0 0 0 0 0 0 0 76
77 19270.0 1 0 0 0 0 1 0 0 0 0 0 0 77
78 19422.1 1 0 0 0 0 0 1 0 0 0 0 0 78
79 17579.4 1 0 0 0 0 0 0 1 0 0 0 0 79
80 18637.3 1 0 0 0 0 0 0 0 1 0 0 0 80
81 18076.7 1 0 0 0 0 0 0 0 0 1 0 0 81
82 20438.6 1 0 0 0 0 0 0 0 0 0 1 0 82
83 18075.2 1 0 0 0 0 0 0 0 0 0 0 1 83
84 19563.0 1 0 0 0 0 0 0 0 0 0 0 0 84
85 19899.2 1 1 0 0 0 0 0 0 0 0 0 0 85
86 19227.5 1 0 1 0 0 0 0 0 0 0 0 0 86
87 17789.6 1 0 0 1 0 0 0 0 0 0 0 0 87
88 19220.8 1 0 0 0 1 0 0 0 0 0 0 0 88
89 22058.6 1 0 0 0 0 1 0 0 0 0 0 0 89
90 21230.8 1 0 0 0 0 0 1 0 0 0 0 0 90
91 19504.4 1 0 0 0 0 0 0 1 0 0 0 0 91
92 23913.1 1 0 0 0 0 0 0 0 1 0 0 0 92
93 23165.7 1 0 0 0 0 0 0 0 0 1 0 0 93
94 23574.3 1 0 0 0 0 0 0 0 0 0 1 0 94
95 25002.0 1 0 0 0 0 0 0 0 0 0 0 1 95
96 22603.9 1 0 0 0 0 0 0 0 0 0 0 0 96
97 23408.6 1 1 0 0 0 0 0 0 0 0 0 0 97
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) x M1 M2 M3 M4
12329.80 310.79 776.24 -1332.97 -1780.45 22.80
M5 M6 M7 M8 M9 M10
700.31 467.29 -145.50 183.02 -13.25 1718.07
M11 t
148.15 80.28
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-2206.8 -849.3 -263.8 747.6 4586.8
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 12329.80 582.34 21.173 < 2e-16 ***
x 310.79 574.68 0.541 0.5901
M1 776.24 684.94 1.133 0.2604
M2 -1332.97 706.36 -1.887 0.0626 .
M3 -1780.45 705.85 -2.522 0.0136 *
M4 22.80 705.49 0.032 0.9743
M5 700.31 705.28 0.993 0.3236
M6 467.29 705.21 0.663 0.5094
M7 -145.50 705.30 -0.206 0.8371
M8 183.02 705.53 0.259 0.7960
M9 -13.25 705.91 -0.019 0.9851
M10 1718.07 704.59 2.438 0.0169 *
M11 148.15 704.36 0.210 0.8339
t 80.28 10.24 7.843 1.34e-11 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1409 on 83 degrees of freedom
Multiple R-squared: 0.7977, Adjusted R-squared: 0.766
F-statistic: 25.17 on 13 and 83 DF, p-value: < 2.2e-16
> postscript(file="/var/www/html/rcomp/tmp/1t2wf1227454851.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/2geik1227454851.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/3y2ir1227454851.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/43sw71227454851.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/59aha1227454851.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> qqnorm(mysum$resid, main='Residual Normal Q-Q Plot')
> grid()
> dev.off()
null device
1
> (myerror <- as.ts(mysum$resid))
Time Series:
Start = 1
End = 97
Frequency = 1
1 2 3 4 5 6
511.975368 1320.208809 2349.508809 1858.483809 1735.496309 2792.333809
7 8 9 10 11 12
2081.246309 1926.946309 1733.646309 1232.344561 911.482061 1930.157061
13 14 15 16 17 18
747.635140 941.868580 850.268580 -7.256420 26.256080 -263.806420
19 20 21 22 23 24
-781.693920 -688.493920 -532.293920 -1028.195668 -32.358168 52.316832
25 26 27 28 29 30
-1100.005089 156.828352 -563.471648 -310.696648 311.015852 -1035.946648
31 32 33 34 35 36
-103.134148 -612.634148 -700.834148 -1456.435896 -854.198396 -1528.323396
37 38 39 40 41 42
-1882.245317 -528.211877 -1822.311877 -947.736877 -678.124377 -2091.586877
43 44 45 46 47 48
-748.774377 -1885.174377 -1286.074377 -858.962142 -1175.724642 -2206.849642
49 50 51 52 53 54
176.128437 -824.538123 -556.038123 -208.263123 -924.950623 -828.113123
55 56 57 58 59 60
-41.200623 -1655.300623 -843.300623 -567.102371 -636.164871 -952.289871
61 62 63 64 65 66
7.088208 -1232.778352 -217.878352 334.096648 -1790.290852 476.746648
67 68 69 70 71 72
1147.259148 -180.640852 -390.440852 1511.757401 -1423.205099 269.369901
73 74 75 76 77 78
-324.052020 -849.318580 14.681420 -211.543580 -252.331080 52.506420
79 80 81 82 83 84
-1257.681080 -608.581080 -1053.181080 -502.882828 -1376.645328 179.029672
85 86 87 88 89 90
-341.292249 1015.941191 -54.758809 -507.083809 1572.928691 897.866191
91 92 93 94 95 96
-296.021309 3703.878691 3072.478691 1669.476943 4586.814443 2256.589443
97
2204.767522
> postscript(file="/var/www/html/rcomp/tmp/687qw1227454851.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 = 97
Frequency = 1
lag(myerror, k = 1) myerror
0 511.975368 NA
1 1320.208809 511.975368
2 2349.508809 1320.208809
3 1858.483809 2349.508809
4 1735.496309 1858.483809
5 2792.333809 1735.496309
6 2081.246309 2792.333809
7 1926.946309 2081.246309
8 1733.646309 1926.946309
9 1232.344561 1733.646309
10 911.482061 1232.344561
11 1930.157061 911.482061
12 747.635140 1930.157061
13 941.868580 747.635140
14 850.268580 941.868580
15 -7.256420 850.268580
16 26.256080 -7.256420
17 -263.806420 26.256080
18 -781.693920 -263.806420
19 -688.493920 -781.693920
20 -532.293920 -688.493920
21 -1028.195668 -532.293920
22 -32.358168 -1028.195668
23 52.316832 -32.358168
24 -1100.005089 52.316832
25 156.828352 -1100.005089
26 -563.471648 156.828352
27 -310.696648 -563.471648
28 311.015852 -310.696648
29 -1035.946648 311.015852
30 -103.134148 -1035.946648
31 -612.634148 -103.134148
32 -700.834148 -612.634148
33 -1456.435896 -700.834148
34 -854.198396 -1456.435896
35 -1528.323396 -854.198396
36 -1882.245317 -1528.323396
37 -528.211877 -1882.245317
38 -1822.311877 -528.211877
39 -947.736877 -1822.311877
40 -678.124377 -947.736877
41 -2091.586877 -678.124377
42 -748.774377 -2091.586877
43 -1885.174377 -748.774377
44 -1286.074377 -1885.174377
45 -858.962142 -1286.074377
46 -1175.724642 -858.962142
47 -2206.849642 -1175.724642
48 176.128437 -2206.849642
49 -824.538123 176.128437
50 -556.038123 -824.538123
51 -208.263123 -556.038123
52 -924.950623 -208.263123
53 -828.113123 -924.950623
54 -41.200623 -828.113123
55 -1655.300623 -41.200623
56 -843.300623 -1655.300623
57 -567.102371 -843.300623
58 -636.164871 -567.102371
59 -952.289871 -636.164871
60 7.088208 -952.289871
61 -1232.778352 7.088208
62 -217.878352 -1232.778352
63 334.096648 -217.878352
64 -1790.290852 334.096648
65 476.746648 -1790.290852
66 1147.259148 476.746648
67 -180.640852 1147.259148
68 -390.440852 -180.640852
69 1511.757401 -390.440852
70 -1423.205099 1511.757401
71 269.369901 -1423.205099
72 -324.052020 269.369901
73 -849.318580 -324.052020
74 14.681420 -849.318580
75 -211.543580 14.681420
76 -252.331080 -211.543580
77 52.506420 -252.331080
78 -1257.681080 52.506420
79 -608.581080 -1257.681080
80 -1053.181080 -608.581080
81 -502.882828 -1053.181080
82 -1376.645328 -502.882828
83 179.029672 -1376.645328
84 -341.292249 179.029672
85 1015.941191 -341.292249
86 -54.758809 1015.941191
87 -507.083809 -54.758809
88 1572.928691 -507.083809
89 897.866191 1572.928691
90 -296.021309 897.866191
91 3703.878691 -296.021309
92 3072.478691 3703.878691
93 1669.476943 3072.478691
94 4586.814443 1669.476943
95 2256.589443 4586.814443
96 2204.767522 2256.589443
97 NA 2204.767522
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 1320.208809 511.975368
[2,] 2349.508809 1320.208809
[3,] 1858.483809 2349.508809
[4,] 1735.496309 1858.483809
[5,] 2792.333809 1735.496309
[6,] 2081.246309 2792.333809
[7,] 1926.946309 2081.246309
[8,] 1733.646309 1926.946309
[9,] 1232.344561 1733.646309
[10,] 911.482061 1232.344561
[11,] 1930.157061 911.482061
[12,] 747.635140 1930.157061
[13,] 941.868580 747.635140
[14,] 850.268580 941.868580
[15,] -7.256420 850.268580
[16,] 26.256080 -7.256420
[17,] -263.806420 26.256080
[18,] -781.693920 -263.806420
[19,] -688.493920 -781.693920
[20,] -532.293920 -688.493920
[21,] -1028.195668 -532.293920
[22,] -32.358168 -1028.195668
[23,] 52.316832 -32.358168
[24,] -1100.005089 52.316832
[25,] 156.828352 -1100.005089
[26,] -563.471648 156.828352
[27,] -310.696648 -563.471648
[28,] 311.015852 -310.696648
[29,] -1035.946648 311.015852
[30,] -103.134148 -1035.946648
[31,] -612.634148 -103.134148
[32,] -700.834148 -612.634148
[33,] -1456.435896 -700.834148
[34,] -854.198396 -1456.435896
[35,] -1528.323396 -854.198396
[36,] -1882.245317 -1528.323396
[37,] -528.211877 -1882.245317
[38,] -1822.311877 -528.211877
[39,] -947.736877 -1822.311877
[40,] -678.124377 -947.736877
[41,] -2091.586877 -678.124377
[42,] -748.774377 -2091.586877
[43,] -1885.174377 -748.774377
[44,] -1286.074377 -1885.174377
[45,] -858.962142 -1286.074377
[46,] -1175.724642 -858.962142
[47,] -2206.849642 -1175.724642
[48,] 176.128437 -2206.849642
[49,] -824.538123 176.128437
[50,] -556.038123 -824.538123
[51,] -208.263123 -556.038123
[52,] -924.950623 -208.263123
[53,] -828.113123 -924.950623
[54,] -41.200623 -828.113123
[55,] -1655.300623 -41.200623
[56,] -843.300623 -1655.300623
[57,] -567.102371 -843.300623
[58,] -636.164871 -567.102371
[59,] -952.289871 -636.164871
[60,] 7.088208 -952.289871
[61,] -1232.778352 7.088208
[62,] -217.878352 -1232.778352
[63,] 334.096648 -217.878352
[64,] -1790.290852 334.096648
[65,] 476.746648 -1790.290852
[66,] 1147.259148 476.746648
[67,] -180.640852 1147.259148
[68,] -390.440852 -180.640852
[69,] 1511.757401 -390.440852
[70,] -1423.205099 1511.757401
[71,] 269.369901 -1423.205099
[72,] -324.052020 269.369901
[73,] -849.318580 -324.052020
[74,] 14.681420 -849.318580
[75,] -211.543580 14.681420
[76,] -252.331080 -211.543580
[77,] 52.506420 -252.331080
[78,] -1257.681080 52.506420
[79,] -608.581080 -1257.681080
[80,] -1053.181080 -608.581080
[81,] -502.882828 -1053.181080
[82,] -1376.645328 -502.882828
[83,] 179.029672 -1376.645328
[84,] -341.292249 179.029672
[85,] 1015.941191 -341.292249
[86,] -54.758809 1015.941191
[87,] -507.083809 -54.758809
[88,] 1572.928691 -507.083809
[89,] 897.866191 1572.928691
[90,] -296.021309 897.866191
[91,] 3703.878691 -296.021309
[92,] 3072.478691 3703.878691
[93,] 1669.476943 3072.478691
[94,] 4586.814443 1669.476943
[95,] 2256.589443 4586.814443
[96,] 2204.767522 2256.589443
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 1320.208809 511.975368
2 2349.508809 1320.208809
3 1858.483809 2349.508809
4 1735.496309 1858.483809
5 2792.333809 1735.496309
6 2081.246309 2792.333809
7 1926.946309 2081.246309
8 1733.646309 1926.946309
9 1232.344561 1733.646309
10 911.482061 1232.344561
11 1930.157061 911.482061
12 747.635140 1930.157061
13 941.868580 747.635140
14 850.268580 941.868580
15 -7.256420 850.268580
16 26.256080 -7.256420
17 -263.806420 26.256080
18 -781.693920 -263.806420
19 -688.493920 -781.693920
20 -532.293920 -688.493920
21 -1028.195668 -532.293920
22 -32.358168 -1028.195668
23 52.316832 -32.358168
24 -1100.005089 52.316832
25 156.828352 -1100.005089
26 -563.471648 156.828352
27 -310.696648 -563.471648
28 311.015852 -310.696648
29 -1035.946648 311.015852
30 -103.134148 -1035.946648
31 -612.634148 -103.134148
32 -700.834148 -612.634148
33 -1456.435896 -700.834148
34 -854.198396 -1456.435896
35 -1528.323396 -854.198396
36 -1882.245317 -1528.323396
37 -528.211877 -1882.245317
38 -1822.311877 -528.211877
39 -947.736877 -1822.311877
40 -678.124377 -947.736877
41 -2091.586877 -678.124377
42 -748.774377 -2091.586877
43 -1885.174377 -748.774377
44 -1286.074377 -1885.174377
45 -858.962142 -1286.074377
46 -1175.724642 -858.962142
47 -2206.849642 -1175.724642
48 176.128437 -2206.849642
49 -824.538123 176.128437
50 -556.038123 -824.538123
51 -208.263123 -556.038123
52 -924.950623 -208.263123
53 -828.113123 -924.950623
54 -41.200623 -828.113123
55 -1655.300623 -41.200623
56 -843.300623 -1655.300623
57 -567.102371 -843.300623
58 -636.164871 -567.102371
59 -952.289871 -636.164871
60 7.088208 -952.289871
61 -1232.778352 7.088208
62 -217.878352 -1232.778352
63 334.096648 -217.878352
64 -1790.290852 334.096648
65 476.746648 -1790.290852
66 1147.259148 476.746648
67 -180.640852 1147.259148
68 -390.440852 -180.640852
69 1511.757401 -390.440852
70 -1423.205099 1511.757401
71 269.369901 -1423.205099
72 -324.052020 269.369901
73 -849.318580 -324.052020
74 14.681420 -849.318580
75 -211.543580 14.681420
76 -252.331080 -211.543580
77 52.506420 -252.331080
78 -1257.681080 52.506420
79 -608.581080 -1257.681080
80 -1053.181080 -608.581080
81 -502.882828 -1053.181080
82 -1376.645328 -502.882828
83 179.029672 -1376.645328
84 -341.292249 179.029672
85 1015.941191 -341.292249
86 -54.758809 1015.941191
87 -507.083809 -54.758809
88 1572.928691 -507.083809
89 897.866191 1572.928691
90 -296.021309 897.866191
91 3703.878691 -296.021309
92 3072.478691 3703.878691
93 1669.476943 3072.478691
94 4586.814443 1669.476943
95 2256.589443 4586.814443
96 2204.767522 2256.589443
> 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/7xgy61227454851.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/8zm9x1227454851.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/99zlt1227454851.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
>
> #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/10m5u31227454851.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/114w941227454851.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/128qr91227454851.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/13wxuc1227454851.tab")
>
> system("convert tmp/1t2wf1227454851.ps tmp/1t2wf1227454851.png")
> system("convert tmp/2geik1227454851.ps tmp/2geik1227454851.png")
> system("convert tmp/3y2ir1227454851.ps tmp/3y2ir1227454851.png")
> system("convert tmp/43sw71227454851.ps tmp/43sw71227454851.png")
> system("convert tmp/59aha1227454851.ps tmp/59aha1227454851.png")
> system("convert tmp/687qw1227454851.ps tmp/687qw1227454851.png")
> system("convert tmp/7xgy61227454851.ps tmp/7xgy61227454851.png")
> system("convert tmp/8zm9x1227454851.ps tmp/8zm9x1227454851.png")
> system("convert tmp/99zlt1227454851.ps tmp/99zlt1227454851.png")
>
>
> proc.time()
user system elapsed
4.231 2.598 4.534