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(1221.53
+ ,2617.2
+ ,10168.52
+ ,6957.61
+ ,23448.78
+ ,1180.55
+ ,2506.13
+ ,9937.04
+ ,6688.49
+ ,23007.99
+ ,1183.26
+ ,2679.07
+ ,9202.45
+ ,6601.37
+ ,23096.32
+ ,1141.2
+ ,2589.73
+ ,9369.35
+ ,6229.02
+ ,22358.17
+ ,1049.33
+ ,2457.46
+ ,8824.06
+ ,5925.22
+ ,20536.49
+ ,1101.6
+ ,2517.3
+ ,9537.3
+ ,6147.97
+ ,21029.81
+ ,1030.71
+ ,2386.53
+ ,9382.64
+ ,5965.52
+ ,20128.99
+ ,1089.41
+ ,2453.37
+ ,9768.7
+ ,5964.33
+ ,19765.19
+ ,1186.69
+ ,2529.66
+ ,11057.4
+ ,6135.7
+ ,21108.59
+ ,1169.43
+ ,2475.14
+ ,11089.94
+ ,6153.55
+ ,21239.35
+ ,1104.49
+ ,2525.93
+ ,10126.03
+ ,5598.46
+ ,20608.7
+ ,1073.87
+ ,2480.93
+ ,10198.04
+ ,5608.79
+ ,20121.99
+ ,1115.1
+ ,2229.85
+ ,10546.44
+ ,5957.43
+ ,21872.5
+ ,1095.63
+ ,2169.14
+ ,9345.55
+ ,5625.95
+ ,21821.5
+ ,1036.19
+ ,2030.98
+ ,10034.74
+ ,5414.96
+ ,21752.87
+ ,1057.08
+ ,2071.37
+ ,10133.23
+ ,5675.16
+ ,20955.25
+ ,1020.62
+ ,1953.35
+ ,10492.53
+ ,5458.04
+ ,19724.19
+ ,987.48
+ ,1748.74
+ ,10356.83
+ ,5332.14
+ ,20573.33
+ ,919.32
+ ,1696.58
+ ,9958.44
+ ,4808.64
+ ,18378.73
+ ,919.14
+ ,1900.09
+ ,9522.5
+ ,4940.82
+ ,18171
+ ,872.81
+ ,1908.64
+ ,8828.26
+ ,4769.45
+ ,15520.99
+ ,797.87
+ ,1881.46
+ ,8109.53
+ ,4084.76
+ ,13576.02
+ ,735.09
+ ,2100.18
+ ,7568.42
+ ,3843.74
+ ,12811.57
+ ,825.88
+ ,2672.2
+ ,7994.05
+ ,4338.35
+ ,13278.21
+ ,903.25
+ ,3136
+ ,8859.56
+ ,4810.2
+ ,14387.48
+ ,896.24
+ ,2994.38
+ ,8512.27
+ ,4669.44
+ ,13888.24
+ ,968.75
+ ,3168.22
+ ,8576.98
+ ,4987.97
+ ,13968.67
+ ,1166.36
+ ,3751.41
+ ,11259.86
+ ,5831.02
+ ,18016.21
+ ,1282.83
+ ,3925.43
+ ,13072.87
+ ,6422.3
+ ,21261.89
+ ,1267.38
+ ,3719.52
+ ,13376.81
+ ,6479.56
+ ,22731.1
+ ,1280
+ ,3757.12
+ ,13481.38
+ ,6418.32
+ ,22102.01
+ ,1400.38
+ ,3722.23
+ ,14338.54
+ ,7096.79
+ ,24533.12
+ ,1385.59
+ ,4127.47
+ ,13849.99
+ ,6948.82
+ ,25755.35
+ ,1322.7
+ ,4162.5
+ ,12525.54
+ ,6534.97
+ ,22849.2
+ ,1330.63
+ ,4441.82
+ ,13603.02
+ ,6748.13
+ ,24331.67
+ ,1378.55
+ ,4325.29
+ ,13592.47
+ ,6851.75
+ ,23455.74
+ ,1468.36
+ ,4350.83
+ ,15307.78
+ ,8067.32
+ ,27812.65
+ ,1481.14
+ ,4384.47
+ ,15680.67
+ ,7870.52
+ ,28643.61
+ ,1549.38
+ ,4639.4
+ ,16737.63
+ ,8019.22
+ ,31352.58
+ ,1526.75
+ ,4697.86
+ ,16785.69
+ ,7861.51
+ ,27142.47
+ ,1473.99
+ ,4614.76
+ ,16569.09
+ ,7638.17
+ ,23984.14
+ ,1455.27
+ ,4471.65
+ ,17248.89
+ ,7584.14
+ ,23184.94
+ ,1503.35
+ ,4305.23
+ ,18138.36
+ ,8007.32
+ ,21772.73
+ ,1530.62
+ ,4433.57
+ ,17875.75
+ ,7883.04
+ ,20634.47
+ ,1482.37
+ ,4388.53
+ ,17400.41
+ ,7408.87
+ ,20318.98
+ ,1420.86
+ ,4140.3
+ ,17287.65
+ ,6917.03
+ ,19800.93
+ ,1406.82
+ ,4144.38
+ ,17604.12
+ ,6715.44
+ ,19651.51
+ ,1438.24
+ ,4070.78
+ ,17383.42
+ ,6789.11
+ ,20106.42
+ ,1418.3
+ ,3906.01
+ ,17225.83
+ ,6596.92
+ ,19964.72
+ ,1400.63
+ ,3795.91
+ ,16274.33
+ ,6309.19
+ ,18960.48
+ ,1377.94
+ ,3703.32
+ ,16399.39
+ ,6268.92
+ ,18324.35
+ ,1335.85
+ ,3675.8
+ ,16127.58
+ ,6004.33
+ ,17543.05
+ ,1303.82
+ ,3911.06
+ ,16140.76
+ ,5859.57
+ ,17392.27
+ ,1276.66
+ ,3912.28
+ ,15456.81
+ ,5681.97
+ ,16971.34
+ ,1270.2
+ ,3839.25
+ ,15505.18
+ ,5683.31
+ ,16267.62
+ ,1270.09
+ ,3744.63
+ ,15467.33
+ ,5692.86
+ ,15857.89
+ ,1310.61
+ ,3549.25
+ ,16906.23
+ ,6009.89
+ ,16661.3
+ ,1294.87
+ ,3394.14
+ ,17059.66
+ ,5970.08
+ ,15805.04
+ ,1280.66
+ ,3264.26
+ ,16205.43
+ ,5796.04
+ ,15918.48
+ ,1280.08
+ ,3328.8
+ ,16649.82
+ ,5674.15
+ ,15753.14
+ ,1248.29
+ ,3223.98
+ ,16111.43
+ ,5408.26
+ ,14876.43
+ ,1249.48
+ ,3228.01
+ ,14872.15
+ ,5193.4
+ ,14937.14
+ ,1207.01
+ ,3112.83
+ ,13606.5
+ ,4929.07
+ ,14386.37
+ ,1228.81
+ ,3051.67
+ ,13574.3
+ ,5044.12
+ ,15428.52
+ ,1220.33
+ ,3039.71
+ ,12413.6
+ ,4829.69
+ ,14903.55
+ ,1234.18
+ ,3125.67
+ ,11899.6
+ ,4886.5
+ ,14880.98
+ ,1191.33
+ ,3106.54
+ ,11584.01
+ ,4586.28
+ ,14201.06
+ ,1191.5
+ ,11276.59
+ ,4460.63
+ ,13867.07
+ ,1156.85
+ ,11008.9
+ ,4184.84
+ ,13908.97
+ ,1180.59
+ ,11668.95
+ ,4348.77
+ ,13516.88
+ ,1203.6
+ ,11740.6
+ ,4350.49
+ ,14195.35
+ ,1181.27
+ ,11387.59
+ ,4254.85
+ ,13721.69)
+ ,dim=c(5
+ ,72)
+ ,dimnames=list(c('S&P'
+ ,'Bel20'
+ ,'Nikkei'
+ ,'DAX'
+ ,'HangSeng')
+ ,1:72))
> y <- array(NA,dim=c(5,72),dimnames=list(c('S&P','Bel20','Nikkei','DAX','HangSeng'),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
S&P Bel20 Nikkei DAX HangSeng
1 1221.53 2617.20 10168.52 6957.61 23448.78
2 1180.55 2506.13 9937.04 6688.49 23007.99
3 1183.26 2679.07 9202.45 6601.37 23096.32
4 1141.20 2589.73 9369.35 6229.02 22358.17
5 1049.33 2457.46 8824.06 5925.22 20536.49
6 1101.60 2517.30 9537.30 6147.97 21029.81
7 1030.71 2386.53 9382.64 5965.52 20128.99
8 1089.41 2453.37 9768.70 5964.33 19765.19
9 1186.69 2529.66 11057.40 6135.70 21108.59
10 1169.43 2475.14 11089.94 6153.55 21239.35
11 1104.49 2525.93 10126.03 5598.46 20608.70
12 1073.87 2480.93 10198.04 5608.79 20121.99
13 1115.10 2229.85 10546.44 5957.43 21872.50
14 1095.63 2169.14 9345.55 5625.95 21821.50
15 1036.19 2030.98 10034.74 5414.96 21752.87
16 1057.08 2071.37 10133.23 5675.16 20955.25
17 1020.62 1953.35 10492.53 5458.04 19724.19
18 987.48 1748.74 10356.83 5332.14 20573.33
19 919.32 1696.58 9958.44 4808.64 18378.73
20 919.14 1900.09 9522.50 4940.82 18171.00
21 872.81 1908.64 8828.26 4769.45 15520.99
22 797.87 1881.46 8109.53 4084.76 13576.02
23 735.09 2100.18 7568.42 3843.74 12811.57
24 825.88 2672.20 7994.05 4338.35 13278.21
25 903.25 3136.00 8859.56 4810.20 14387.48
26 896.24 2994.38 8512.27 4669.44 13888.24
27 968.75 3168.22 8576.98 4987.97 13968.67
28 1166.36 3751.41 11259.86 5831.02 18016.21
29 1282.83 3925.43 13072.87 6422.30 21261.89
30 1267.38 3719.52 13376.81 6479.56 22731.10
31 1280.00 3757.12 13481.38 6418.32 22102.01
32 1400.38 3722.23 14338.54 7096.79 24533.12
33 1385.59 4127.47 13849.99 6948.82 25755.35
34 1322.70 4162.50 12525.54 6534.97 22849.20
35 1330.63 4441.82 13603.02 6748.13 24331.67
36 1378.55 4325.29 13592.47 6851.75 23455.74
37 1468.36 4350.83 15307.78 8067.32 27812.65
38 1481.14 4384.47 15680.67 7870.52 28643.61
39 1549.38 4639.40 16737.63 8019.22 31352.58
40 1526.75 4697.86 16785.69 7861.51 27142.47
41 1473.99 4614.76 16569.09 7638.17 23984.14
42 1455.27 4471.65 17248.89 7584.14 23184.94
43 1503.35 4305.23 18138.36 8007.32 21772.73
44 1530.62 4433.57 17875.75 7883.04 20634.47
45 1482.37 4388.53 17400.41 7408.87 20318.98
46 1420.86 4140.30 17287.65 6917.03 19800.93
47 1406.82 4144.38 17604.12 6715.44 19651.51
48 1438.24 4070.78 17383.42 6789.11 20106.42
49 1418.30 3906.01 17225.83 6596.92 19964.72
50 1400.63 3795.91 16274.33 6309.19 18960.48
51 1377.94 3703.32 16399.39 6268.92 18324.35
52 1335.85 3675.80 16127.58 6004.33 17543.05
53 1303.82 3911.06 16140.76 5859.57 17392.27
54 1276.66 3912.28 15456.81 5681.97 16971.34
55 1270.20 3839.25 15505.18 5683.31 16267.62
56 1270.09 3744.63 15467.33 5692.86 15857.89
57 1310.61 3549.25 16906.23 6009.89 16661.30
58 1294.87 3394.14 17059.66 5970.08 15805.04
59 1280.66 3264.26 16205.43 5796.04 15918.48
60 1280.08 3328.80 16649.82 5674.15 15753.14
61 1248.29 3223.98 16111.43 5408.26 14876.43
62 1249.48 3228.01 14872.15 5193.40 14937.14
63 1207.01 3112.83 13606.50 4929.07 14386.37
64 1228.81 3051.67 13574.30 5044.12 15428.52
65 1220.33 3039.71 12413.60 4829.69 14903.55
66 1234.18 3125.67 11899.60 4886.50 14880.98
67 1191.33 3106.54 11584.01 4586.28 14201.06
68 1191.50 11276.59 4460.63 13867.07 1156.85
69 11008.90 4184.84 13908.97 1180.59 11668.95
70 4348.77 13516.88 1203.60 11740.60 4350.49
71 14195.35 1181.27 11387.59 4254.85 13721.69
72 1221.53 2617.20 10168.52 6957.61 23448.78
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Bel20 Nikkei DAX HangSeng
3.902e+03 6.792e-01 9.317e-03 -8.159e-01 1.112e-02
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-1669.8 -702.9 -232.1 251.0 12704.3
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.902e+03 1.240e+03 3.145 0.00247 **
Bel20 6.792e-01 2.630e-01 2.583 0.01199 *
Nikkei 9.317e-03 6.719e-02 0.139 0.89013
DAX -8.159e-01 2.688e-01 -3.035 0.00342 **
HangSeng 1.112e-02 6.397e-02 0.174 0.86256
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1813 on 67 degrees of freedom
Multiple R-squared: 0.1762, Adjusted R-squared: 0.127
F-statistic: 3.583 on 4 and 67 DF, p-value: 0.01044
> 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,] 5.362733e-07 1.072547e-06 0.9999995
[2,] 3.155278e-09 6.310556e-09 1.0000000
[3,] 1.794940e-11 3.589881e-11 1.0000000
[4,] 1.041310e-13 2.082620e-13 1.0000000
[5,] 1.691555e-15 3.383109e-15 1.0000000
[6,] 2.705269e-17 5.410538e-17 1.0000000
[7,] 5.565783e-18 1.113157e-17 1.0000000
[8,] 1.844689e-19 3.689379e-19 1.0000000
[9,] 1.630588e-21 3.261175e-21 1.0000000
[10,] 1.386875e-23 2.773751e-23 1.0000000
[11,] 1.358553e-25 2.717106e-25 1.0000000
[12,] 1.325442e-27 2.650885e-27 1.0000000
[13,] 9.821959e-30 1.964392e-29 1.0000000
[14,] 4.448040e-31 8.896081e-31 1.0000000
[15,] 1.048737e-32 2.097474e-32 1.0000000
[16,] 2.627507e-34 5.255014e-34 1.0000000
[17,] 9.171115e-36 1.834223e-35 1.0000000
[18,] 1.058024e-36 2.116047e-36 1.0000000
[19,] 1.413650e-38 2.827300e-38 1.0000000
[20,] 6.278030e-40 1.255606e-39 1.0000000
[21,] 9.199733e-42 1.839947e-41 1.0000000
[22,] 2.792693e-43 5.585386e-43 1.0000000
[23,] 2.268699e-44 4.537398e-44 1.0000000
[24,] 2.790795e-46 5.581590e-46 1.0000000
[25,] 2.774673e-48 5.549345e-48 1.0000000
[26,] 4.019439e-50 8.038877e-50 1.0000000
[27,] 5.440000e-52 1.088000e-51 1.0000000
[28,] 3.810086e-53 7.620171e-53 1.0000000
[29,] 4.540850e-55 9.081699e-55 1.0000000
[30,] 1.961343e-54 3.922685e-54 1.0000000
[31,] 1.884622e-55 3.769243e-55 1.0000000
[32,] 8.513555e-57 1.702711e-56 1.0000000
[33,] 1.211381e-58 2.422761e-58 1.0000000
[34,] 1.542763e-60 3.085525e-60 1.0000000
[35,] 2.263655e-62 4.527310e-62 1.0000000
[36,] 4.597059e-64 9.194117e-64 1.0000000
[37,] 3.059018e-65 6.118036e-65 1.0000000
[38,] 1.639205e-66 3.278410e-66 1.0000000
[39,] 6.341946e-68 1.268389e-67 1.0000000
[40,] 2.131396e-69 4.262792e-69 1.0000000
[41,] 2.832178e-70 5.664356e-70 1.0000000
[42,] 4.606906e-71 9.213812e-71 1.0000000
[43,] 4.269537e-71 8.539074e-71 1.0000000
[44,] 5.970329e-72 1.194066e-71 1.0000000
[45,] 3.090157e-73 6.180314e-73 1.0000000
[46,] 6.431413e-75 1.286283e-74 1.0000000
[47,] 1.150448e-76 2.300896e-76 1.0000000
[48,] 1.666358e-78 3.332716e-78 1.0000000
[49,] 2.376178e-80 4.752356e-80 1.0000000
[50,] 6.660091e-82 1.332018e-81 1.0000000
[51,] 2.294558e-83 4.589115e-83 1.0000000
[52,] 7.232105e-85 1.446421e-84 1.0000000
[53,] 1.540571e-85 3.081142e-85 1.0000000
[54,] 8.466895e-86 1.693379e-85 1.0000000
[55,] 3.367485e-85 6.734971e-85 1.0000000
[56,] 1.165344e-85 2.330688e-85 1.0000000
[57,] 7.020471e-86 1.404094e-85 1.0000000
> postscript(file="/var/www/html/rcomp/tmp/1u0rc1293227494.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/www/html/rcomp/tmp/2u0rc1293227494.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/www/html/rcomp/tmp/3u0rc1293227494.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/www/html/rcomp/tmp/4n98f1293227494.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/www/html/rcomp/tmp/5n98f1293227494.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 = 72
Frequency = 1
1 2 3 4 5 6
863.57488 685.51543 505.54517 227.01514 2.44333 183.68259
7 8 9 10 11 12
64.20554 76.98393 235.32822 267.90542 -268.43818 -255.32637
13 14 15 16 17 18
218.18814 -18.74657 -162.15294 51.55016 -71.56147 -76.62661
19 20 21 22 23 24
-508.37576 -532.56378 -688.59557 -1275.39650 -1669.84063 -1573.16840
25 26 27 28 29 30
-1446.22469 -1463.10684 -1250.27769 -830.91644 -403.18704 -251.22785
31 32 33 34 35 36
-308.09276 354.53759 -65.25587 -444.95350 -479.33921 -257.89289
37 38 39 40 41 42
741.94035 558.59156 535.04584 390.38670 248.97193 285.91993
43 44 45 46 47 48
799.71668 653.51762 256.91928 -30.47531 -213.05129 -74.53527
49 50 51 52 53 54
-136.32753 -293.94795 -280.70073 -508.76058 -817.13611 -978.97717
55 56 57 58 59 60
-927.36969 -850.51444 -440.96402 -375.74507 -437.04178 -583.20994
61 62 63 64 65 66
-745.98383 -911.96394 -1073.95590 -928.03132 -1086.69174 -1079.83491
67 68 69 70 71 72
-1344.14285 890.47747 4968.82655 785.99455 12704.26999 863.57488
> postscript(file="/var/www/html/rcomp/tmp/6n98f1293227494.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 = 72
Frequency = 1
lag(myerror, k = 1) myerror
0 863.57488 NA
1 685.51543 863.57488
2 505.54517 685.51543
3 227.01514 505.54517
4 2.44333 227.01514
5 183.68259 2.44333
6 64.20554 183.68259
7 76.98393 64.20554
8 235.32822 76.98393
9 267.90542 235.32822
10 -268.43818 267.90542
11 -255.32637 -268.43818
12 218.18814 -255.32637
13 -18.74657 218.18814
14 -162.15294 -18.74657
15 51.55016 -162.15294
16 -71.56147 51.55016
17 -76.62661 -71.56147
18 -508.37576 -76.62661
19 -532.56378 -508.37576
20 -688.59557 -532.56378
21 -1275.39650 -688.59557
22 -1669.84063 -1275.39650
23 -1573.16840 -1669.84063
24 -1446.22469 -1573.16840
25 -1463.10684 -1446.22469
26 -1250.27769 -1463.10684
27 -830.91644 -1250.27769
28 -403.18704 -830.91644
29 -251.22785 -403.18704
30 -308.09276 -251.22785
31 354.53759 -308.09276
32 -65.25587 354.53759
33 -444.95350 -65.25587
34 -479.33921 -444.95350
35 -257.89289 -479.33921
36 741.94035 -257.89289
37 558.59156 741.94035
38 535.04584 558.59156
39 390.38670 535.04584
40 248.97193 390.38670
41 285.91993 248.97193
42 799.71668 285.91993
43 653.51762 799.71668
44 256.91928 653.51762
45 -30.47531 256.91928
46 -213.05129 -30.47531
47 -74.53527 -213.05129
48 -136.32753 -74.53527
49 -293.94795 -136.32753
50 -280.70073 -293.94795
51 -508.76058 -280.70073
52 -817.13611 -508.76058
53 -978.97717 -817.13611
54 -927.36969 -978.97717
55 -850.51444 -927.36969
56 -440.96402 -850.51444
57 -375.74507 -440.96402
58 -437.04178 -375.74507
59 -583.20994 -437.04178
60 -745.98383 -583.20994
61 -911.96394 -745.98383
62 -1073.95590 -911.96394
63 -928.03132 -1073.95590
64 -1086.69174 -928.03132
65 -1079.83491 -1086.69174
66 -1344.14285 -1079.83491
67 890.47747 -1344.14285
68 4968.82655 890.47747
69 785.99455 4968.82655
70 12704.26999 785.99455
71 863.57488 12704.26999
72 NA 863.57488
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 685.51543 863.57488
[2,] 505.54517 685.51543
[3,] 227.01514 505.54517
[4,] 2.44333 227.01514
[5,] 183.68259 2.44333
[6,] 64.20554 183.68259
[7,] 76.98393 64.20554
[8,] 235.32822 76.98393
[9,] 267.90542 235.32822
[10,] -268.43818 267.90542
[11,] -255.32637 -268.43818
[12,] 218.18814 -255.32637
[13,] -18.74657 218.18814
[14,] -162.15294 -18.74657
[15,] 51.55016 -162.15294
[16,] -71.56147 51.55016
[17,] -76.62661 -71.56147
[18,] -508.37576 -76.62661
[19,] -532.56378 -508.37576
[20,] -688.59557 -532.56378
[21,] -1275.39650 -688.59557
[22,] -1669.84063 -1275.39650
[23,] -1573.16840 -1669.84063
[24,] -1446.22469 -1573.16840
[25,] -1463.10684 -1446.22469
[26,] -1250.27769 -1463.10684
[27,] -830.91644 -1250.27769
[28,] -403.18704 -830.91644
[29,] -251.22785 -403.18704
[30,] -308.09276 -251.22785
[31,] 354.53759 -308.09276
[32,] -65.25587 354.53759
[33,] -444.95350 -65.25587
[34,] -479.33921 -444.95350
[35,] -257.89289 -479.33921
[36,] 741.94035 -257.89289
[37,] 558.59156 741.94035
[38,] 535.04584 558.59156
[39,] 390.38670 535.04584
[40,] 248.97193 390.38670
[41,] 285.91993 248.97193
[42,] 799.71668 285.91993
[43,] 653.51762 799.71668
[44,] 256.91928 653.51762
[45,] -30.47531 256.91928
[46,] -213.05129 -30.47531
[47,] -74.53527 -213.05129
[48,] -136.32753 -74.53527
[49,] -293.94795 -136.32753
[50,] -280.70073 -293.94795
[51,] -508.76058 -280.70073
[52,] -817.13611 -508.76058
[53,] -978.97717 -817.13611
[54,] -927.36969 -978.97717
[55,] -850.51444 -927.36969
[56,] -440.96402 -850.51444
[57,] -375.74507 -440.96402
[58,] -437.04178 -375.74507
[59,] -583.20994 -437.04178
[60,] -745.98383 -583.20994
[61,] -911.96394 -745.98383
[62,] -1073.95590 -911.96394
[63,] -928.03132 -1073.95590
[64,] -1086.69174 -928.03132
[65,] -1079.83491 -1086.69174
[66,] -1344.14285 -1079.83491
[67,] 890.47747 -1344.14285
[68,] 4968.82655 890.47747
[69,] 785.99455 4968.82655
[70,] 12704.26999 785.99455
[71,] 863.57488 12704.26999
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 685.51543 863.57488
2 505.54517 685.51543
3 227.01514 505.54517
4 2.44333 227.01514
5 183.68259 2.44333
6 64.20554 183.68259
7 76.98393 64.20554
8 235.32822 76.98393
9 267.90542 235.32822
10 -268.43818 267.90542
11 -255.32637 -268.43818
12 218.18814 -255.32637
13 -18.74657 218.18814
14 -162.15294 -18.74657
15 51.55016 -162.15294
16 -71.56147 51.55016
17 -76.62661 -71.56147
18 -508.37576 -76.62661
19 -532.56378 -508.37576
20 -688.59557 -532.56378
21 -1275.39650 -688.59557
22 -1669.84063 -1275.39650
23 -1573.16840 -1669.84063
24 -1446.22469 -1573.16840
25 -1463.10684 -1446.22469
26 -1250.27769 -1463.10684
27 -830.91644 -1250.27769
28 -403.18704 -830.91644
29 -251.22785 -403.18704
30 -308.09276 -251.22785
31 354.53759 -308.09276
32 -65.25587 354.53759
33 -444.95350 -65.25587
34 -479.33921 -444.95350
35 -257.89289 -479.33921
36 741.94035 -257.89289
37 558.59156 741.94035
38 535.04584 558.59156
39 390.38670 535.04584
40 248.97193 390.38670
41 285.91993 248.97193
42 799.71668 285.91993
43 653.51762 799.71668
44 256.91928 653.51762
45 -30.47531 256.91928
46 -213.05129 -30.47531
47 -74.53527 -213.05129
48 -136.32753 -74.53527
49 -293.94795 -136.32753
50 -280.70073 -293.94795
51 -508.76058 -280.70073
52 -817.13611 -508.76058
53 -978.97717 -817.13611
54 -927.36969 -978.97717
55 -850.51444 -927.36969
56 -440.96402 -850.51444
57 -375.74507 -440.96402
58 -437.04178 -375.74507
59 -583.20994 -437.04178
60 -745.98383 -583.20994
61 -911.96394 -745.98383
62 -1073.95590 -911.96394
63 -928.03132 -1073.95590
64 -1086.69174 -928.03132
65 -1079.83491 -1086.69174
66 -1344.14285 -1079.83491
67 890.47747 -1344.14285
68 4968.82655 890.47747
69 785.99455 4968.82655
70 12704.26999 785.99455
71 863.57488 12704.26999
> 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/7gjqi1293227494.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/www/html/rcomp/tmp/88s6l1293227494.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/www/html/rcomp/tmp/98s6l1293227494.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/www/html/rcomp/tmp/108s6l1293227494.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/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/11n14t1293227494.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/12xbmw1293227494.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/134c1q1293227494.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/14xlit1293227494.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/15i4gz1293227494.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/16wve81293227494.tab")
+ }
>
> try(system("convert tmp/1u0rc1293227494.ps tmp/1u0rc1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/2u0rc1293227494.ps tmp/2u0rc1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/3u0rc1293227494.ps tmp/3u0rc1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/4n98f1293227494.ps tmp/4n98f1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/5n98f1293227494.ps tmp/5n98f1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/6n98f1293227494.ps tmp/6n98f1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/7gjqi1293227494.ps tmp/7gjqi1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/88s6l1293227494.ps tmp/88s6l1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/98s6l1293227494.ps tmp/98s6l1293227494.png",intern=TRUE))
character(0)
> try(system("convert tmp/108s6l1293227494.ps tmp/108s6l1293227494.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.718 1.709 6.309