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(101.5,0,99.2,0,107.8,0,92.3,0,99.2,0,101.6,0,87,0,71.4,0,104.7,0,115.1,0,102.5,0,75.3,0,96.7,1,94.6,1,98.6,1,99.5,1,92,1,93.6,1,89.3,1,66.9,1,108.8,1,113.2,1,105.5,1,77.8,1,102.1,1,97,1,95.5,1,99.3,1,86.4,1,92.4,1,85.7,1,61.9,1,104.9,1,107.9,1,95.6,1,79.8,1,94.8,1,93.7,1,108.1,1,96.9,1,88.8,1,106.7,1,86.8,1,69.8,1,110.9,1,105.4,1,99.2,1,84.4,1,87.2,1,91.9,1,97.9,1,94.5,1,85,1,100.3,1,78.7,1,65.8,1,104.8,1,96,1,103.3,1,82.9,1,91.4,1,94.5,1,109.3,1,92.1,1,99.3,1,109.6,1,87.5,1,73.1,1,110.7,1,111.6,1,110.7,1,84,1,101.6,1,102.1,1,113.9,1,99,1,100.4,1,109.5,1,93,1,76.8,1,105.3,1),dim=c(2,81),dimnames=list(c('Y','X'),1:81))
> y <- array(NA,dim=c(2,81),dimnames=list(c('Y','X'),1:81))
> 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
Y X
1 101.5 0
2 99.2 0
3 107.8 0
4 92.3 0
5 99.2 0
6 101.6 0
7 87.0 0
8 71.4 0
9 104.7 0
10 115.1 0
11 102.5 0
12 75.3 0
13 96.7 1
14 94.6 1
15 98.6 1
16 99.5 1
17 92.0 1
18 93.6 1
19 89.3 1
20 66.9 1
21 108.8 1
22 113.2 1
23 105.5 1
24 77.8 1
25 102.1 1
26 97.0 1
27 95.5 1
28 99.3 1
29 86.4 1
30 92.4 1
31 85.7 1
32 61.9 1
33 104.9 1
34 107.9 1
35 95.6 1
36 79.8 1
37 94.8 1
38 93.7 1
39 108.1 1
40 96.9 1
41 88.8 1
42 106.7 1
43 86.8 1
44 69.8 1
45 110.9 1
46 105.4 1
47 99.2 1
48 84.4 1
49 87.2 1
50 91.9 1
51 97.9 1
52 94.5 1
53 85.0 1
54 100.3 1
55 78.7 1
56 65.8 1
57 104.8 1
58 96.0 1
59 103.3 1
60 82.9 1
61 91.4 1
62 94.5 1
63 109.3 1
64 92.1 1
65 99.3 1
66 109.6 1
67 87.5 1
68 73.1 1
69 110.7 1
70 111.6 1
71 110.7 1
72 84.0 1
73 101.6 1
74 102.1 1
75 113.9 1
76 99.0 1
77 100.4 1
78 109.5 1
79 93.0 1
80 76.8 1
81 105.3 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X
96.467 -1.478
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-33.088 -7.488 1.912 8.312 18.912
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 96.467 3.509 27.491 <2e-16 ***
X -1.478 3.802 -0.389 0.698
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 12.16 on 79 degrees of freedom
Multiple R-squared: 0.00191, Adjusted R-squared: -0.01072
F-statistic: 0.1512 on 1 and 79 DF, p-value: 0.6985
> 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.15071808 0.3014362 0.8492819
[2,] 0.06341651 0.1268330 0.9365835
[3,] 0.12028435 0.2405687 0.8797157
[4,] 0.57359708 0.8528058 0.4264029
[5,] 0.51401230 0.9719754 0.4859877
[6,] 0.63696544 0.7260691 0.3630346
[7,] 0.57388818 0.8522236 0.4261118
[8,] 0.72065138 0.5586972 0.2793486
[9,] 0.63640051 0.7271990 0.3635995
[10,] 0.54827863 0.9034427 0.4517214
[11,] 0.46290456 0.9258091 0.5370954
[12,] 0.38205898 0.7641180 0.6179410
[13,] 0.31398187 0.6279637 0.6860181
[14,] 0.24590103 0.4918021 0.7540990
[15,] 0.20019070 0.4003814 0.7998093
[16,] 0.47331648 0.9466330 0.5266835
[17,] 0.52096015 0.9580797 0.4790399
[18,] 0.61127290 0.7774542 0.3887271
[19,] 0.58517995 0.8296401 0.4148200
[20,] 0.65013116 0.6997377 0.3498688
[21,] 0.60423378 0.7915324 0.3957662
[22,] 0.53633083 0.9273383 0.4636692
[23,] 0.46649132 0.9329826 0.5335087
[24,] 0.40519798 0.8103960 0.5948020
[25,] 0.37292612 0.7458522 0.6270739
[26,] 0.31282092 0.6256418 0.6871791
[27,] 0.28659996 0.5731999 0.7134000
[28,] 0.66582454 0.6683509 0.3341755
[29,] 0.65018357 0.6996329 0.3498164
[30,] 0.65972470 0.6805506 0.3402753
[31,] 0.59810467 0.8037907 0.4018953
[32,] 0.62522017 0.7495597 0.3747798
[33,] 0.56234054 0.8753189 0.4376595
[34,] 0.49829474 0.9965895 0.5017053
[35,] 0.50935157 0.9812969 0.4906484
[36,] 0.44649054 0.8929811 0.5535095
[37,] 0.39811182 0.7962236 0.6018882
[38,] 0.39247141 0.7849428 0.6075286
[39,] 0.35727262 0.7145452 0.6427274
[40,] 0.55858680 0.8828264 0.4414132
[41,] 0.60041881 0.7991624 0.3995812
[42,] 0.58114467 0.8377107 0.4188553
[43,] 0.52317714 0.9536457 0.4768229
[44,] 0.50691639 0.9861672 0.4930836
[45,] 0.46885101 0.9377020 0.5311490
[46,] 0.40845339 0.8169068 0.5915466
[47,] 0.34773592 0.6954718 0.6522641
[48,] 0.28827825 0.5765565 0.7117218
[49,] 0.27183507 0.5436701 0.7281649
[50,] 0.22583704 0.4516741 0.7741630
[51,] 0.27304556 0.5460911 0.7269544
[52,] 0.62786754 0.7442649 0.3721325
[53,] 0.59036391 0.8192722 0.4096361
[54,] 0.52092122 0.9581576 0.4790788
[55,] 0.46908806 0.9381761 0.5309119
[56,] 0.49879785 0.9975957 0.5012022
[57,] 0.44551870 0.8910374 0.5544813
[58,] 0.37894660 0.7578932 0.6210534
[59,] 0.37073904 0.7414781 0.6292610
[60,] 0.31483819 0.6296764 0.6851618
[61,] 0.24856113 0.4971223 0.7514389
[62,] 0.23990996 0.4798199 0.7600900
[63,] 0.21926087 0.4385217 0.7807391
[64,] 0.50293239 0.9941352 0.4970676
[65,] 0.48442384 0.9688477 0.5155762
[66,] 0.48568685 0.9713737 0.5143131
[67,] 0.48286440 0.9657288 0.5171356
[68,] 0.53886511 0.9222698 0.4611349
[69,] 0.42253400 0.8450680 0.5774660
[70,] 0.30786624 0.6157325 0.6921338
[71,] 0.36102148 0.7220430 0.6389785
[72,] 0.22615495 0.4523099 0.7738450
> postscript(file="/var/www/html/rcomp/tmp/1laco1260960560.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/2qajl1260960560.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/31nov1260960560.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/4bjs71260960560.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/59nsl1260960560.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 = 81
Frequency = 1
1 2 3 4 5 6
5.0333333 2.7333333 11.3333333 -4.1666667 2.7333333 5.1333333
7 8 9 10 11 12
-9.4666667 -25.0666667 8.2333333 18.6333333 6.0333333 -21.1666667
13 14 15 16 17 18
1.7115942 -0.3884058 3.6115942 4.5115942 -2.9884058 -1.3884058
19 20 21 22 23 24
-5.6884058 -28.0884058 13.8115942 18.2115942 10.5115942 -17.1884058
25 26 27 28 29 30
7.1115942 2.0115942 0.5115942 4.3115942 -8.5884058 -2.5884058
31 32 33 34 35 36
-9.2884058 -33.0884058 9.9115942 12.9115942 0.6115942 -15.1884058
37 38 39 40 41 42
-0.1884058 -1.2884058 13.1115942 1.9115942 -6.1884058 11.7115942
43 44 45 46 47 48
-8.1884058 -25.1884058 15.9115942 10.4115942 4.2115942 -10.5884058
49 50 51 52 53 54
-7.7884058 -3.0884058 2.9115942 -0.4884058 -9.9884058 5.3115942
55 56 57 58 59 60
-16.2884058 -29.1884058 9.8115942 1.0115942 8.3115942 -12.0884058
61 62 63 64 65 66
-3.5884058 -0.4884058 14.3115942 -2.8884058 4.3115942 14.6115942
67 68 69 70 71 72
-7.4884058 -21.8884058 15.7115942 16.6115942 15.7115942 -10.9884058
73 74 75 76 77 78
6.6115942 7.1115942 18.9115942 4.0115942 5.4115942 14.5115942
79 80 81
-1.9884058 -18.1884058 10.3115942
> postscript(file="/var/www/html/rcomp/tmp/6rtdr1260960560.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 = 81
Frequency = 1
lag(myerror, k = 1) myerror
0 5.0333333 NA
1 2.7333333 5.0333333
2 11.3333333 2.7333333
3 -4.1666667 11.3333333
4 2.7333333 -4.1666667
5 5.1333333 2.7333333
6 -9.4666667 5.1333333
7 -25.0666667 -9.4666667
8 8.2333333 -25.0666667
9 18.6333333 8.2333333
10 6.0333333 18.6333333
11 -21.1666667 6.0333333
12 1.7115942 -21.1666667
13 -0.3884058 1.7115942
14 3.6115942 -0.3884058
15 4.5115942 3.6115942
16 -2.9884058 4.5115942
17 -1.3884058 -2.9884058
18 -5.6884058 -1.3884058
19 -28.0884058 -5.6884058
20 13.8115942 -28.0884058
21 18.2115942 13.8115942
22 10.5115942 18.2115942
23 -17.1884058 10.5115942
24 7.1115942 -17.1884058
25 2.0115942 7.1115942
26 0.5115942 2.0115942
27 4.3115942 0.5115942
28 -8.5884058 4.3115942
29 -2.5884058 -8.5884058
30 -9.2884058 -2.5884058
31 -33.0884058 -9.2884058
32 9.9115942 -33.0884058
33 12.9115942 9.9115942
34 0.6115942 12.9115942
35 -15.1884058 0.6115942
36 -0.1884058 -15.1884058
37 -1.2884058 -0.1884058
38 13.1115942 -1.2884058
39 1.9115942 13.1115942
40 -6.1884058 1.9115942
41 11.7115942 -6.1884058
42 -8.1884058 11.7115942
43 -25.1884058 -8.1884058
44 15.9115942 -25.1884058
45 10.4115942 15.9115942
46 4.2115942 10.4115942
47 -10.5884058 4.2115942
48 -7.7884058 -10.5884058
49 -3.0884058 -7.7884058
50 2.9115942 -3.0884058
51 -0.4884058 2.9115942
52 -9.9884058 -0.4884058
53 5.3115942 -9.9884058
54 -16.2884058 5.3115942
55 -29.1884058 -16.2884058
56 9.8115942 -29.1884058
57 1.0115942 9.8115942
58 8.3115942 1.0115942
59 -12.0884058 8.3115942
60 -3.5884058 -12.0884058
61 -0.4884058 -3.5884058
62 14.3115942 -0.4884058
63 -2.8884058 14.3115942
64 4.3115942 -2.8884058
65 14.6115942 4.3115942
66 -7.4884058 14.6115942
67 -21.8884058 -7.4884058
68 15.7115942 -21.8884058
69 16.6115942 15.7115942
70 15.7115942 16.6115942
71 -10.9884058 15.7115942
72 6.6115942 -10.9884058
73 7.1115942 6.6115942
74 18.9115942 7.1115942
75 4.0115942 18.9115942
76 5.4115942 4.0115942
77 14.5115942 5.4115942
78 -1.9884058 14.5115942
79 -18.1884058 -1.9884058
80 10.3115942 -18.1884058
81 NA 10.3115942
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 2.7333333 5.0333333
[2,] 11.3333333 2.7333333
[3,] -4.1666667 11.3333333
[4,] 2.7333333 -4.1666667
[5,] 5.1333333 2.7333333
[6,] -9.4666667 5.1333333
[7,] -25.0666667 -9.4666667
[8,] 8.2333333 -25.0666667
[9,] 18.6333333 8.2333333
[10,] 6.0333333 18.6333333
[11,] -21.1666667 6.0333333
[12,] 1.7115942 -21.1666667
[13,] -0.3884058 1.7115942
[14,] 3.6115942 -0.3884058
[15,] 4.5115942 3.6115942
[16,] -2.9884058 4.5115942
[17,] -1.3884058 -2.9884058
[18,] -5.6884058 -1.3884058
[19,] -28.0884058 -5.6884058
[20,] 13.8115942 -28.0884058
[21,] 18.2115942 13.8115942
[22,] 10.5115942 18.2115942
[23,] -17.1884058 10.5115942
[24,] 7.1115942 -17.1884058
[25,] 2.0115942 7.1115942
[26,] 0.5115942 2.0115942
[27,] 4.3115942 0.5115942
[28,] -8.5884058 4.3115942
[29,] -2.5884058 -8.5884058
[30,] -9.2884058 -2.5884058
[31,] -33.0884058 -9.2884058
[32,] 9.9115942 -33.0884058
[33,] 12.9115942 9.9115942
[34,] 0.6115942 12.9115942
[35,] -15.1884058 0.6115942
[36,] -0.1884058 -15.1884058
[37,] -1.2884058 -0.1884058
[38,] 13.1115942 -1.2884058
[39,] 1.9115942 13.1115942
[40,] -6.1884058 1.9115942
[41,] 11.7115942 -6.1884058
[42,] -8.1884058 11.7115942
[43,] -25.1884058 -8.1884058
[44,] 15.9115942 -25.1884058
[45,] 10.4115942 15.9115942
[46,] 4.2115942 10.4115942
[47,] -10.5884058 4.2115942
[48,] -7.7884058 -10.5884058
[49,] -3.0884058 -7.7884058
[50,] 2.9115942 -3.0884058
[51,] -0.4884058 2.9115942
[52,] -9.9884058 -0.4884058
[53,] 5.3115942 -9.9884058
[54,] -16.2884058 5.3115942
[55,] -29.1884058 -16.2884058
[56,] 9.8115942 -29.1884058
[57,] 1.0115942 9.8115942
[58,] 8.3115942 1.0115942
[59,] -12.0884058 8.3115942
[60,] -3.5884058 -12.0884058
[61,] -0.4884058 -3.5884058
[62,] 14.3115942 -0.4884058
[63,] -2.8884058 14.3115942
[64,] 4.3115942 -2.8884058
[65,] 14.6115942 4.3115942
[66,] -7.4884058 14.6115942
[67,] -21.8884058 -7.4884058
[68,] 15.7115942 -21.8884058
[69,] 16.6115942 15.7115942
[70,] 15.7115942 16.6115942
[71,] -10.9884058 15.7115942
[72,] 6.6115942 -10.9884058
[73,] 7.1115942 6.6115942
[74,] 18.9115942 7.1115942
[75,] 4.0115942 18.9115942
[76,] 5.4115942 4.0115942
[77,] 14.5115942 5.4115942
[78,] -1.9884058 14.5115942
[79,] -18.1884058 -1.9884058
[80,] 10.3115942 -18.1884058
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 2.7333333 5.0333333
2 11.3333333 2.7333333
3 -4.1666667 11.3333333
4 2.7333333 -4.1666667
5 5.1333333 2.7333333
6 -9.4666667 5.1333333
7 -25.0666667 -9.4666667
8 8.2333333 -25.0666667
9 18.6333333 8.2333333
10 6.0333333 18.6333333
11 -21.1666667 6.0333333
12 1.7115942 -21.1666667
13 -0.3884058 1.7115942
14 3.6115942 -0.3884058
15 4.5115942 3.6115942
16 -2.9884058 4.5115942
17 -1.3884058 -2.9884058
18 -5.6884058 -1.3884058
19 -28.0884058 -5.6884058
20 13.8115942 -28.0884058
21 18.2115942 13.8115942
22 10.5115942 18.2115942
23 -17.1884058 10.5115942
24 7.1115942 -17.1884058
25 2.0115942 7.1115942
26 0.5115942 2.0115942
27 4.3115942 0.5115942
28 -8.5884058 4.3115942
29 -2.5884058 -8.5884058
30 -9.2884058 -2.5884058
31 -33.0884058 -9.2884058
32 9.9115942 -33.0884058
33 12.9115942 9.9115942
34 0.6115942 12.9115942
35 -15.1884058 0.6115942
36 -0.1884058 -15.1884058
37 -1.2884058 -0.1884058
38 13.1115942 -1.2884058
39 1.9115942 13.1115942
40 -6.1884058 1.9115942
41 11.7115942 -6.1884058
42 -8.1884058 11.7115942
43 -25.1884058 -8.1884058
44 15.9115942 -25.1884058
45 10.4115942 15.9115942
46 4.2115942 10.4115942
47 -10.5884058 4.2115942
48 -7.7884058 -10.5884058
49 -3.0884058 -7.7884058
50 2.9115942 -3.0884058
51 -0.4884058 2.9115942
52 -9.9884058 -0.4884058
53 5.3115942 -9.9884058
54 -16.2884058 5.3115942
55 -29.1884058 -16.2884058
56 9.8115942 -29.1884058
57 1.0115942 9.8115942
58 8.3115942 1.0115942
59 -12.0884058 8.3115942
60 -3.5884058 -12.0884058
61 -0.4884058 -3.5884058
62 14.3115942 -0.4884058
63 -2.8884058 14.3115942
64 4.3115942 -2.8884058
65 14.6115942 4.3115942
66 -7.4884058 14.6115942
67 -21.8884058 -7.4884058
68 15.7115942 -21.8884058
69 16.6115942 15.7115942
70 15.7115942 16.6115942
71 -10.9884058 15.7115942
72 6.6115942 -10.9884058
73 7.1115942 6.6115942
74 18.9115942 7.1115942
75 4.0115942 18.9115942
76 5.4115942 4.0115942
77 14.5115942 5.4115942
78 -1.9884058 14.5115942
79 -18.1884058 -1.9884058
80 10.3115942 -18.1884058
> 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/7mnk01260960560.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/8gwws1260960560.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/9ejkj1260960560.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/10g5to1260960560.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/11k48o1260960560.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/12fv291260960560.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/13yjiu1260960560.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/14rgsh1260960560.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/15dutd1260960560.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/16hq4c1260960560.tab")
+ }
>
> try(system("convert tmp/1laco1260960560.ps tmp/1laco1260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/2qajl1260960560.ps tmp/2qajl1260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/31nov1260960560.ps tmp/31nov1260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/4bjs71260960560.ps tmp/4bjs71260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/59nsl1260960560.ps tmp/59nsl1260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/6rtdr1260960560.ps tmp/6rtdr1260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/7mnk01260960560.ps tmp/7mnk01260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/8gwws1260960560.ps tmp/8gwws1260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/9ejkj1260960560.ps tmp/9ejkj1260960560.png",intern=TRUE))
character(0)
> try(system("convert tmp/10g5to1260960560.ps tmp/10g5to1260960560.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.629 1.578 4.114