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(127,2.75,123,2.75,118,2.55,114,2.5,108,2.5,111,2.1,151,2,159,2,158,2,148,2,138,2,137,2,136,2,133,2,126,2,120,2,114,2,116,2,153,2,162,2,161,2,149,2,139,2,135,2,130,2,127,2,122,2,117,2,112,2,113,2,149,2,157,2,157,2,147,2,137,2,132,2.21,125,2.25,123,2.25,117,2.45,114,2.5,111,2.5,112,2.64,144,2.75,150,2.93,149,3,134,3.17,123,3.25,116,3.39,117,3.5,111,3.5,105,3.65,102,3.75,95,3.75,93,3.9,124,4,130,4,124,4,115,4,106,4,105,4,105,4,101,4,95,4,93,4,84,4,87,4,116,4.18,120,4.25,117,4.25,109,3.97,105,3.42,107,2.75),dim=c(2,72),dimnames=list(c('Werkloosheid','Rente'),1:72))
> y <- array(NA,dim=c(2,72),dimnames=list(c('Werkloosheid','Rente'),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
Werkloosheid Rente
1 127 2.75
2 123 2.75
3 118 2.55
4 114 2.50
5 108 2.50
6 111 2.10
7 151 2.00
8 159 2.00
9 158 2.00
10 148 2.00
11 138 2.00
12 137 2.00
13 136 2.00
14 133 2.00
15 126 2.00
16 120 2.00
17 114 2.00
18 116 2.00
19 153 2.00
20 162 2.00
21 161 2.00
22 149 2.00
23 139 2.00
24 135 2.00
25 130 2.00
26 127 2.00
27 122 2.00
28 117 2.00
29 112 2.00
30 113 2.00
31 149 2.00
32 157 2.00
33 157 2.00
34 147 2.00
35 137 2.00
36 132 2.21
37 125 2.25
38 123 2.25
39 117 2.45
40 114 2.50
41 111 2.50
42 112 2.64
43 144 2.75
44 150 2.93
45 149 3.00
46 134 3.17
47 123 3.25
48 116 3.39
49 117 3.50
50 111 3.50
51 105 3.65
52 102 3.75
53 95 3.75
54 93 3.90
55 124 4.00
56 130 4.00
57 124 4.00
58 115 4.00
59 106 4.00
60 105 4.00
61 105 4.00
62 101 4.00
63 95 4.00
64 93 4.00
65 84 4.00
66 87 4.00
67 116 4.18
68 120 4.25
69 117 4.25
70 109 3.97
71 105 3.42
72 107 2.75
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Rente
163.54 -14.18
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-23.176 -13.336 -1.178 13.753 28.012
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 163.536 6.162 26.54 < 2e-16 ***
Rente -14.180 2.110 -6.72 4.02e-09 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 15.1 on 70 degrees of freedom
Multiple R-squared: 0.3922, Adjusted R-squared: 0.3835
F-statistic: 45.16 on 1 and 70 DF, p-value: 4.023e-09
> 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.01745967 0.03491933 0.982540333
[2,] 0.02732049 0.05464098 0.972679509
[3,] 0.55655530 0.88688940 0.443444700
[4,] 0.70628264 0.58743472 0.293717362
[5,] 0.71324575 0.57350849 0.286754246
[6,] 0.62537203 0.74925593 0.374627966
[7,] 0.53835941 0.92328118 0.461640589
[8,] 0.45253295 0.90506591 0.547467046
[9,] 0.37210848 0.74421697 0.627891516
[10,] 0.30980018 0.61960036 0.690199819
[11,] 0.29876009 0.59752019 0.701239907
[12,] 0.33889385 0.67778769 0.661106154
[13,] 0.44188702 0.88377403 0.558112984
[14,] 0.49192245 0.98384491 0.508077546
[15,] 0.52007041 0.95985918 0.479929591
[16,] 0.65749891 0.68500219 0.342501094
[17,] 0.75147017 0.49705966 0.248529831
[18,] 0.72809576 0.54380847 0.271904237
[19,] 0.66663409 0.66673183 0.333365914
[20,] 0.60111973 0.79776055 0.398880275
[21,] 0.54582381 0.90835239 0.454176193
[22,] 0.50326608 0.99346784 0.496733919
[23,] 0.49462208 0.98924416 0.505377918
[24,] 0.53288479 0.93423042 0.467115209
[25,] 0.63082100 0.73835799 0.369178996
[26,] 0.71263705 0.57472590 0.287362952
[27,] 0.69332320 0.61335361 0.306676803
[28,] 0.74337294 0.51325411 0.256627057
[29,] 0.79647956 0.40704088 0.203520442
[30,] 0.78238998 0.43522004 0.217610019
[31,] 0.73349057 0.53301885 0.266509426
[32,] 0.67674577 0.64650846 0.323254231
[33,] 0.61548962 0.76902076 0.384510381
[34,] 0.55630583 0.88738834 0.443694168
[35,] 0.50564345 0.98871310 0.494356549
[36,] 0.47242363 0.94484726 0.527576369
[37,] 0.48061637 0.96123273 0.519383634
[38,] 0.48913269 0.97826539 0.510867305
[39,] 0.56896854 0.86206293 0.431031464
[40,] 0.75181353 0.49637294 0.248186470
[41,] 0.89923891 0.20152218 0.100761088
[42,] 0.92952053 0.14095893 0.070479467
[43,] 0.92661778 0.14676444 0.073382219
[44,] 0.90969028 0.18061944 0.090309718
[45,] 0.89506864 0.20986272 0.104931361
[46,] 0.86543074 0.26913851 0.134569256
[47,] 0.82282122 0.35435757 0.177178784
[48,] 0.77430164 0.45139671 0.225698356
[49,] 0.75319396 0.49361207 0.246806036
[50,] 0.74957928 0.50084145 0.250420723
[51,] 0.77431373 0.45137255 0.225686273
[52,] 0.87399467 0.25201067 0.126005334
[53,] 0.91030964 0.17938072 0.089690362
[54,] 0.89387912 0.21224177 0.106120884
[55,] 0.84277417 0.31445165 0.157225827
[56,] 0.77444722 0.45110555 0.225552777
[57,] 0.68882947 0.62234105 0.311170526
[58,] 0.59142861 0.81714278 0.408571391
[59,] 0.52471706 0.95056588 0.475282941
[60,] 0.48706617 0.97413234 0.512933828
[61,] 0.70856833 0.58286333 0.291431667
[62,] 0.99003600 0.01992801 0.009964004
[63,] 0.96134091 0.07731819 0.038659095
> postscript(file="/var/www/html/rcomp/tmp/14irx1258710612.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/2qxmx1258710612.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/33s001258710612.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/47me81258710612.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/5tqab1258710612.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> qqnorm(mysum$resid, main='Residual Normal Q-Q Plot')
> qqline(mysum$resid)
> grid()
> dev.off()
null device
1
> (myerror <- as.ts(mysum$resid))
Time Series:
Start = 1
End = 72
Frequency = 1
1 2 3 4 5 6
2.4592172 -1.5407828 -9.3768168 -14.0858253 -20.0858253 -22.7578932
7 8 9 10 11 12
15.8240898 23.8240898 22.8240898 12.8240898 2.8240898 1.8240898
13 14 15 16 17 18
0.8240898 -2.1759102 -9.1759102 -15.1759102 -21.1759102 -19.1759102
19 20 21 22 23 24
17.8240898 26.8240898 25.8240898 13.8240898 3.8240898 -0.1759102
25 26 27 28 29 30
-5.1759102 -8.1759102 -13.1759102 -18.1759102 -23.1759102 -22.1759102
31 32 33 34 35 36
13.8240898 21.8240898 21.8240898 11.8240898 1.8240898 -0.1980745
37 38 39 40 41 42
-6.6308677 -8.6308677 -11.7948337 -14.0858253 -17.0858253 -14.1006015
43 44 45 46 47 48
19.4592172 28.0116478 28.0042596 15.4148885 5.5493021 0.5345259
49 50 51 52 53 54
3.0943445 -2.9056555 -6.7786300 -8.3606130 -15.3606130 -15.2335875
55 56 57 58 59 60
17.1844294 23.1844294 17.1844294 8.1844294 -0.8155706 -1.8155706
61 62 63 64 65 66
-1.8155706 -5.8155706 -11.8155706 -13.8155706 -22.8155706 -19.8155706
67 68 69 70 71 72
11.7368600 16.7294719 13.7294719 1.7590243 -10.0400690 -17.5407828
> postscript(file="/var/www/html/rcomp/tmp/6vd731258710612.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> dum <- cbind(lag(myerror,k=1),myerror)
> dum
Time Series:
Start = 0
End = 72
Frequency = 1
lag(myerror, k = 1) myerror
0 2.4592172 NA
1 -1.5407828 2.4592172
2 -9.3768168 -1.5407828
3 -14.0858253 -9.3768168
4 -20.0858253 -14.0858253
5 -22.7578932 -20.0858253
6 15.8240898 -22.7578932
7 23.8240898 15.8240898
8 22.8240898 23.8240898
9 12.8240898 22.8240898
10 2.8240898 12.8240898
11 1.8240898 2.8240898
12 0.8240898 1.8240898
13 -2.1759102 0.8240898
14 -9.1759102 -2.1759102
15 -15.1759102 -9.1759102
16 -21.1759102 -15.1759102
17 -19.1759102 -21.1759102
18 17.8240898 -19.1759102
19 26.8240898 17.8240898
20 25.8240898 26.8240898
21 13.8240898 25.8240898
22 3.8240898 13.8240898
23 -0.1759102 3.8240898
24 -5.1759102 -0.1759102
25 -8.1759102 -5.1759102
26 -13.1759102 -8.1759102
27 -18.1759102 -13.1759102
28 -23.1759102 -18.1759102
29 -22.1759102 -23.1759102
30 13.8240898 -22.1759102
31 21.8240898 13.8240898
32 21.8240898 21.8240898
33 11.8240898 21.8240898
34 1.8240898 11.8240898
35 -0.1980745 1.8240898
36 -6.6308677 -0.1980745
37 -8.6308677 -6.6308677
38 -11.7948337 -8.6308677
39 -14.0858253 -11.7948337
40 -17.0858253 -14.0858253
41 -14.1006015 -17.0858253
42 19.4592172 -14.1006015
43 28.0116478 19.4592172
44 28.0042596 28.0116478
45 15.4148885 28.0042596
46 5.5493021 15.4148885
47 0.5345259 5.5493021
48 3.0943445 0.5345259
49 -2.9056555 3.0943445
50 -6.7786300 -2.9056555
51 -8.3606130 -6.7786300
52 -15.3606130 -8.3606130
53 -15.2335875 -15.3606130
54 17.1844294 -15.2335875
55 23.1844294 17.1844294
56 17.1844294 23.1844294
57 8.1844294 17.1844294
58 -0.8155706 8.1844294
59 -1.8155706 -0.8155706
60 -1.8155706 -1.8155706
61 -5.8155706 -1.8155706
62 -11.8155706 -5.8155706
63 -13.8155706 -11.8155706
64 -22.8155706 -13.8155706
65 -19.8155706 -22.8155706
66 11.7368600 -19.8155706
67 16.7294719 11.7368600
68 13.7294719 16.7294719
69 1.7590243 13.7294719
70 -10.0400690 1.7590243
71 -17.5407828 -10.0400690
72 NA -17.5407828
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -1.5407828 2.4592172
[2,] -9.3768168 -1.5407828
[3,] -14.0858253 -9.3768168
[4,] -20.0858253 -14.0858253
[5,] -22.7578932 -20.0858253
[6,] 15.8240898 -22.7578932
[7,] 23.8240898 15.8240898
[8,] 22.8240898 23.8240898
[9,] 12.8240898 22.8240898
[10,] 2.8240898 12.8240898
[11,] 1.8240898 2.8240898
[12,] 0.8240898 1.8240898
[13,] -2.1759102 0.8240898
[14,] -9.1759102 -2.1759102
[15,] -15.1759102 -9.1759102
[16,] -21.1759102 -15.1759102
[17,] -19.1759102 -21.1759102
[18,] 17.8240898 -19.1759102
[19,] 26.8240898 17.8240898
[20,] 25.8240898 26.8240898
[21,] 13.8240898 25.8240898
[22,] 3.8240898 13.8240898
[23,] -0.1759102 3.8240898
[24,] -5.1759102 -0.1759102
[25,] -8.1759102 -5.1759102
[26,] -13.1759102 -8.1759102
[27,] -18.1759102 -13.1759102
[28,] -23.1759102 -18.1759102
[29,] -22.1759102 -23.1759102
[30,] 13.8240898 -22.1759102
[31,] 21.8240898 13.8240898
[32,] 21.8240898 21.8240898
[33,] 11.8240898 21.8240898
[34,] 1.8240898 11.8240898
[35,] -0.1980745 1.8240898
[36,] -6.6308677 -0.1980745
[37,] -8.6308677 -6.6308677
[38,] -11.7948337 -8.6308677
[39,] -14.0858253 -11.7948337
[40,] -17.0858253 -14.0858253
[41,] -14.1006015 -17.0858253
[42,] 19.4592172 -14.1006015
[43,] 28.0116478 19.4592172
[44,] 28.0042596 28.0116478
[45,] 15.4148885 28.0042596
[46,] 5.5493021 15.4148885
[47,] 0.5345259 5.5493021
[48,] 3.0943445 0.5345259
[49,] -2.9056555 3.0943445
[50,] -6.7786300 -2.9056555
[51,] -8.3606130 -6.7786300
[52,] -15.3606130 -8.3606130
[53,] -15.2335875 -15.3606130
[54,] 17.1844294 -15.2335875
[55,] 23.1844294 17.1844294
[56,] 17.1844294 23.1844294
[57,] 8.1844294 17.1844294
[58,] -0.8155706 8.1844294
[59,] -1.8155706 -0.8155706
[60,] -1.8155706 -1.8155706
[61,] -5.8155706 -1.8155706
[62,] -11.8155706 -5.8155706
[63,] -13.8155706 -11.8155706
[64,] -22.8155706 -13.8155706
[65,] -19.8155706 -22.8155706
[66,] 11.7368600 -19.8155706
[67,] 16.7294719 11.7368600
[68,] 13.7294719 16.7294719
[69,] 1.7590243 13.7294719
[70,] -10.0400690 1.7590243
[71,] -17.5407828 -10.0400690
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -1.5407828 2.4592172
2 -9.3768168 -1.5407828
3 -14.0858253 -9.3768168
4 -20.0858253 -14.0858253
5 -22.7578932 -20.0858253
6 15.8240898 -22.7578932
7 23.8240898 15.8240898
8 22.8240898 23.8240898
9 12.8240898 22.8240898
10 2.8240898 12.8240898
11 1.8240898 2.8240898
12 0.8240898 1.8240898
13 -2.1759102 0.8240898
14 -9.1759102 -2.1759102
15 -15.1759102 -9.1759102
16 -21.1759102 -15.1759102
17 -19.1759102 -21.1759102
18 17.8240898 -19.1759102
19 26.8240898 17.8240898
20 25.8240898 26.8240898
21 13.8240898 25.8240898
22 3.8240898 13.8240898
23 -0.1759102 3.8240898
24 -5.1759102 -0.1759102
25 -8.1759102 -5.1759102
26 -13.1759102 -8.1759102
27 -18.1759102 -13.1759102
28 -23.1759102 -18.1759102
29 -22.1759102 -23.1759102
30 13.8240898 -22.1759102
31 21.8240898 13.8240898
32 21.8240898 21.8240898
33 11.8240898 21.8240898
34 1.8240898 11.8240898
35 -0.1980745 1.8240898
36 -6.6308677 -0.1980745
37 -8.6308677 -6.6308677
38 -11.7948337 -8.6308677
39 -14.0858253 -11.7948337
40 -17.0858253 -14.0858253
41 -14.1006015 -17.0858253
42 19.4592172 -14.1006015
43 28.0116478 19.4592172
44 28.0042596 28.0116478
45 15.4148885 28.0042596
46 5.5493021 15.4148885
47 0.5345259 5.5493021
48 3.0943445 0.5345259
49 -2.9056555 3.0943445
50 -6.7786300 -2.9056555
51 -8.3606130 -6.7786300
52 -15.3606130 -8.3606130
53 -15.2335875 -15.3606130
54 17.1844294 -15.2335875
55 23.1844294 17.1844294
56 17.1844294 23.1844294
57 8.1844294 17.1844294
58 -0.8155706 8.1844294
59 -1.8155706 -0.8155706
60 -1.8155706 -1.8155706
61 -5.8155706 -1.8155706
62 -11.8155706 -5.8155706
63 -13.8155706 -11.8155706
64 -22.8155706 -13.8155706
65 -19.8155706 -22.8155706
66 11.7368600 -19.8155706
67 16.7294719 11.7368600
68 13.7294719 16.7294719
69 1.7590243 13.7294719
70 -10.0400690 1.7590243
71 -17.5407828 -10.0400690
> 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/7r8171258710612.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/8yz121258710612.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/91yqs1258710612.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/10bpww1258710612.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/11pkln1258710612.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/12vw801258710612.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/13o9cr1258710612.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/14ji681258710612.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/15jt4k1258710612.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/16ajci1258710612.tab")
+ }
>
> system("convert tmp/14irx1258710612.ps tmp/14irx1258710612.png")
> system("convert tmp/2qxmx1258710612.ps tmp/2qxmx1258710612.png")
> system("convert tmp/33s001258710612.ps tmp/33s001258710612.png")
> system("convert tmp/47me81258710612.ps tmp/47me81258710612.png")
> system("convert tmp/5tqab1258710612.ps tmp/5tqab1258710612.png")
> system("convert tmp/6vd731258710612.ps tmp/6vd731258710612.png")
> system("convert tmp/7r8171258710612.ps tmp/7r8171258710612.png")
> system("convert tmp/8yz121258710612.ps tmp/8yz121258710612.png")
> system("convert tmp/91yqs1258710612.ps tmp/91yqs1258710612.png")
> system("convert tmp/10bpww1258710612.ps tmp/10bpww1258710612.png")
>
>
> proc.time()
user system elapsed
2.613 1.587 3.475