R version 2.8.0 (2008-10-20)
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.
Natural language support but running in an English locale
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(376.974,0,377.632,0,378.205,0,370.861,0,369.167,0,371.551,0,382.842,0,381.903,0,384.502,0,392.058,0,384.359,0,388.884,0,386.586,0,387.495,0,385.705,0,378.67,0,377.367,0,376.911,0,389.827,0,387.82,0,387.267,0,380.575,0,372.402,0,376.74,0,377.795,0,376.126,0,370.804,0,367.98,0,367.866,0,366.121,0,379.421,0,378.519,0,372.423,0,355.072,0,344.693,0,342.892,0,344.178,0,337.606,0,327.103,0,323.953,0,316.532,0,306.307,0,327.225,0,329.573,0,313.761,0,307.836,0,300.074,0,304.198,0,306.122,0,300.414,0,292.133,0,290.616,0,280.244,1,285.179,1,305.486,1,305.957,1,293.886,1,289.441,1,288.776,1,299.149,1,306.532,1,309.914,1,313.468,1,314.901,1,309.16,1,316.15,1,336.544,1,339.196,1,326.738,1,320.838,1,318.62,1,331.533,1,335.378,1),dim=c(2,73),dimnames=list(c('Maandelijksewerkloosheid','x'),1:73))
> y <- array(NA,dim=c(2,73),dimnames=list(c('Maandelijksewerkloosheid','x'),1:73))
> 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
Maandelijksewerkloosheid x
1 376.974 0
2 377.632 0
3 378.205 0
4 370.861 0
5 369.167 0
6 371.551 0
7 382.842 0
8 381.903 0
9 384.502 0
10 392.058 0
11 384.359 0
12 388.884 0
13 386.586 0
14 387.495 0
15 385.705 0
16 378.670 0
17 377.367 0
18 376.911 0
19 389.827 0
20 387.820 0
21 387.267 0
22 380.575 0
23 372.402 0
24 376.740 0
25 377.795 0
26 376.126 0
27 370.804 0
28 367.980 0
29 367.866 0
30 366.121 0
31 379.421 0
32 378.519 0
33 372.423 0
34 355.072 0
35 344.693 0
36 342.892 0
37 344.178 0
38 337.606 0
39 327.103 0
40 323.953 0
41 316.532 0
42 306.307 0
43 327.225 0
44 329.573 0
45 313.761 0
46 307.836 0
47 300.074 0
48 304.198 0
49 306.122 0
50 300.414 0
51 292.133 0
52 290.616 0
53 280.244 1
54 285.179 1
55 305.486 1
56 305.957 1
57 293.886 1
58 289.441 1
59 288.776 1
60 299.149 1
61 306.532 1
62 309.914 1
63 313.468 1
64 314.901 1
65 309.160 1
66 316.150 1
67 336.544 1
68 339.196 1
69 326.738 1
70 320.838 1
71 318.620 1
72 331.533 1
73 335.378 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) x
357.19 -46.37
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-66.57 -19.58 10.68 21.33 34.87
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 357.186 3.924 91.018 < 2e-16 ***
x -46.372 7.317 -6.338 1.88e-08 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 28.3 on 71 degrees of freedom
Multiple R-squared: 0.3613, Adjusted R-squared: 0.3523
F-statistic: 40.17 on 1 and 71 DF, p-value: 1.877e-08
> 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,] 6.535779e-03 1.307156e-02 0.9934642214
[2,] 1.089992e-03 2.179984e-03 0.9989100082
[3,] 6.434099e-04 1.286820e-03 0.9993565901
[4,] 2.115500e-04 4.230999e-04 0.9997884500
[5,] 1.014368e-04 2.028736e-04 0.9998985632
[6,] 2.063869e-04 4.127738e-04 0.9997936131
[7,] 6.974292e-05 1.394858e-04 0.9999302571
[8,] 4.076710e-05 8.153420e-05 0.9999592329
[9,] 1.629129e-05 3.258258e-05 0.9999837087
[10,] 7.094749e-06 1.418950e-05 0.9999929053
[11,] 2.529948e-06 5.059896e-06 0.9999974701
[12,] 7.464055e-07 1.492811e-06 0.9999992536
[13,] 2.331801e-07 4.663602e-07 0.9999997668
[14,] 7.474405e-08 1.494881e-07 0.9999999253
[15,] 5.936767e-08 1.187353e-07 0.9999999406
[16,] 3.377807e-08 6.755614e-08 0.9999999662
[17,] 1.884536e-08 3.769071e-08 0.9999999812
[18,] 7.235157e-09 1.447031e-08 0.9999999928
[19,] 5.852693e-09 1.170539e-08 0.9999999941
[20,] 2.928030e-09 5.856059e-09 0.9999999971
[21,] 1.509682e-09 3.019365e-09 0.9999999985
[22,] 9.509685e-10 1.901937e-09 0.9999999990
[23,] 1.244361e-09 2.488723e-09 0.9999999988
[24,] 2.657251e-09 5.314501e-09 0.9999999973
[25,] 5.451304e-09 1.090261e-08 0.9999999945
[26,] 1.452498e-08 2.904995e-08 0.9999999855
[27,] 2.601123e-08 5.202247e-08 0.9999999740
[28,] 7.696402e-08 1.539280e-07 0.9999999230
[29,] 3.790141e-07 7.580281e-07 0.9999996210
[30,] 1.627365e-05 3.254730e-05 0.9999837264
[31,] 9.961002e-04 1.992200e-03 0.9990038998
[32,] 1.290087e-02 2.580173e-02 0.9870991329
[33,] 6.332037e-02 1.266407e-01 0.9366796329
[34,] 2.116224e-01 4.232448e-01 0.7883776128
[35,] 4.791511e-01 9.583023e-01 0.5208488605
[36,] 7.057105e-01 5.885791e-01 0.2942895310
[37,] 8.567121e-01 2.865759e-01 0.1432879403
[38,] 9.429256e-01 1.141488e-01 0.0570743978
[39,] 9.646471e-01 7.070573e-02 0.0353528625
[40,] 9.811019e-01 3.779620e-02 0.0188980996
[41,] 9.884632e-01 2.307363e-02 0.0115368164
[42,] 9.924901e-01 1.501979e-02 0.0075098928
[43,] 9.950292e-01 9.941513e-03 0.0049707567
[44,] 9.958611e-01 8.277823e-03 0.0041389114
[45,] 9.962566e-01 7.486865e-03 0.0037434325
[46,] 9.964722e-01 7.055535e-03 0.0035277677
[47,] 9.966409e-01 6.718115e-03 0.0033590576
[48,] 9.964970e-01 7.005984e-03 0.0035029921
[49,] 9.977789e-01 4.442136e-03 0.0022210681
[50,] 9.984078e-01 3.184474e-03 0.0015922368
[51,] 9.971420e-01 5.716024e-03 0.0028580119
[52,] 9.948675e-01 1.026501e-02 0.0051325058
[53,] 9.944549e-01 1.109015e-02 0.0055450739
[54,] 9.963780e-01 7.243938e-03 0.0036219690
[55,] 9.986908e-01 2.618351e-03 0.0013091755
[56,] 9.990562e-01 1.887677e-03 0.0009438384
[57,] 9.988224e-01 2.355109e-03 0.0011775547
[58,] 9.982331e-01 3.533890e-03 0.0017669450
[59,] 9.967507e-01 6.498575e-03 0.0032492874
[60,] 9.938407e-01 1.231865e-02 0.0061593265
[61,] 9.950436e-01 9.912868e-03 0.0049564340
[62,] 9.930491e-01 1.390171e-02 0.0069508565
[63,] 9.822026e-01 3.559472e-02 0.0177973603
[64,] 9.710239e-01 5.795227e-02 0.0289761337
> postscript(file="/var/www/html/freestat/rcomp/tmp/19fif1291051166.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/freestat/rcomp/tmp/2koii1291051166.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/freestat/rcomp/tmp/3koii1291051166.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/freestat/rcomp/tmp/4koii1291051166.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/freestat/rcomp/tmp/5koii1291051166.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 = 73
Frequency = 1
1 2 3 4 5 6
19.7885000 20.4465000 21.0195000 13.6755000 11.9815000 14.3655000
7 8 9 10 11 12
25.6565000 24.7175000 27.3165000 34.8725000 27.1735000 31.6985000
13 14 15 16 17 18
29.4005000 30.3095000 28.5195000 21.4845000 20.1815000 19.7255000
19 20 21 22 23 24
32.6415000 30.6345000 30.0815000 23.3895000 15.2165000 19.5545000
25 26 27 28 29 30
20.6095000 18.9405000 13.6185000 10.7945000 10.6805000 8.9355000
31 32 33 34 35 36
22.2355000 21.3335000 15.2375000 -2.1135000 -12.4925000 -14.2935000
37 38 39 40 41 42
-13.0075000 -19.5795000 -30.0825000 -33.2325000 -40.6535000 -50.8785000
43 44 45 46 47 48
-29.9605000 -27.6125000 -43.4245000 -49.3495000 -57.1115000 -52.9875000
49 50 51 52 53 54
-51.0635000 -56.7715000 -65.0525000 -66.5695000 -30.5698095 -25.6348095
55 56 57 58 59 60
-5.3278095 -4.8568095 -16.9278095 -21.3728095 -22.0378095 -11.6648095
61 62 63 64 65 66
-4.2818095 -0.8998095 2.6541905 4.0871905 -1.6538095 5.3361905
67 68 69 70 71 72
25.7301905 28.3821905 15.9241905 10.0241905 7.8061905 20.7191905
73
24.5641905
> postscript(file="/var/www/html/freestat/rcomp/tmp/6cfz31291051166.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 = 73
Frequency = 1
lag(myerror, k = 1) myerror
0 19.7885000 NA
1 20.4465000 19.7885000
2 21.0195000 20.4465000
3 13.6755000 21.0195000
4 11.9815000 13.6755000
5 14.3655000 11.9815000
6 25.6565000 14.3655000
7 24.7175000 25.6565000
8 27.3165000 24.7175000
9 34.8725000 27.3165000
10 27.1735000 34.8725000
11 31.6985000 27.1735000
12 29.4005000 31.6985000
13 30.3095000 29.4005000
14 28.5195000 30.3095000
15 21.4845000 28.5195000
16 20.1815000 21.4845000
17 19.7255000 20.1815000
18 32.6415000 19.7255000
19 30.6345000 32.6415000
20 30.0815000 30.6345000
21 23.3895000 30.0815000
22 15.2165000 23.3895000
23 19.5545000 15.2165000
24 20.6095000 19.5545000
25 18.9405000 20.6095000
26 13.6185000 18.9405000
27 10.7945000 13.6185000
28 10.6805000 10.7945000
29 8.9355000 10.6805000
30 22.2355000 8.9355000
31 21.3335000 22.2355000
32 15.2375000 21.3335000
33 -2.1135000 15.2375000
34 -12.4925000 -2.1135000
35 -14.2935000 -12.4925000
36 -13.0075000 -14.2935000
37 -19.5795000 -13.0075000
38 -30.0825000 -19.5795000
39 -33.2325000 -30.0825000
40 -40.6535000 -33.2325000
41 -50.8785000 -40.6535000
42 -29.9605000 -50.8785000
43 -27.6125000 -29.9605000
44 -43.4245000 -27.6125000
45 -49.3495000 -43.4245000
46 -57.1115000 -49.3495000
47 -52.9875000 -57.1115000
48 -51.0635000 -52.9875000
49 -56.7715000 -51.0635000
50 -65.0525000 -56.7715000
51 -66.5695000 -65.0525000
52 -30.5698095 -66.5695000
53 -25.6348095 -30.5698095
54 -5.3278095 -25.6348095
55 -4.8568095 -5.3278095
56 -16.9278095 -4.8568095
57 -21.3728095 -16.9278095
58 -22.0378095 -21.3728095
59 -11.6648095 -22.0378095
60 -4.2818095 -11.6648095
61 -0.8998095 -4.2818095
62 2.6541905 -0.8998095
63 4.0871905 2.6541905
64 -1.6538095 4.0871905
65 5.3361905 -1.6538095
66 25.7301905 5.3361905
67 28.3821905 25.7301905
68 15.9241905 28.3821905
69 10.0241905 15.9241905
70 7.8061905 10.0241905
71 20.7191905 7.8061905
72 24.5641905 20.7191905
73 NA 24.5641905
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 20.4465000 19.7885000
[2,] 21.0195000 20.4465000
[3,] 13.6755000 21.0195000
[4,] 11.9815000 13.6755000
[5,] 14.3655000 11.9815000
[6,] 25.6565000 14.3655000
[7,] 24.7175000 25.6565000
[8,] 27.3165000 24.7175000
[9,] 34.8725000 27.3165000
[10,] 27.1735000 34.8725000
[11,] 31.6985000 27.1735000
[12,] 29.4005000 31.6985000
[13,] 30.3095000 29.4005000
[14,] 28.5195000 30.3095000
[15,] 21.4845000 28.5195000
[16,] 20.1815000 21.4845000
[17,] 19.7255000 20.1815000
[18,] 32.6415000 19.7255000
[19,] 30.6345000 32.6415000
[20,] 30.0815000 30.6345000
[21,] 23.3895000 30.0815000
[22,] 15.2165000 23.3895000
[23,] 19.5545000 15.2165000
[24,] 20.6095000 19.5545000
[25,] 18.9405000 20.6095000
[26,] 13.6185000 18.9405000
[27,] 10.7945000 13.6185000
[28,] 10.6805000 10.7945000
[29,] 8.9355000 10.6805000
[30,] 22.2355000 8.9355000
[31,] 21.3335000 22.2355000
[32,] 15.2375000 21.3335000
[33,] -2.1135000 15.2375000
[34,] -12.4925000 -2.1135000
[35,] -14.2935000 -12.4925000
[36,] -13.0075000 -14.2935000
[37,] -19.5795000 -13.0075000
[38,] -30.0825000 -19.5795000
[39,] -33.2325000 -30.0825000
[40,] -40.6535000 -33.2325000
[41,] -50.8785000 -40.6535000
[42,] -29.9605000 -50.8785000
[43,] -27.6125000 -29.9605000
[44,] -43.4245000 -27.6125000
[45,] -49.3495000 -43.4245000
[46,] -57.1115000 -49.3495000
[47,] -52.9875000 -57.1115000
[48,] -51.0635000 -52.9875000
[49,] -56.7715000 -51.0635000
[50,] -65.0525000 -56.7715000
[51,] -66.5695000 -65.0525000
[52,] -30.5698095 -66.5695000
[53,] -25.6348095 -30.5698095
[54,] -5.3278095 -25.6348095
[55,] -4.8568095 -5.3278095
[56,] -16.9278095 -4.8568095
[57,] -21.3728095 -16.9278095
[58,] -22.0378095 -21.3728095
[59,] -11.6648095 -22.0378095
[60,] -4.2818095 -11.6648095
[61,] -0.8998095 -4.2818095
[62,] 2.6541905 -0.8998095
[63,] 4.0871905 2.6541905
[64,] -1.6538095 4.0871905
[65,] 5.3361905 -1.6538095
[66,] 25.7301905 5.3361905
[67,] 28.3821905 25.7301905
[68,] 15.9241905 28.3821905
[69,] 10.0241905 15.9241905
[70,] 7.8061905 10.0241905
[71,] 20.7191905 7.8061905
[72,] 24.5641905 20.7191905
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 20.4465000 19.7885000
2 21.0195000 20.4465000
3 13.6755000 21.0195000
4 11.9815000 13.6755000
5 14.3655000 11.9815000
6 25.6565000 14.3655000
7 24.7175000 25.6565000
8 27.3165000 24.7175000
9 34.8725000 27.3165000
10 27.1735000 34.8725000
11 31.6985000 27.1735000
12 29.4005000 31.6985000
13 30.3095000 29.4005000
14 28.5195000 30.3095000
15 21.4845000 28.5195000
16 20.1815000 21.4845000
17 19.7255000 20.1815000
18 32.6415000 19.7255000
19 30.6345000 32.6415000
20 30.0815000 30.6345000
21 23.3895000 30.0815000
22 15.2165000 23.3895000
23 19.5545000 15.2165000
24 20.6095000 19.5545000
25 18.9405000 20.6095000
26 13.6185000 18.9405000
27 10.7945000 13.6185000
28 10.6805000 10.7945000
29 8.9355000 10.6805000
30 22.2355000 8.9355000
31 21.3335000 22.2355000
32 15.2375000 21.3335000
33 -2.1135000 15.2375000
34 -12.4925000 -2.1135000
35 -14.2935000 -12.4925000
36 -13.0075000 -14.2935000
37 -19.5795000 -13.0075000
38 -30.0825000 -19.5795000
39 -33.2325000 -30.0825000
40 -40.6535000 -33.2325000
41 -50.8785000 -40.6535000
42 -29.9605000 -50.8785000
43 -27.6125000 -29.9605000
44 -43.4245000 -27.6125000
45 -49.3495000 -43.4245000
46 -57.1115000 -49.3495000
47 -52.9875000 -57.1115000
48 -51.0635000 -52.9875000
49 -56.7715000 -51.0635000
50 -65.0525000 -56.7715000
51 -66.5695000 -65.0525000
52 -30.5698095 -66.5695000
53 -25.6348095 -30.5698095
54 -5.3278095 -25.6348095
55 -4.8568095 -5.3278095
56 -16.9278095 -4.8568095
57 -21.3728095 -16.9278095
58 -22.0378095 -21.3728095
59 -11.6648095 -22.0378095
60 -4.2818095 -11.6648095
61 -0.8998095 -4.2818095
62 2.6541905 -0.8998095
63 4.0871905 2.6541905
64 -1.6538095 4.0871905
65 5.3361905 -1.6538095
66 25.7301905 5.3361905
67 28.3821905 25.7301905
68 15.9241905 28.3821905
69 10.0241905 15.9241905
70 7.8061905 10.0241905
71 20.7191905 7.8061905
72 24.5641905 20.7191905
> 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/freestat/rcomp/tmp/75pgo1291051166.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/freestat/rcomp/tmp/85pgo1291051166.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/freestat/rcomp/tmp/95pgo1291051166.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/freestat/rcomp/tmp/10yyfq1291051166.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/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/html/freestat/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/freestat/rcomp/tmp/11jyew1291051166.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/freestat/rcomp/tmp/125zc21291051166.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/freestat/rcomp/tmp/1319at1291051166.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/freestat/rcomp/tmp/144rrz1291051166.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/freestat/rcomp/tmp/15797n1291051166.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/freestat/rcomp/tmp/16bsob1291051166.tab")
+ }
>
> try(system("convert tmp/19fif1291051166.ps tmp/19fif1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/2koii1291051166.ps tmp/2koii1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/3koii1291051166.ps tmp/3koii1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/4koii1291051166.ps tmp/4koii1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/5koii1291051166.ps tmp/5koii1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/6cfz31291051166.ps tmp/6cfz31291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/75pgo1291051166.ps tmp/75pgo1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/85pgo1291051166.ps tmp/85pgo1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/95pgo1291051166.ps tmp/95pgo1291051166.png",intern=TRUE))
character(0)
> try(system("convert tmp/10yyfq1291051166.ps tmp/10yyfq1291051166.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.950 2.488 4.405