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(0.61328
+ ,1.5334
+ ,0.62168
+ ,0.62915
+ ,0.634
+ ,0.6348
+ ,0.6089
+ ,1.5225
+ ,0.61328
+ ,0.62168
+ ,0.62915
+ ,0.634
+ ,0.60857
+ ,1.5135
+ ,0.6089
+ ,0.61328
+ ,0.62168
+ ,0.62915
+ ,0.62672
+ ,1.5144
+ ,0.60857
+ ,0.6089
+ ,0.61328
+ ,0.62168
+ ,0.62291
+ ,1.4913
+ ,0.62672
+ ,0.60857
+ ,0.6089
+ ,0.61328
+ ,0.62393
+ ,1.4793
+ ,0.62291
+ ,0.62672
+ ,0.60857
+ ,0.6089
+ ,0.61838
+ ,1.4663
+ ,0.62393
+ ,0.62291
+ ,0.62672
+ ,0.60857
+ ,0.62012
+ ,1.4749
+ ,0.61838
+ ,0.62393
+ ,0.62291
+ ,0.62672
+ ,0.61659
+ ,1.4745
+ ,0.62012
+ ,0.61838
+ ,0.62393
+ ,0.62291
+ ,0.6116
+ ,1.4775
+ ,0.61659
+ ,0.62012
+ ,0.61838
+ ,0.62393
+ ,0.61573
+ ,1.4678
+ ,0.6116
+ ,0.61659
+ ,0.62012
+ ,0.61838
+ ,0.61407
+ ,1.4658
+ ,0.61573
+ ,0.6116
+ ,0.61659
+ ,0.62012
+ ,0.62823
+ ,1.4572
+ ,0.61407
+ ,0.61573
+ ,0.6116
+ ,0.61659
+ ,0.64405
+ ,1.4721
+ ,0.62823
+ ,0.61407
+ ,0.61573
+ ,0.6116
+ ,0.6387
+ ,1.4624
+ ,0.64405
+ ,0.62823
+ ,0.61407
+ ,0.61573
+ ,0.63633
+ ,1.4636
+ ,0.6387
+ ,0.64405
+ ,0.62823
+ ,0.61407
+ ,0.63059
+ ,1.4649
+ ,0.63633
+ ,0.6387
+ ,0.64405
+ ,0.62823
+ ,0.62994
+ ,1.465
+ ,0.63059
+ ,0.63633
+ ,0.6387
+ ,0.64405
+ ,0.63709
+ ,1.4673
+ ,0.62994
+ ,0.63059
+ ,0.63633
+ ,0.6387
+ ,0.64217
+ ,1.4679
+ ,0.63709
+ ,0.62994
+ ,0.63059
+ ,0.63633
+ ,0.65711
+ ,1.4621
+ ,0.64217
+ ,0.63709
+ ,0.62994
+ ,0.63059
+ ,0.66977
+ ,1.4674
+ ,0.65711
+ ,0.64217
+ ,0.63709
+ ,0.62994
+ ,0.68255
+ ,1.4695
+ ,0.66977
+ ,0.65711
+ ,0.64217
+ ,0.63709
+ ,0.68902
+ ,1.4964
+ ,0.68255
+ ,0.66977
+ ,0.65711
+ ,0.64217
+ ,0.71322
+ ,1.5155
+ ,0.68902
+ ,0.68255
+ ,0.66977
+ ,0.65711
+ ,0.70224
+ ,1.5411
+ ,0.71322
+ ,0.68902
+ ,0.68255
+ ,0.66977
+ ,0.70045
+ ,1.5476
+ ,0.70224
+ ,0.71322
+ ,0.68902
+ ,0.68255
+ ,0.69919
+ ,1.54
+ ,0.70045
+ ,0.70224
+ ,0.71322
+ ,0.68902
+ ,0.69693
+ ,1.5474
+ ,0.69919
+ ,0.70045
+ ,0.70224
+ ,0.71322
+ ,0.69763
+ ,1.5485
+ ,0.69693
+ ,0.69919
+ ,0.70045
+ ,0.70224
+ ,0.69278
+ ,1.559
+ ,0.69763
+ ,0.69693
+ ,0.69919
+ ,0.70045
+ ,0.70196
+ ,1.5544
+ ,0.69278
+ ,0.69763
+ ,0.69693
+ ,0.69919
+ ,0.69215
+ ,1.5657
+ ,0.70196
+ ,0.69278
+ ,0.69763
+ ,0.69693
+ ,0.6769
+ ,1.5734
+ ,0.69215
+ ,0.70196
+ ,0.69278
+ ,0.69763
+ ,0.67124
+ ,1.567
+ ,0.6769
+ ,0.69215
+ ,0.70196
+ ,0.69278
+ ,0.66532
+ ,1.5547
+ ,0.67124
+ ,0.6769
+ ,0.69215
+ ,0.70196
+ ,0.67157
+ ,1.54
+ ,0.66532
+ ,0.67124
+ ,0.6769
+ ,0.69215
+ ,0.66428
+ ,1.5192
+ ,0.67157
+ ,0.66532
+ ,0.67124
+ ,0.6769
+ ,0.66576
+ ,1.527
+ ,0.66428
+ ,0.67157
+ ,0.66532
+ ,0.67124
+ ,0.66942
+ ,1.5387
+ ,0.66576
+ ,0.66428
+ ,0.67157
+ ,0.66532
+ ,0.6813
+ ,1.5431
+ ,0.66942
+ ,0.66576
+ ,0.66428
+ ,0.67157
+ ,0.69144
+ ,1.5426
+ ,0.6813
+ ,0.66942
+ ,0.66576
+ ,0.66428
+ ,0.69862
+ ,1.5216
+ ,0.69144
+ ,0.6813
+ ,0.66942
+ ,0.66576
+ ,0.695
+ ,1.5364
+ ,0.69862
+ ,0.69144
+ ,0.6813
+ ,0.66942
+ ,0.69867
+ ,1.5469
+ ,0.695
+ ,0.69862
+ ,0.69144
+ ,0.6813
+ ,0.68968
+ ,1.5501
+ ,0.69867
+ ,0.695
+ ,0.69862
+ ,0.69144
+ ,0.69233
+ ,1.5494
+ ,0.68968
+ ,0.69867
+ ,0.695
+ ,0.69862
+ ,0.68293
+ ,1.5475
+ ,0.69233
+ ,0.68968
+ ,0.69867
+ ,0.695
+ ,0.68399
+ ,1.5448
+ ,0.68293
+ ,0.69233
+ ,0.68968
+ ,0.69867
+ ,0.66895
+ ,1.5391
+ ,0.68399
+ ,0.68293
+ ,0.69233
+ ,0.68968
+ ,0.68756
+ ,1.5578
+ ,0.66895
+ ,0.68399
+ ,0.68293
+ ,0.69233
+ ,0.68527
+ ,1.5528
+ ,0.68756
+ ,0.66895
+ ,0.68399
+ ,0.68293
+ ,0.6776
+ ,1.5496
+ ,0.68527
+ ,0.68756
+ ,0.66895
+ ,0.68399
+ ,0.68137
+ ,1.549
+ ,0.6776
+ ,0.68527
+ ,0.68756
+ ,0.66895
+ ,0.67933
+ ,1.5449
+ ,0.68137
+ ,0.6776
+ ,0.68527
+ ,0.68756
+ ,0.67922
+ ,1.5479
+ ,0.67933
+ ,0.68137
+ ,0.6776
+ ,0.68527)
+ ,dim=c(6
+ ,56)
+ ,dimnames=list(c('Britse_pond'
+ ,'Zwitserse_frank'
+ ,'Britse_pond_-1'
+ ,'Britse_pond_-2'
+ ,'Britse_pond_-3'
+ ,'Britse_pond_-4')
+ ,1:56))
> y <- array(NA,dim=c(6,56),dimnames=list(c('Britse_pond','Zwitserse_frank','Britse_pond_-1','Britse_pond_-2','Britse_pond_-3','Britse_pond_-4'),1:56))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = 'Linear Trend'
> par2 = 'Include Monthly Dummies'
> par1 = '1'
> #'GNU S' R Code compiled by R2WASP v. 1.0.44 ()
> #Author: Prof. Dr. P. Wessa
> #To cite this work: AUTHOR(S), (YEAR), YOUR SOFTWARE TITLE (vNUMBER) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_YOURPAGE.wasp/
> #Source of accompanying publication: Office for Research, Development, and Education
> #Technical description: Write here your technical program description (don't use hard returns!)
> library(lattice)
> 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
Britse_pond Zwitserse_frank Britse_pond_-1 Britse_pond_-2 Britse_pond_-3
1 0.61328 1.5334 0.62168 0.62915 0.63400
2 0.60890 1.5225 0.61328 0.62168 0.62915
3 0.60857 1.5135 0.60890 0.61328 0.62168
4 0.62672 1.5144 0.60857 0.60890 0.61328
5 0.62291 1.4913 0.62672 0.60857 0.60890
6 0.62393 1.4793 0.62291 0.62672 0.60857
7 0.61838 1.4663 0.62393 0.62291 0.62672
8 0.62012 1.4749 0.61838 0.62393 0.62291
9 0.61659 1.4745 0.62012 0.61838 0.62393
10 0.61160 1.4775 0.61659 0.62012 0.61838
11 0.61573 1.4678 0.61160 0.61659 0.62012
12 0.61407 1.4658 0.61573 0.61160 0.61659
13 0.62823 1.4572 0.61407 0.61573 0.61160
14 0.64405 1.4721 0.62823 0.61407 0.61573
15 0.63870 1.4624 0.64405 0.62823 0.61407
16 0.63633 1.4636 0.63870 0.64405 0.62823
17 0.63059 1.4649 0.63633 0.63870 0.64405
18 0.62994 1.4650 0.63059 0.63633 0.63870
19 0.63709 1.4673 0.62994 0.63059 0.63633
20 0.64217 1.4679 0.63709 0.62994 0.63059
21 0.65711 1.4621 0.64217 0.63709 0.62994
22 0.66977 1.4674 0.65711 0.64217 0.63709
23 0.68255 1.4695 0.66977 0.65711 0.64217
24 0.68902 1.4964 0.68255 0.66977 0.65711
25 0.71322 1.5155 0.68902 0.68255 0.66977
26 0.70224 1.5411 0.71322 0.68902 0.68255
27 0.70045 1.5476 0.70224 0.71322 0.68902
28 0.69919 1.5400 0.70045 0.70224 0.71322
29 0.69693 1.5474 0.69919 0.70045 0.70224
30 0.69763 1.5485 0.69693 0.69919 0.70045
31 0.69278 1.5590 0.69763 0.69693 0.69919
32 0.70196 1.5544 0.69278 0.69763 0.69693
33 0.69215 1.5657 0.70196 0.69278 0.69763
34 0.67690 1.5734 0.69215 0.70196 0.69278
35 0.67124 1.5670 0.67690 0.69215 0.70196
36 0.66532 1.5547 0.67124 0.67690 0.69215
37 0.67157 1.5400 0.66532 0.67124 0.67690
38 0.66428 1.5192 0.67157 0.66532 0.67124
39 0.66576 1.5270 0.66428 0.67157 0.66532
40 0.66942 1.5387 0.66576 0.66428 0.67157
41 0.68130 1.5431 0.66942 0.66576 0.66428
42 0.69144 1.5426 0.68130 0.66942 0.66576
43 0.69862 1.5216 0.69144 0.68130 0.66942
44 0.69500 1.5364 0.69862 0.69144 0.68130
45 0.69867 1.5469 0.69500 0.69862 0.69144
46 0.68968 1.5501 0.69867 0.69500 0.69862
47 0.69233 1.5494 0.68968 0.69867 0.69500
48 0.68293 1.5475 0.69233 0.68968 0.69867
49 0.68399 1.5448 0.68293 0.69233 0.68968
50 0.66895 1.5391 0.68399 0.68293 0.69233
51 0.68756 1.5578 0.66895 0.68399 0.68293
52 0.68527 1.5528 0.68756 0.66895 0.68399
53 0.67760 1.5496 0.68527 0.68756 0.66895
54 0.68137 1.5490 0.67760 0.68527 0.68756
55 0.67933 1.5449 0.68137 0.67760 0.68527
56 0.67922 1.5479 0.67933 0.68137 0.67760
Britse_pond_-4 M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 t
1 0.63480 1 0 0 0 0 0 0 0 0 0 0 1
2 0.63400 0 1 0 0 0 0 0 0 0 0 0 2
3 0.62915 0 0 1 0 0 0 0 0 0 0 0 3
4 0.62168 0 0 0 1 0 0 0 0 0 0 0 4
5 0.61328 0 0 0 0 1 0 0 0 0 0 0 5
6 0.60890 0 0 0 0 0 1 0 0 0 0 0 6
7 0.60857 0 0 0 0 0 0 1 0 0 0 0 7
8 0.62672 0 0 0 0 0 0 0 1 0 0 0 8
9 0.62291 0 0 0 0 0 0 0 0 1 0 0 9
10 0.62393 0 0 0 0 0 0 0 0 0 1 0 10
11 0.61838 0 0 0 0 0 0 0 0 0 0 1 11
12 0.62012 0 0 0 0 0 0 0 0 0 0 0 12
13 0.61659 1 0 0 0 0 0 0 0 0 0 0 13
14 0.61160 0 1 0 0 0 0 0 0 0 0 0 14
15 0.61573 0 0 1 0 0 0 0 0 0 0 0 15
16 0.61407 0 0 0 1 0 0 0 0 0 0 0 16
17 0.62823 0 0 0 0 1 0 0 0 0 0 0 17
18 0.64405 0 0 0 0 0 1 0 0 0 0 0 18
19 0.63870 0 0 0 0 0 0 1 0 0 0 0 19
20 0.63633 0 0 0 0 0 0 0 1 0 0 0 20
21 0.63059 0 0 0 0 0 0 0 0 1 0 0 21
22 0.62994 0 0 0 0 0 0 0 0 0 1 0 22
23 0.63709 0 0 0 0 0 0 0 0 0 0 1 23
24 0.64217 0 0 0 0 0 0 0 0 0 0 0 24
25 0.65711 1 0 0 0 0 0 0 0 0 0 0 25
26 0.66977 0 1 0 0 0 0 0 0 0 0 0 26
27 0.68255 0 0 1 0 0 0 0 0 0 0 0 27
28 0.68902 0 0 0 1 0 0 0 0 0 0 0 28
29 0.71322 0 0 0 0 1 0 0 0 0 0 0 29
30 0.70224 0 0 0 0 0 1 0 0 0 0 0 30
31 0.70045 0 0 0 0 0 0 1 0 0 0 0 31
32 0.69919 0 0 0 0 0 0 0 1 0 0 0 32
33 0.69693 0 0 0 0 0 0 0 0 1 0 0 33
34 0.69763 0 0 0 0 0 0 0 0 0 1 0 34
35 0.69278 0 0 0 0 0 0 0 0 0 0 1 35
36 0.70196 0 0 0 0 0 0 0 0 0 0 0 36
37 0.69215 1 0 0 0 0 0 0 0 0 0 0 37
38 0.67690 0 1 0 0 0 0 0 0 0 0 0 38
39 0.67124 0 0 1 0 0 0 0 0 0 0 0 39
40 0.66532 0 0 0 1 0 0 0 0 0 0 0 40
41 0.67157 0 0 0 0 1 0 0 0 0 0 0 41
42 0.66428 0 0 0 0 0 1 0 0 0 0 0 42
43 0.66576 0 0 0 0 0 0 1 0 0 0 0 43
44 0.66942 0 0 0 0 0 0 0 1 0 0 0 44
45 0.68130 0 0 0 0 0 0 0 0 1 0 0 45
46 0.69144 0 0 0 0 0 0 0 0 0 1 0 46
47 0.69862 0 0 0 0 0 0 0 0 0 0 1 47
48 0.69500 0 0 0 0 0 0 0 0 0 0 0 48
49 0.69867 1 0 0 0 0 0 0 0 0 0 0 49
50 0.68968 0 1 0 0 0 0 0 0 0 0 0 50
51 0.69233 0 0 1 0 0 0 0 0 0 0 0 51
52 0.68293 0 0 0 1 0 0 0 0 0 0 0 52
53 0.68399 0 0 0 0 1 0 0 0 0 0 0 53
54 0.66895 0 0 0 0 0 1 0 0 0 0 0 54
55 0.68756 0 0 0 0 0 0 1 0 0 0 0 55
56 0.68527 0 0 0 0 0 0 0 1 0 0 0 56
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) Zwitserse_frank `Britse_pond_-1` `Britse_pond_-2`
0.1146370 -0.0565950 1.1244479 -0.0005468
`Britse_pond_-3` `Britse_pond_-4` M1 M2
-0.2590966 0.0820676 0.0101173 -0.0019884
M3 M4 M5 M6
0.0044024 0.0068426 -0.0002812 0.0051666
M7 M8 M9 M10
0.0023676 0.0038598 0.0033893 -0.0020740
M11 t
0.0064357 0.0001373
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-0.014919 -0.005367 -0.001137 0.003580 0.017985
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.1146370 0.0744975 1.539 0.132
Zwitserse_frank -0.0565950 0.0750536 -0.754 0.455
`Britse_pond_-1` 1.1244479 0.1617883 6.950 2.87e-08 ***
`Britse_pond_-2` -0.0005468 0.2393492 -0.002 0.998
`Britse_pond_-3` -0.2590966 0.2380478 -1.088 0.283
`Britse_pond_-4` 0.0820676 0.1762383 0.466 0.644
M1 0.0101173 0.0061411 1.647 0.108
M2 -0.0019884 0.0059576 -0.334 0.740
M3 0.0044024 0.0064515 0.682 0.499
M4 0.0068426 0.0060656 1.128 0.266
M5 -0.0002812 0.0060614 -0.046 0.963
M6 0.0051666 0.0060878 0.849 0.401
M7 0.0023676 0.0058835 0.402 0.690
M8 0.0038598 0.0060282 0.640 0.526
M9 0.0033893 0.0062308 0.544 0.590
M10 -0.0020740 0.0063087 -0.329 0.744
M11 0.0064357 0.0063824 1.008 0.320
t 0.0001373 0.0001304 1.053 0.299
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.008727 on 38 degrees of freedom
Multiple R-squared: 0.9478, Adjusted R-squared: 0.9244
F-statistic: 40.57 on 17 and 38 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.5173171 0.9653658 0.4826829
[2,] 0.7478903 0.5042194 0.2521097
[3,] 0.7378753 0.5242495 0.2621247
[4,] 0.6201589 0.7596821 0.3798411
[5,] 0.6562765 0.6874470 0.3437235
[6,] 0.8226998 0.3546005 0.1773002
[7,] 0.7262350 0.5475301 0.2737650
[8,] 0.6208881 0.7582238 0.3791119
[9,] 0.5185929 0.9628142 0.4814071
[10,] 0.4479718 0.8959436 0.5520282
[11,] 0.4589546 0.9179093 0.5410454
[12,] 0.4786123 0.9572246 0.5213877
[13,] 0.4721933 0.9443866 0.5278067
[14,] 0.4076860 0.8153720 0.5923140
[15,] 0.6624719 0.6750562 0.3375281
> postscript(file="/var/www/html/rcomp/tmp/17fge1258713237.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/2eynp1258713237.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/3yyjs1258713237.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/4b9qr1258713237.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/5z0ih1258713237.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 = 56
Frequency = 1
1 2 3 4 5
-0.0113608609 0.0038610116 -0.0001232975 0.0143054383 -0.0046797450
6 7 8 9 10
-0.0053559856 -0.0053992774 -0.0010375065 -0.0056395210 -0.0026851348
11 12 13 14 15
-0.0012357700 -0.0024145977 0.0018697201 0.0160579260 -0.0149191273
16 17 18 19 20
-0.0099692318 -0.0030502641 -0.0055111718 0.0049835291 -0.0008648912
21 22 23 24 25
0.0086744689 0.0120698767 0.0028237804 0.0062051333 0.0160173675
26 27 28 29 30
-0.0064811791 -0.0014441283 0.0020341627 0.0037644211 0.0019194321
31 32 33 34 35
-0.0006424533 0.0116195069 -0.0071759717 -0.0069423586 -0.0016925971
36 37 38 39 40
0.0010506721 -0.0002783915 -0.0040231154 -0.0014984000 0.0006832804
41 42 43 44 45
0.0132824380 0.0054343102 0.0035190205 -0.0061831694 0.0041410238
46 47 48 49 50
-0.0024423834 0.0001045867 -0.0048412078 -0.0062478351 -0.0094146431
51 52 53 54 55
0.0179849532 -0.0070536496 -0.0093168500 0.0035134150 -0.0024608189
56
-0.0035339398
> postscript(file="/var/www/html/rcomp/tmp/6u58e1258713237.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 = 56
Frequency = 1
lag(myerror, k = 1) myerror
0 -0.0113608609 NA
1 0.0038610116 -0.0113608609
2 -0.0001232975 0.0038610116
3 0.0143054383 -0.0001232975
4 -0.0046797450 0.0143054383
5 -0.0053559856 -0.0046797450
6 -0.0053992774 -0.0053559856
7 -0.0010375065 -0.0053992774
8 -0.0056395210 -0.0010375065
9 -0.0026851348 -0.0056395210
10 -0.0012357700 -0.0026851348
11 -0.0024145977 -0.0012357700
12 0.0018697201 -0.0024145977
13 0.0160579260 0.0018697201
14 -0.0149191273 0.0160579260
15 -0.0099692318 -0.0149191273
16 -0.0030502641 -0.0099692318
17 -0.0055111718 -0.0030502641
18 0.0049835291 -0.0055111718
19 -0.0008648912 0.0049835291
20 0.0086744689 -0.0008648912
21 0.0120698767 0.0086744689
22 0.0028237804 0.0120698767
23 0.0062051333 0.0028237804
24 0.0160173675 0.0062051333
25 -0.0064811791 0.0160173675
26 -0.0014441283 -0.0064811791
27 0.0020341627 -0.0014441283
28 0.0037644211 0.0020341627
29 0.0019194321 0.0037644211
30 -0.0006424533 0.0019194321
31 0.0116195069 -0.0006424533
32 -0.0071759717 0.0116195069
33 -0.0069423586 -0.0071759717
34 -0.0016925971 -0.0069423586
35 0.0010506721 -0.0016925971
36 -0.0002783915 0.0010506721
37 -0.0040231154 -0.0002783915
38 -0.0014984000 -0.0040231154
39 0.0006832804 -0.0014984000
40 0.0132824380 0.0006832804
41 0.0054343102 0.0132824380
42 0.0035190205 0.0054343102
43 -0.0061831694 0.0035190205
44 0.0041410238 -0.0061831694
45 -0.0024423834 0.0041410238
46 0.0001045867 -0.0024423834
47 -0.0048412078 0.0001045867
48 -0.0062478351 -0.0048412078
49 -0.0094146431 -0.0062478351
50 0.0179849532 -0.0094146431
51 -0.0070536496 0.0179849532
52 -0.0093168500 -0.0070536496
53 0.0035134150 -0.0093168500
54 -0.0024608189 0.0035134150
55 -0.0035339398 -0.0024608189
56 NA -0.0035339398
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 0.0038610116 -0.0113608609
[2,] -0.0001232975 0.0038610116
[3,] 0.0143054383 -0.0001232975
[4,] -0.0046797450 0.0143054383
[5,] -0.0053559856 -0.0046797450
[6,] -0.0053992774 -0.0053559856
[7,] -0.0010375065 -0.0053992774
[8,] -0.0056395210 -0.0010375065
[9,] -0.0026851348 -0.0056395210
[10,] -0.0012357700 -0.0026851348
[11,] -0.0024145977 -0.0012357700
[12,] 0.0018697201 -0.0024145977
[13,] 0.0160579260 0.0018697201
[14,] -0.0149191273 0.0160579260
[15,] -0.0099692318 -0.0149191273
[16,] -0.0030502641 -0.0099692318
[17,] -0.0055111718 -0.0030502641
[18,] 0.0049835291 -0.0055111718
[19,] -0.0008648912 0.0049835291
[20,] 0.0086744689 -0.0008648912
[21,] 0.0120698767 0.0086744689
[22,] 0.0028237804 0.0120698767
[23,] 0.0062051333 0.0028237804
[24,] 0.0160173675 0.0062051333
[25,] -0.0064811791 0.0160173675
[26,] -0.0014441283 -0.0064811791
[27,] 0.0020341627 -0.0014441283
[28,] 0.0037644211 0.0020341627
[29,] 0.0019194321 0.0037644211
[30,] -0.0006424533 0.0019194321
[31,] 0.0116195069 -0.0006424533
[32,] -0.0071759717 0.0116195069
[33,] -0.0069423586 -0.0071759717
[34,] -0.0016925971 -0.0069423586
[35,] 0.0010506721 -0.0016925971
[36,] -0.0002783915 0.0010506721
[37,] -0.0040231154 -0.0002783915
[38,] -0.0014984000 -0.0040231154
[39,] 0.0006832804 -0.0014984000
[40,] 0.0132824380 0.0006832804
[41,] 0.0054343102 0.0132824380
[42,] 0.0035190205 0.0054343102
[43,] -0.0061831694 0.0035190205
[44,] 0.0041410238 -0.0061831694
[45,] -0.0024423834 0.0041410238
[46,] 0.0001045867 -0.0024423834
[47,] -0.0048412078 0.0001045867
[48,] -0.0062478351 -0.0048412078
[49,] -0.0094146431 -0.0062478351
[50,] 0.0179849532 -0.0094146431
[51,] -0.0070536496 0.0179849532
[52,] -0.0093168500 -0.0070536496
[53,] 0.0035134150 -0.0093168500
[54,] -0.0024608189 0.0035134150
[55,] -0.0035339398 -0.0024608189
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 0.0038610116 -0.0113608609
2 -0.0001232975 0.0038610116
3 0.0143054383 -0.0001232975
4 -0.0046797450 0.0143054383
5 -0.0053559856 -0.0046797450
6 -0.0053992774 -0.0053559856
7 -0.0010375065 -0.0053992774
8 -0.0056395210 -0.0010375065
9 -0.0026851348 -0.0056395210
10 -0.0012357700 -0.0026851348
11 -0.0024145977 -0.0012357700
12 0.0018697201 -0.0024145977
13 0.0160579260 0.0018697201
14 -0.0149191273 0.0160579260
15 -0.0099692318 -0.0149191273
16 -0.0030502641 -0.0099692318
17 -0.0055111718 -0.0030502641
18 0.0049835291 -0.0055111718
19 -0.0008648912 0.0049835291
20 0.0086744689 -0.0008648912
21 0.0120698767 0.0086744689
22 0.0028237804 0.0120698767
23 0.0062051333 0.0028237804
24 0.0160173675 0.0062051333
25 -0.0064811791 0.0160173675
26 -0.0014441283 -0.0064811791
27 0.0020341627 -0.0014441283
28 0.0037644211 0.0020341627
29 0.0019194321 0.0037644211
30 -0.0006424533 0.0019194321
31 0.0116195069 -0.0006424533
32 -0.0071759717 0.0116195069
33 -0.0069423586 -0.0071759717
34 -0.0016925971 -0.0069423586
35 0.0010506721 -0.0016925971
36 -0.0002783915 0.0010506721
37 -0.0040231154 -0.0002783915
38 -0.0014984000 -0.0040231154
39 0.0006832804 -0.0014984000
40 0.0132824380 0.0006832804
41 0.0054343102 0.0132824380
42 0.0035190205 0.0054343102
43 -0.0061831694 0.0035190205
44 0.0041410238 -0.0061831694
45 -0.0024423834 0.0041410238
46 0.0001045867 -0.0024423834
47 -0.0048412078 0.0001045867
48 -0.0062478351 -0.0048412078
49 -0.0094146431 -0.0062478351
50 0.0179849532 -0.0094146431
51 -0.0070536496 0.0179849532
52 -0.0093168500 -0.0070536496
53 0.0035134150 -0.0093168500
54 -0.0024608189 0.0035134150
55 -0.0035339398 -0.0024608189
> 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/71vgv1258713237.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/8snvv1258713237.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/9agj71258713237.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/10l2u81258713237.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/11znym1258713237.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/1226bi1258713237.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/13d3641258713237.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/146dkt1258713237.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/15an0n1258713237.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/16zamk1258713237.tab")
+ }
>
> system("convert tmp/17fge1258713237.ps tmp/17fge1258713237.png")
> system("convert tmp/2eynp1258713237.ps tmp/2eynp1258713237.png")
> system("convert tmp/3yyjs1258713237.ps tmp/3yyjs1258713237.png")
> system("convert tmp/4b9qr1258713237.ps tmp/4b9qr1258713237.png")
> system("convert tmp/5z0ih1258713237.ps tmp/5z0ih1258713237.png")
> system("convert tmp/6u58e1258713237.ps tmp/6u58e1258713237.png")
> system("convert tmp/71vgv1258713237.ps tmp/71vgv1258713237.png")
> system("convert tmp/8snvv1258713237.ps tmp/8snvv1258713237.png")
> system("convert tmp/9agj71258713237.ps tmp/9agj71258713237.png")
> system("convert tmp/10l2u81258713237.ps tmp/10l2u81258713237.png")
>
>
> proc.time()
user system elapsed
2.307 1.554 2.994