R version 2.15.2 (2012-10-26) -- "Trick or Treat"
Copyright (C) 2012 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: i686-pc-linux-gnu (32-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> x <- array(list(12.008,9.169,8.788,8.417,8.247,8.197,8.236,8.253,7.733,8.366,8.626,8.863,10.102,8.463,9.114,8.563,8.872,8.301,8.301,8.278,7.736,7.973,8.268,9.476,11.100,8.962,9.173,8.738,8.459,8.078,8.411,8.291,7.810,8.616,8.312,9.692,9.911,8.915,9.452,9.112,8.472,8.230,8.384,8.625,8.221,8.649,8.625,10.443,10.357,8.586,8.892,8.329,8.101,7.922,8.120,7.838,7.735,8.406,8.209,9.451,10.041,9.411,10.405,8.467,8.464,8.102,7.627,7.513,7.510,8.291,8.064,9.383,9.706,8.579,9.474,8.318,8.213,8.059,9.111,7.708,7.680,8.014,8.007,8.718,9.486,9.113,9.025,8.476,7.952,7.759,7.835,7.600,7.651,8.319,8.812,8.630),dim=c(1,96),dimnames=list(c('y'),1:96))
> y <- array(NA,dim=c(1,96),dimnames=list(c('y'),1:96))
> 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'
> 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)
> library(lmtest)
Loading required package: zoo
Attaching package: 'zoo'
The following object(s) are masked from 'package:base':
as.Date, 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 M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 t
1 12.008 1 0 0 0 0 0 0 0 0 0 0 1
2 9.169 0 1 0 0 0 0 0 0 0 0 0 2
3 8.788 0 0 1 0 0 0 0 0 0 0 0 3
4 8.417 0 0 0 1 0 0 0 0 0 0 0 4
5 8.247 0 0 0 0 1 0 0 0 0 0 0 5
6 8.197 0 0 0 0 0 1 0 0 0 0 0 6
7 8.236 0 0 0 0 0 0 1 0 0 0 0 7
8 8.253 0 0 0 0 0 0 0 1 0 0 0 8
9 7.733 0 0 0 0 0 0 0 0 1 0 0 9
10 8.366 0 0 0 0 0 0 0 0 0 1 0 10
11 8.626 0 0 0 0 0 0 0 0 0 0 1 11
12 8.863 0 0 0 0 0 0 0 0 0 0 0 12
13 10.102 1 0 0 0 0 0 0 0 0 0 0 13
14 8.463 0 1 0 0 0 0 0 0 0 0 0 14
15 9.114 0 0 1 0 0 0 0 0 0 0 0 15
16 8.563 0 0 0 1 0 0 0 0 0 0 0 16
17 8.872 0 0 0 0 1 0 0 0 0 0 0 17
18 8.301 0 0 0 0 0 1 0 0 0 0 0 18
19 8.301 0 0 0 0 0 0 1 0 0 0 0 19
20 8.278 0 0 0 0 0 0 0 1 0 0 0 20
21 7.736 0 0 0 0 0 0 0 0 1 0 0 21
22 7.973 0 0 0 0 0 0 0 0 0 1 0 22
23 8.268 0 0 0 0 0 0 0 0 0 0 1 23
24 9.476 0 0 0 0 0 0 0 0 0 0 0 24
25 11.100 1 0 0 0 0 0 0 0 0 0 0 25
26 8.962 0 1 0 0 0 0 0 0 0 0 0 26
27 9.173 0 0 1 0 0 0 0 0 0 0 0 27
28 8.738 0 0 0 1 0 0 0 0 0 0 0 28
29 8.459 0 0 0 0 1 0 0 0 0 0 0 29
30 8.078 0 0 0 0 0 1 0 0 0 0 0 30
31 8.411 0 0 0 0 0 0 1 0 0 0 0 31
32 8.291 0 0 0 0 0 0 0 1 0 0 0 32
33 7.810 0 0 0 0 0 0 0 0 1 0 0 33
34 8.616 0 0 0 0 0 0 0 0 0 1 0 34
35 8.312 0 0 0 0 0 0 0 0 0 0 1 35
36 9.692 0 0 0 0 0 0 0 0 0 0 0 36
37 9.911 1 0 0 0 0 0 0 0 0 0 0 37
38 8.915 0 1 0 0 0 0 0 0 0 0 0 38
39 9.452 0 0 1 0 0 0 0 0 0 0 0 39
40 9.112 0 0 0 1 0 0 0 0 0 0 0 40
41 8.472 0 0 0 0 1 0 0 0 0 0 0 41
42 8.230 0 0 0 0 0 1 0 0 0 0 0 42
43 8.384 0 0 0 0 0 0 1 0 0 0 0 43
44 8.625 0 0 0 0 0 0 0 1 0 0 0 44
45 8.221 0 0 0 0 0 0 0 0 1 0 0 45
46 8.649 0 0 0 0 0 0 0 0 0 1 0 46
47 8.625 0 0 0 0 0 0 0 0 0 0 1 47
48 10.443 0 0 0 0 0 0 0 0 0 0 0 48
49 10.357 1 0 0 0 0 0 0 0 0 0 0 49
50 8.586 0 1 0 0 0 0 0 0 0 0 0 50
51 8.892 0 0 1 0 0 0 0 0 0 0 0 51
52 8.329 0 0 0 1 0 0 0 0 0 0 0 52
53 8.101 0 0 0 0 1 0 0 0 0 0 0 53
54 7.922 0 0 0 0 0 1 0 0 0 0 0 54
55 8.120 0 0 0 0 0 0 1 0 0 0 0 55
56 7.838 0 0 0 0 0 0 0 1 0 0 0 56
57 7.735 0 0 0 0 0 0 0 0 1 0 0 57
58 8.406 0 0 0 0 0 0 0 0 0 1 0 58
59 8.209 0 0 0 0 0 0 0 0 0 0 1 59
60 9.451 0 0 0 0 0 0 0 0 0 0 0 60
61 10.041 1 0 0 0 0 0 0 0 0 0 0 61
62 9.411 0 1 0 0 0 0 0 0 0 0 0 62
63 10.405 0 0 1 0 0 0 0 0 0 0 0 63
64 8.467 0 0 0 1 0 0 0 0 0 0 0 64
65 8.464 0 0 0 0 1 0 0 0 0 0 0 65
66 8.102 0 0 0 0 0 1 0 0 0 0 0 66
67 7.627 0 0 0 0 0 0 1 0 0 0 0 67
68 7.513 0 0 0 0 0 0 0 1 0 0 0 68
69 7.510 0 0 0 0 0 0 0 0 1 0 0 69
70 8.291 0 0 0 0 0 0 0 0 0 1 0 70
71 8.064 0 0 0 0 0 0 0 0 0 0 1 71
72 9.383 0 0 0 0 0 0 0 0 0 0 0 72
73 9.706 1 0 0 0 0 0 0 0 0 0 0 73
74 8.579 0 1 0 0 0 0 0 0 0 0 0 74
75 9.474 0 0 1 0 0 0 0 0 0 0 0 75
76 8.318 0 0 0 1 0 0 0 0 0 0 0 76
77 8.213 0 0 0 0 1 0 0 0 0 0 0 77
78 8.059 0 0 0 0 0 1 0 0 0 0 0 78
79 9.111 0 0 0 0 0 0 1 0 0 0 0 79
80 7.708 0 0 0 0 0 0 0 1 0 0 0 80
81 7.680 0 0 0 0 0 0 0 0 1 0 0 81
82 8.014 0 0 0 0 0 0 0 0 0 1 0 82
83 8.007 0 0 0 0 0 0 0 0 0 0 1 83
84 8.718 0 0 0 0 0 0 0 0 0 0 0 84
85 9.486 1 0 0 0 0 0 0 0 0 0 0 85
86 9.113 0 1 0 0 0 0 0 0 0 0 0 86
87 9.025 0 0 1 0 0 0 0 0 0 0 0 87
88 8.476 0 0 0 1 0 0 0 0 0 0 0 88
89 7.952 0 0 0 0 1 0 0 0 0 0 0 89
90 7.759 0 0 0 0 0 1 0 0 0 0 0 90
91 7.835 0 0 0 0 0 0 1 0 0 0 0 91
92 7.600 0 0 0 0 0 0 0 1 0 0 0 92
93 7.651 0 0 0 0 0 0 0 0 1 0 0 93
94 8.319 0 0 0 0 0 0 0 0 0 1 0 94
95 8.812 0 0 0 0 0 0 0 0 0 0 1 95
96 8.630 0 0 0 0 0 0 0 0 0 0 0 96
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) M1 M2 M3 M4 M5
9.560571 0.960314 -0.474578 -0.079720 -0.813362 -1.014130
M6 M7 M8 M9 M10 M11
-1.276397 -1.100039 -1.335681 -1.585198 -1.011216 -0.970858
t
-0.004233
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.68015 -0.21857 -0.01982 0.12450 1.49135
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 9.560571 0.164960 57.957 < 2e-16 ***
M1 0.960314 0.203618 4.716 9.59e-06 ***
M2 -0.474578 0.203501 -2.332 0.022120 *
M3 -0.079720 0.203395 -0.392 0.696101
M4 -0.813362 0.203300 -4.001 0.000136 ***
M5 -1.014130 0.203216 -4.990 3.27e-06 ***
M6 -1.276397 0.203144 -6.283 1.46e-08 ***
M7 -1.100039 0.203082 -5.417 5.80e-07 ***
M8 -1.335681 0.203032 -6.579 3.99e-09 ***
M9 -1.585198 0.202993 -7.809 1.56e-11 ***
M10 -1.011216 0.202965 -4.982 3.38e-06 ***
M11 -0.970858 0.202948 -4.784 7.38e-06 ***
t -0.004233 0.001507 -2.809 0.006187 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.4059 on 83 degrees of freedom
Multiple R-squared: 0.7757, Adjusted R-squared: 0.7433
F-statistic: 23.93 on 12 and 83 DF, p-value: < 2.2e-16
> 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.9870200 0.02595997 0.012979986
[2,] 0.9932155 0.01356905 0.006784524
[3,] 0.9871334 0.02573330 0.012866648
[4,] 0.9767903 0.04641946 0.023209728
[5,] 0.9594030 0.08119398 0.040596990
[6,] 0.9354202 0.12915959 0.064579796
[7,] 0.9201474 0.15970527 0.079852635
[8,] 0.8884826 0.22303482 0.111517411
[9,] 0.9009656 0.19806887 0.099034435
[10,] 0.9058265 0.18834699 0.094173493
[11,] 0.8759878 0.24802450 0.124012248
[12,] 0.8605860 0.27882794 0.139413968
[13,] 0.8228778 0.35424443 0.177122217
[14,] 0.7674498 0.46510047 0.232550234
[15,] 0.7124055 0.57518896 0.287594482
[16,] 0.6500217 0.69995658 0.349978292
[17,] 0.5786927 0.84261470 0.421307348
[18,] 0.5126374 0.97472523 0.487362615
[19,] 0.4848100 0.96961997 0.515190015
[20,] 0.4323147 0.86462943 0.567685286
[21,] 0.4170191 0.83403821 0.582980897
[22,] 0.6568587 0.68628253 0.343141266
[23,] 0.6032830 0.79343397 0.396716984
[24,] 0.5892885 0.82142298 0.410711490
[25,] 0.5977951 0.80440971 0.402204853
[26,] 0.5301900 0.93962008 0.469810040
[27,] 0.4608798 0.92175965 0.539120173
[28,] 0.3941805 0.78836095 0.605819524
[29,] 0.4108589 0.82171771 0.589141143
[30,] 0.3897491 0.77949815 0.610250927
[31,] 0.3400590 0.68011798 0.659941008
[32,] 0.2847275 0.56945508 0.715272460
[33,] 0.6277139 0.74457222 0.372286108
[34,] 0.6649383 0.67012346 0.335061729
[35,] 0.6756024 0.64879515 0.324397577
[36,] 0.7619573 0.47608533 0.238042663
[37,] 0.7414211 0.51715778 0.258578888
[38,] 0.7229456 0.55410874 0.277054372
[39,] 0.6844361 0.63112777 0.315563884
[40,] 0.6420587 0.71588262 0.357941310
[41,] 0.6094836 0.78103285 0.390516424
[42,] 0.5436814 0.91263729 0.456318645
[43,] 0.4730491 0.94609821 0.526950893
[44,] 0.4315821 0.86316416 0.568417920
[45,] 0.3842568 0.76851353 0.615743237
[46,] 0.3750225 0.75004504 0.624977481
[47,] 0.3971447 0.79428932 0.602855340
[48,] 0.7663056 0.46738885 0.233694426
[49,] 0.7052513 0.58949734 0.294748672
[50,] 0.6657682 0.66846361 0.334231807
[51,] 0.6004337 0.79913261 0.399566305
[52,] 0.7542784 0.49144324 0.245721619
[53,] 0.7311383 0.53772333 0.268861665
[54,] 0.6811621 0.63767578 0.318837892
[55,] 0.5979370 0.80412606 0.402063032
[56,] 0.6017243 0.79655139 0.398275696
[57,] 0.6045651 0.79086984 0.395434920
[58,] 0.5534951 0.89300985 0.446504924
[59,] 0.5612802 0.87743969 0.438719845
[60,] 0.4929242 0.98584833 0.507075837
[61,] 0.4019149 0.80382980 0.598085101
[62,] 0.2982534 0.59650687 0.701746566
[63,] 0.2083213 0.41664252 0.791678739
[64,] 0.7786636 0.44267286 0.221336432
[65,] 0.6886518 0.62269642 0.311348211
> postscript(file="/var/wessaorg/rcomp/tmp/13tsy1354050199.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index')
> points(x[,1]-mysum$resid)
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/27sek1354050199.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/3atct1354050199.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/4so4l1354050199.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals')
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/5lnnu1354050199.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 = 96
Frequency = 1
1 2 3 4 5
1.4913472222 0.0914722222 -0.6801527778 -0.3132777778 -0.2782777778
6 7 8 9 10
-0.0617777778 -0.1949027778 0.0619722222 -0.2042777778 -0.1410277778
11 12 13 14 15
0.0828472222 -0.6467777778 -0.3638591270 -0.5637341270 -0.3033591270
16 17 18 19 20
-0.1164841270 0.3975158730 0.0930158730 -0.0791091270 0.1377658730
21 22 23 24 25
-0.1504841270 -0.4832341270 -0.2243591270 0.0170158730 0.6849345238
26 27 28 29 30
-0.0139404762 -0.1935654762 0.1093095238 0.0353095238 -0.0791904762
31 32 33 34 35
0.0816845238 0.2015595238 -0.0256904762 0.2105595238 -0.1295654762
36 37 38 39 40
0.2838095238 -0.4532718254 -0.0101468254 0.1362281746 0.5341031746
41 42 43 44 45
0.0991031746 0.1236031746 0.1054781746 0.5863531746 0.4361031746
46 47 48 49 50
0.2943531746 0.2342281746 1.0856031746 0.0435218254 -0.2883531746
51 52 53 54 55
-0.3729781746 -0.1981031746 -0.2211031746 -0.1336031746 -0.1077281746
56 57 58 59 60
-0.1498531746 0.0008968254 0.1021468254 -0.1309781746 0.1443968254
61 62 63 64 65
-0.2216845238 0.5874404762 1.1908154762 -0.0093095238 0.1926904762
66 67 68 69 70
0.0971904762 -0.5499345238 -0.4240595238 -0.1733095238 0.0379404762
71 72 73 74 75
-0.2251845238 0.1271904762 -0.5058908730 -0.1937658730 0.3106091270
76 77 78 79 80
-0.1075158730 -0.0075158730 0.1049841270 0.9848591270 -0.1782658730
81 82 83 84 85
0.0474841270 -0.1882658730 -0.2313908730 -0.4870158730 -0.6750972222
86 87 88 89 90
0.3910277778 -0.0875972222 0.1012777778 -0.2177222222 -0.1442222222
91 92 93 94 95
-0.2403472222 -0.2354722222 0.0692777778 0.1675277778 0.6244027778
96
-0.5242222222
> postscript(file="/var/wessaorg/rcomp/tmp/6zrzc1354050199.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 = 96
Frequency = 1
lag(myerror, k = 1) myerror
0 1.4913472222 NA
1 0.0914722222 1.4913472222
2 -0.6801527778 0.0914722222
3 -0.3132777778 -0.6801527778
4 -0.2782777778 -0.3132777778
5 -0.0617777778 -0.2782777778
6 -0.1949027778 -0.0617777778
7 0.0619722222 -0.1949027778
8 -0.2042777778 0.0619722222
9 -0.1410277778 -0.2042777778
10 0.0828472222 -0.1410277778
11 -0.6467777778 0.0828472222
12 -0.3638591270 -0.6467777778
13 -0.5637341270 -0.3638591270
14 -0.3033591270 -0.5637341270
15 -0.1164841270 -0.3033591270
16 0.3975158730 -0.1164841270
17 0.0930158730 0.3975158730
18 -0.0791091270 0.0930158730
19 0.1377658730 -0.0791091270
20 -0.1504841270 0.1377658730
21 -0.4832341270 -0.1504841270
22 -0.2243591270 -0.4832341270
23 0.0170158730 -0.2243591270
24 0.6849345238 0.0170158730
25 -0.0139404762 0.6849345238
26 -0.1935654762 -0.0139404762
27 0.1093095238 -0.1935654762
28 0.0353095238 0.1093095238
29 -0.0791904762 0.0353095238
30 0.0816845238 -0.0791904762
31 0.2015595238 0.0816845238
32 -0.0256904762 0.2015595238
33 0.2105595238 -0.0256904762
34 -0.1295654762 0.2105595238
35 0.2838095238 -0.1295654762
36 -0.4532718254 0.2838095238
37 -0.0101468254 -0.4532718254
38 0.1362281746 -0.0101468254
39 0.5341031746 0.1362281746
40 0.0991031746 0.5341031746
41 0.1236031746 0.0991031746
42 0.1054781746 0.1236031746
43 0.5863531746 0.1054781746
44 0.4361031746 0.5863531746
45 0.2943531746 0.4361031746
46 0.2342281746 0.2943531746
47 1.0856031746 0.2342281746
48 0.0435218254 1.0856031746
49 -0.2883531746 0.0435218254
50 -0.3729781746 -0.2883531746
51 -0.1981031746 -0.3729781746
52 -0.2211031746 -0.1981031746
53 -0.1336031746 -0.2211031746
54 -0.1077281746 -0.1336031746
55 -0.1498531746 -0.1077281746
56 0.0008968254 -0.1498531746
57 0.1021468254 0.0008968254
58 -0.1309781746 0.1021468254
59 0.1443968254 -0.1309781746
60 -0.2216845238 0.1443968254
61 0.5874404762 -0.2216845238
62 1.1908154762 0.5874404762
63 -0.0093095238 1.1908154762
64 0.1926904762 -0.0093095238
65 0.0971904762 0.1926904762
66 -0.5499345238 0.0971904762
67 -0.4240595238 -0.5499345238
68 -0.1733095238 -0.4240595238
69 0.0379404762 -0.1733095238
70 -0.2251845238 0.0379404762
71 0.1271904762 -0.2251845238
72 -0.5058908730 0.1271904762
73 -0.1937658730 -0.5058908730
74 0.3106091270 -0.1937658730
75 -0.1075158730 0.3106091270
76 -0.0075158730 -0.1075158730
77 0.1049841270 -0.0075158730
78 0.9848591270 0.1049841270
79 -0.1782658730 0.9848591270
80 0.0474841270 -0.1782658730
81 -0.1882658730 0.0474841270
82 -0.2313908730 -0.1882658730
83 -0.4870158730 -0.2313908730
84 -0.6750972222 -0.4870158730
85 0.3910277778 -0.6750972222
86 -0.0875972222 0.3910277778
87 0.1012777778 -0.0875972222
88 -0.2177222222 0.1012777778
89 -0.1442222222 -0.2177222222
90 -0.2403472222 -0.1442222222
91 -0.2354722222 -0.2403472222
92 0.0692777778 -0.2354722222
93 0.1675277778 0.0692777778
94 0.6244027778 0.1675277778
95 -0.5242222222 0.6244027778
96 NA -0.5242222222
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 0.0914722222 1.4913472222
[2,] -0.6801527778 0.0914722222
[3,] -0.3132777778 -0.6801527778
[4,] -0.2782777778 -0.3132777778
[5,] -0.0617777778 -0.2782777778
[6,] -0.1949027778 -0.0617777778
[7,] 0.0619722222 -0.1949027778
[8,] -0.2042777778 0.0619722222
[9,] -0.1410277778 -0.2042777778
[10,] 0.0828472222 -0.1410277778
[11,] -0.6467777778 0.0828472222
[12,] -0.3638591270 -0.6467777778
[13,] -0.5637341270 -0.3638591270
[14,] -0.3033591270 -0.5637341270
[15,] -0.1164841270 -0.3033591270
[16,] 0.3975158730 -0.1164841270
[17,] 0.0930158730 0.3975158730
[18,] -0.0791091270 0.0930158730
[19,] 0.1377658730 -0.0791091270
[20,] -0.1504841270 0.1377658730
[21,] -0.4832341270 -0.1504841270
[22,] -0.2243591270 -0.4832341270
[23,] 0.0170158730 -0.2243591270
[24,] 0.6849345238 0.0170158730
[25,] -0.0139404762 0.6849345238
[26,] -0.1935654762 -0.0139404762
[27,] 0.1093095238 -0.1935654762
[28,] 0.0353095238 0.1093095238
[29,] -0.0791904762 0.0353095238
[30,] 0.0816845238 -0.0791904762
[31,] 0.2015595238 0.0816845238
[32,] -0.0256904762 0.2015595238
[33,] 0.2105595238 -0.0256904762
[34,] -0.1295654762 0.2105595238
[35,] 0.2838095238 -0.1295654762
[36,] -0.4532718254 0.2838095238
[37,] -0.0101468254 -0.4532718254
[38,] 0.1362281746 -0.0101468254
[39,] 0.5341031746 0.1362281746
[40,] 0.0991031746 0.5341031746
[41,] 0.1236031746 0.0991031746
[42,] 0.1054781746 0.1236031746
[43,] 0.5863531746 0.1054781746
[44,] 0.4361031746 0.5863531746
[45,] 0.2943531746 0.4361031746
[46,] 0.2342281746 0.2943531746
[47,] 1.0856031746 0.2342281746
[48,] 0.0435218254 1.0856031746
[49,] -0.2883531746 0.0435218254
[50,] -0.3729781746 -0.2883531746
[51,] -0.1981031746 -0.3729781746
[52,] -0.2211031746 -0.1981031746
[53,] -0.1336031746 -0.2211031746
[54,] -0.1077281746 -0.1336031746
[55,] -0.1498531746 -0.1077281746
[56,] 0.0008968254 -0.1498531746
[57,] 0.1021468254 0.0008968254
[58,] -0.1309781746 0.1021468254
[59,] 0.1443968254 -0.1309781746
[60,] -0.2216845238 0.1443968254
[61,] 0.5874404762 -0.2216845238
[62,] 1.1908154762 0.5874404762
[63,] -0.0093095238 1.1908154762
[64,] 0.1926904762 -0.0093095238
[65,] 0.0971904762 0.1926904762
[66,] -0.5499345238 0.0971904762
[67,] -0.4240595238 -0.5499345238
[68,] -0.1733095238 -0.4240595238
[69,] 0.0379404762 -0.1733095238
[70,] -0.2251845238 0.0379404762
[71,] 0.1271904762 -0.2251845238
[72,] -0.5058908730 0.1271904762
[73,] -0.1937658730 -0.5058908730
[74,] 0.3106091270 -0.1937658730
[75,] -0.1075158730 0.3106091270
[76,] -0.0075158730 -0.1075158730
[77,] 0.1049841270 -0.0075158730
[78,] 0.9848591270 0.1049841270
[79,] -0.1782658730 0.9848591270
[80,] 0.0474841270 -0.1782658730
[81,] -0.1882658730 0.0474841270
[82,] -0.2313908730 -0.1882658730
[83,] -0.4870158730 -0.2313908730
[84,] -0.6750972222 -0.4870158730
[85,] 0.3910277778 -0.6750972222
[86,] -0.0875972222 0.3910277778
[87,] 0.1012777778 -0.0875972222
[88,] -0.2177222222 0.1012777778
[89,] -0.1442222222 -0.2177222222
[90,] -0.2403472222 -0.1442222222
[91,] -0.2354722222 -0.2403472222
[92,] 0.0692777778 -0.2354722222
[93,] 0.1675277778 0.0692777778
[94,] 0.6244027778 0.1675277778
[95,] -0.5242222222 0.6244027778
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 0.0914722222 1.4913472222
2 -0.6801527778 0.0914722222
3 -0.3132777778 -0.6801527778
4 -0.2782777778 -0.3132777778
5 -0.0617777778 -0.2782777778
6 -0.1949027778 -0.0617777778
7 0.0619722222 -0.1949027778
8 -0.2042777778 0.0619722222
9 -0.1410277778 -0.2042777778
10 0.0828472222 -0.1410277778
11 -0.6467777778 0.0828472222
12 -0.3638591270 -0.6467777778
13 -0.5637341270 -0.3638591270
14 -0.3033591270 -0.5637341270
15 -0.1164841270 -0.3033591270
16 0.3975158730 -0.1164841270
17 0.0930158730 0.3975158730
18 -0.0791091270 0.0930158730
19 0.1377658730 -0.0791091270
20 -0.1504841270 0.1377658730
21 -0.4832341270 -0.1504841270
22 -0.2243591270 -0.4832341270
23 0.0170158730 -0.2243591270
24 0.6849345238 0.0170158730
25 -0.0139404762 0.6849345238
26 -0.1935654762 -0.0139404762
27 0.1093095238 -0.1935654762
28 0.0353095238 0.1093095238
29 -0.0791904762 0.0353095238
30 0.0816845238 -0.0791904762
31 0.2015595238 0.0816845238
32 -0.0256904762 0.2015595238
33 0.2105595238 -0.0256904762
34 -0.1295654762 0.2105595238
35 0.2838095238 -0.1295654762
36 -0.4532718254 0.2838095238
37 -0.0101468254 -0.4532718254
38 0.1362281746 -0.0101468254
39 0.5341031746 0.1362281746
40 0.0991031746 0.5341031746
41 0.1236031746 0.0991031746
42 0.1054781746 0.1236031746
43 0.5863531746 0.1054781746
44 0.4361031746 0.5863531746
45 0.2943531746 0.4361031746
46 0.2342281746 0.2943531746
47 1.0856031746 0.2342281746
48 0.0435218254 1.0856031746
49 -0.2883531746 0.0435218254
50 -0.3729781746 -0.2883531746
51 -0.1981031746 -0.3729781746
52 -0.2211031746 -0.1981031746
53 -0.1336031746 -0.2211031746
54 -0.1077281746 -0.1336031746
55 -0.1498531746 -0.1077281746
56 0.0008968254 -0.1498531746
57 0.1021468254 0.0008968254
58 -0.1309781746 0.1021468254
59 0.1443968254 -0.1309781746
60 -0.2216845238 0.1443968254
61 0.5874404762 -0.2216845238
62 1.1908154762 0.5874404762
63 -0.0093095238 1.1908154762
64 0.1926904762 -0.0093095238
65 0.0971904762 0.1926904762
66 -0.5499345238 0.0971904762
67 -0.4240595238 -0.5499345238
68 -0.1733095238 -0.4240595238
69 0.0379404762 -0.1733095238
70 -0.2251845238 0.0379404762
71 0.1271904762 -0.2251845238
72 -0.5058908730 0.1271904762
73 -0.1937658730 -0.5058908730
74 0.3106091270 -0.1937658730
75 -0.1075158730 0.3106091270
76 -0.0075158730 -0.1075158730
77 0.1049841270 -0.0075158730
78 0.9848591270 0.1049841270
79 -0.1782658730 0.9848591270
80 0.0474841270 -0.1782658730
81 -0.1882658730 0.0474841270
82 -0.2313908730 -0.1882658730
83 -0.4870158730 -0.2313908730
84 -0.6750972222 -0.4870158730
85 0.3910277778 -0.6750972222
86 -0.0875972222 0.3910277778
87 0.1012777778 -0.0875972222
88 -0.2177222222 0.1012777778
89 -0.1442222222 -0.2177222222
90 -0.2403472222 -0.1442222222
91 -0.2354722222 -0.2403472222
92 0.0692777778 -0.2354722222
93 0.1675277778 0.0692777778
94 0.6244027778 0.1675277778
95 -0.5242222222 0.6244027778
> plot(z,main=paste('Residual Lag plot, lowess, and regression line'), ylab='values of Residuals', xlab='lagged values of Residuals')
> lines(lowess(z))
> abline(lm(z))
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/71en61354050200.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/815ck1354050200.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/9uwqt1354050200.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
> plot(mylm, las = 1, sub='Residual Diagnostics')
> par(opar)
> dev.off()
null device
1
> if (n > n25) {
+ postscript(file="/var/wessaorg/rcomp/tmp/10y7n41354050200.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
+ plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint')
+ grid()
+ dev.off()
+ }
null device
1
>
> #Note: the /var/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/wessaorg/rcomp/createtable")
>
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Estimated Regression Equation', 1, TRUE)
> a<-table.row.end(a)
> myeq <- colnames(x)[1]
> myeq <- paste(myeq, '[t] = ', sep='')
> for (i in 1:k){
+ if (mysum$coefficients[i,1] > 0) myeq <- paste(myeq, '+', '')
+ myeq <- paste(myeq, mysum$coefficients[i,1], sep=' ')
+ if (rownames(mysum$coefficients)[i] != '(Intercept)') {
+ myeq <- paste(myeq, rownames(mysum$coefficients)[i], sep='')
+ if (rownames(mysum$coefficients)[i] != 't') myeq <- paste(myeq, '[t]', sep='')
+ }
+ }
> myeq <- paste(myeq, ' + e[t]')
> a<-table.row.start(a)
> a<-table.element(a, myeq)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/11p0eg1354050200.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,hyperlink('http://www.xycoon.com/ols1.htm','Multiple Linear Regression - Ordinary Least Squares',''), 6, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Variable',header=TRUE)
> a<-table.element(a,'Parameter',header=TRUE)
> a<-table.element(a,'S.D.',header=TRUE)
> a<-table.element(a,'T-STAT
H0: parameter = 0',header=TRUE)
> a<-table.element(a,'2-tail p-value',header=TRUE)
> a<-table.element(a,'1-tail p-value',header=TRUE)
> a<-table.row.end(a)
> for (i in 1:k){
+ a<-table.row.start(a)
+ a<-table.element(a,rownames(mysum$coefficients)[i],header=TRUE)
+ a<-table.element(a,mysum$coefficients[i,1])
+ a<-table.element(a, round(mysum$coefficients[i,2],6))
+ a<-table.element(a, round(mysum$coefficients[i,3],4))
+ a<-table.element(a, round(mysum$coefficients[i,4],6))
+ a<-table.element(a, round(mysum$coefficients[i,4]/2,6))
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/12ws8d1354050200.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Regression Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple R',1,TRUE)
> a<-table.element(a, sqrt(mysum$r.squared))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'R-squared',1,TRUE)
> a<-table.element(a, mysum$r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Adjusted R-squared',1,TRUE)
> a<-table.element(a, mysum$adj.r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (value)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[3])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'p-value',1,TRUE)
> a<-table.element(a, 1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3]))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Residual Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Residual Standard Deviation',1,TRUE)
> a<-table.element(a, mysum$sigma)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Sum Squared Residuals',1,TRUE)
> a<-table.element(a, sum(myerror*myerror))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/13ciit1354050200.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Actuals, Interpolation, and Residuals', 4, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Time or Index', 1, TRUE)
> a<-table.element(a, 'Actuals', 1, TRUE)
> a<-table.element(a, 'Interpolation
Forecast', 1, TRUE)
> a<-table.element(a, 'Residuals
Prediction Error', 1, TRUE)
> a<-table.row.end(a)
> for (i in 1:n) {
+ a<-table.row.start(a)
+ a<-table.element(a,i, 1, TRUE)
+ a<-table.element(a,x[i])
+ a<-table.element(a,x[i]-mysum$resid[i])
+ a<-table.element(a,mysum$resid[i])
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/149wem1354050200.tab")
> if (n > n25) {
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'p-values',header=TRUE)
+ a<-table.element(a,'Alternative Hypothesis',3,header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'breakpoint index',header=TRUE)
+ a<-table.element(a,'greater',header=TRUE)
+ a<-table.element(a,'2-sided',header=TRUE)
+ a<-table.element(a,'less',header=TRUE)
+ a<-table.row.end(a)
+ for (mypoint in kp3:nmkm3) {
+ a<-table.row.start(a)
+ a<-table.element(a,mypoint,header=TRUE)
+ a<-table.element(a,gqarr[mypoint-kp3+1,1])
+ a<-table.element(a,gqarr[mypoint-kp3+1,2])
+ a<-table.element(a,gqarr[mypoint-kp3+1,3])
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/1516f21354050200.tab")
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Meta Analysis of Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'Description',header=TRUE)
+ a<-table.element(a,'# significant tests',header=TRUE)
+ a<-table.element(a,'% significant tests',header=TRUE)
+ a<-table.element(a,'OK/NOK',header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'1% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant1)
+ a<-table.element(a,numsignificant1/numgqtests)
+ if (numsignificant1/numgqtests < 0.01) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'5% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant5)
+ a<-table.element(a,numsignificant5/numgqtests)
+ if (numsignificant5/numgqtests < 0.05) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'10% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant10)
+ a<-table.element(a,numsignificant10/numgqtests)
+ if (numsignificant10/numgqtests < 0.1) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/168qtv1354050200.tab")
+ }
>
> try(system("convert tmp/13tsy1354050199.ps tmp/13tsy1354050199.png",intern=TRUE))
character(0)
> try(system("convert tmp/27sek1354050199.ps tmp/27sek1354050199.png",intern=TRUE))
character(0)
> try(system("convert tmp/3atct1354050199.ps tmp/3atct1354050199.png",intern=TRUE))
character(0)
> try(system("convert tmp/4so4l1354050199.ps tmp/4so4l1354050199.png",intern=TRUE))
character(0)
> try(system("convert tmp/5lnnu1354050199.ps tmp/5lnnu1354050199.png",intern=TRUE))
character(0)
> try(system("convert tmp/6zrzc1354050199.ps tmp/6zrzc1354050199.png",intern=TRUE))
character(0)
> try(system("convert tmp/71en61354050200.ps tmp/71en61354050200.png",intern=TRUE))
character(0)
> try(system("convert tmp/815ck1354050200.ps tmp/815ck1354050200.png",intern=TRUE))
character(0)
> try(system("convert tmp/9uwqt1354050200.ps tmp/9uwqt1354050200.png",intern=TRUE))
character(0)
> try(system("convert tmp/10y7n41354050200.ps tmp/10y7n41354050200.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
8.103 1.218 9.331