R version 2.13.0 (2011-04-13)
Copyright (C) 2011 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: i486-pc-linux-gnu (32-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> x <- array(list(79,30,115,108,30,116,43,26,100,78,38,140,86,44,166,44,30,99,104,40,139,158,47,181,102,30,116,77,31,116,80,30,108,123,34,129,73,31,118,105,33,125,107,33,127,84,36,136,33,14,46,42,17,54,96,32,124,106,30,115,56,35,128,59,28,97,76,34,125,91,39,149,115,39,149,76,29,108,101,44,166,94,21,80,92,28,107,75,28,107,128,38,146,56,32,123,41,29,111,67,27,105,77,40,155,66,40,155,69,28,104,105,34,132,116,33,127,62,33,122,100,35,87,67,29,109,46,20,78,135,37,141,124,33,124,58,29,112,68,28,108,37,21,78,93,41,158,56,20,78,83,30,119,59,22,88,133,42,155,106,32,123,71,36,136,116,31,117,98,33,124,64,40,151,32,38,145,25,24,87,46,43,165,63,31,120,95,40,150,113,37,136,111,31,116,120,39,150,87,32,118,25,18,71,131,39,144,47,30,110,109,37,147,37,32,111,15,17,68,54,12,48,16,13,51,22,17,68,37,17,64,29,20,76,55,17,66,5,17,68,0,17,66,27,22,83,37,15,55,29,12,41,17,17,66),dim=c(3,85),dimnames=list(c('blog','reviews','fdb120'),1:85))
> y <- array(NA,dim=c(3,85),dimnames=list(c('blog','reviews','fdb120'),1:85))
> 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 = '3'
> #'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
> 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
fdb120 blog reviews
1 115 79 30
2 116 108 30
3 100 43 26
4 140 78 38
5 166 86 44
6 99 44 30
7 139 104 40
8 181 158 47
9 116 102 30
10 116 77 31
11 108 80 30
12 129 123 34
13 118 73 31
14 125 105 33
15 127 107 33
16 136 84 36
17 46 33 14
18 54 42 17
19 124 96 32
20 115 106 30
21 128 56 35
22 97 59 28
23 125 76 34
24 149 91 39
25 149 115 39
26 108 76 29
27 166 101 44
28 80 94 21
29 107 92 28
30 107 75 28
31 146 128 38
32 123 56 32
33 111 41 29
34 105 67 27
35 155 77 40
36 155 66 40
37 104 69 28
38 132 105 34
39 127 116 33
40 122 62 33
41 87 100 35
42 109 67 29
43 78 46 20
44 141 135 37
45 124 124 33
46 112 58 29
47 108 68 28
48 78 37 21
49 158 93 41
50 78 56 20
51 119 83 30
52 88 59 22
53 155 133 42
54 123 106 32
55 136 71 36
56 117 116 31
57 124 98 33
58 151 64 40
59 145 32 38
60 87 25 24
61 165 46 43
62 120 63 31
63 150 95 40
64 136 113 37
65 116 111 31
66 150 120 39
67 118 87 32
68 71 25 18
69 144 131 39
70 110 47 30
71 147 109 37
72 111 37 32
73 68 15 17
74 48 54 12
75 51 16 13
76 68 22 17
77 64 37 17
78 76 29 20
79 66 55 17
80 68 5 17
81 66 0 17
82 83 27 22
83 55 37 15
84 41 29 12
85 66 17 17
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) blog reviews
0.71510 0.01008 3.71090
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-44.604 -1.097 1.511 2.930 7.883
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.71510 2.52013 0.284 0.777
blog 0.01008 0.02782 0.362 0.718
reviews 3.71090 0.11355 32.682 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 6.287 on 82 degrees of freedom
Multiple R-squared: 0.9635, Adjusted R-squared: 0.9626
F-statistic: 1082 on 2 and 82 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.6459300942 7.081398e-01 3.540699e-01
[2,] 0.8194579661 3.610841e-01 1.805420e-01
[3,] 0.7328587384 5.342825e-01 2.671413e-01
[4,] 0.6218288663 7.563423e-01 3.781711e-01
[5,] 0.5044449142 9.911102e-01 4.955551e-01
[6,] 0.4417533753 8.835068e-01 5.582466e-01
[7,] 0.3451382191 6.902764e-01 6.548618e-01
[8,] 0.2827918423 5.655837e-01 7.172082e-01
[9,] 0.2041036257 4.082073e-01 7.958964e-01
[10,] 0.1465906717 2.931813e-01 8.534093e-01
[11,] 0.1090870476 2.181741e-01 8.909130e-01
[12,] 0.0985409630 1.970819e-01 9.014590e-01
[13,] 0.1159479741 2.318959e-01 8.840520e-01
[14,] 0.0940571365 1.881143e-01 9.059429e-01
[15,] 0.0640037506 1.280075e-01 9.359962e-01
[16,] 0.0441230553 8.824611e-02 9.558769e-01
[17,] 0.0421208906 8.424178e-02 9.578791e-01
[18,] 0.0275642738 5.512855e-02 9.724357e-01
[19,] 0.0211149409 4.222988e-02 9.788851e-01
[20,] 0.0133179233 2.663585e-02 9.866821e-01
[21,] 0.0081600156 1.632003e-02 9.918400e-01
[22,] 0.0048375249 9.675050e-03 9.951625e-01
[23,] 0.0027622940 5.524588e-03 9.972377e-01
[24,] 0.0016276699 3.255340e-03 9.983723e-01
[25,] 0.0011460340 2.292068e-03 9.988540e-01
[26,] 0.0006527070 1.305414e-03 9.993473e-01
[27,] 0.0008812541 1.762508e-03 9.991187e-01
[28,] 0.0010832469 2.166494e-03 9.989168e-01
[29,] 0.0009975283 1.995057e-03 9.990025e-01
[30,] 0.0010394200 2.078840e-03 9.989606e-01
[31,] 0.0011029190 2.205838e-03 9.988971e-01
[32,] 0.0006201090 1.240218e-03 9.993799e-01
[33,] 0.0004241762 8.483525e-04 9.995758e-01
[34,] 0.0002468482 4.936965e-04 9.997532e-01
[35,] 0.0001324535 2.649070e-04 9.998675e-01
[36,] 0.9999999940 1.192249e-08 5.961243e-09
[37,] 0.9999999850 2.998544e-08 1.499272e-08
[38,] 0.9999999696 6.081764e-08 3.040882e-08
[39,] 0.9999999250 1.500096e-07 7.500482e-08
[40,] 0.9999998177 3.646110e-07 1.823055e-07
[41,] 0.9999996508 6.983814e-07 3.491907e-07
[42,] 0.9999992925 1.415047e-06 7.075235e-07
[43,] 0.9999985581 2.883845e-06 1.441922e-06
[44,] 0.9999978500 4.299919e-06 2.149959e-06
[45,] 0.9999957132 8.573602e-06 4.286801e-06
[46,] 0.9999968265 6.346991e-06 3.173495e-06
[47,] 0.9999966177 6.764623e-06 3.382312e-06
[48,] 0.9999942827 1.143462e-05 5.717312e-06
[49,] 0.9999889559 2.208826e-05 1.104413e-05
[50,] 0.9999751188 4.976250e-05 2.488125e-05
[51,] 0.9999447591 1.104819e-04 5.524093e-05
[52,] 0.9998809942 2.380116e-04 1.190058e-04
[53,] 0.9997530793 4.938414e-04 2.469207e-04
[54,] 0.9995636662 8.726677e-04 4.363338e-04
[55,] 0.9994103927 1.179215e-03 5.896073e-04
[56,] 0.9992923803 1.415239e-03 7.076197e-04
[57,] 0.9990487073 1.902585e-03 9.512927e-04
[58,] 0.9981432107 3.713579e-03 1.856789e-03
[59,] 0.9971232347 5.753531e-03 2.876765e-03
[60,] 0.9949011994 1.019760e-02 5.098801e-03
[61,] 0.9929587820 1.408244e-02 7.041218e-03
[62,] 0.9883322761 2.333545e-02 1.166772e-02
[63,] 0.9819427357 3.611453e-02 1.805726e-02
[64,] 0.9764700740 4.705985e-02 2.352993e-02
[65,] 0.9648654151 7.026917e-02 3.513458e-02
[66,] 0.9968861046 6.227791e-03 3.113895e-03
[67,] 0.9992014910 1.597018e-03 7.985090e-04
[68,] 0.9984562009 3.087598e-03 1.543799e-03
[69,] 0.9984792245 3.041551e-03 1.520775e-03
[70,] 0.9964525136 7.094973e-03 3.547486e-03
[71,] 0.9951700508 9.659898e-03 4.829949e-03
[72,] 0.9844807930 3.103841e-02 1.551921e-02
[73,] 0.9572969817 8.540604e-02 4.270302e-02
[74,] 0.9710857635 5.782847e-02 2.891424e-02
> postscript(file="/var/wessaorg/rcomp/tmp/1jm961324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index')
> points(x[,1]-mysum$resid)
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/243fe1324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/3pefi1324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/49h8a1324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals')
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/5o3341324132041.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 = 85
Frequency = 1
1 2 3 4 5
2.161918705 2.869673050 2.368288831 -2.515168689 1.138838197
6 7 8 9 10
-13.485371021 -11.198973243 4.280575966 2.930137668 -0.528822026
11 12 13 14 15
-4.848158731 0.874929089 1.511487720 0.767218548 2.747063675
16 17 18 19 20
0.846157901 -7.000189559 -10.223573298 3.568811079 1.889827923
21 22 23 24 25
-3.160778277 -8.214741359 -2.651431401 2.642929034 2.401070561
26 27 28 29 30
-1.096953382 0.987676651 0.408817594 1.452703240 1.624019659
31 32 33 34 35
2.980959492 2.971908535 2.255756892 3.415534754 5.073117540
36 37 38 39 40
5.183969340 -1.315515723 4.056322944 2.656366747 -1.799451687
41 42 43 44 45
-44.604185478 -0.006256454 2.603430145 1.621313041 -0.424252744
46 47 48 49 50
3.084440473 2.694561713 -1.016768532 4.200982954 2.502655781
51 52 53 54 55
6.121608960 5.050632264 -2.913010106 2.468036715 0.977164574
56 57 58 59 60
0.078157955 -0.162239397 1.204124213 2.948393385 -3.028526107
61 62 63 64 65
4.252831257 3.612262084 -0.108276315 -3.156983359 -0.871454863
66 67 68 69 70
3.350683379 -2.340491993 3.236847516 -2.760168421 -2.515603330
71 72 73 74 75
7.883326387 -8.836620174 4.048517484 2.209975484 1.882022463
76 77 78 79 80
3.977975429 -0.173186117 0.774746563 1.645420028 4.149291848
81 82 83 84 85
2.199679030 0.373110228 -1.751394909 -4.538088606 2.028362611
> postscript(file="/var/wessaorg/rcomp/tmp/6hqc41324132041.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 = 85
Frequency = 1
lag(myerror, k = 1) myerror
0 2.161918705 NA
1 2.869673050 2.161918705
2 2.368288831 2.869673050
3 -2.515168689 2.368288831
4 1.138838197 -2.515168689
5 -13.485371021 1.138838197
6 -11.198973243 -13.485371021
7 4.280575966 -11.198973243
8 2.930137668 4.280575966
9 -0.528822026 2.930137668
10 -4.848158731 -0.528822026
11 0.874929089 -4.848158731
12 1.511487720 0.874929089
13 0.767218548 1.511487720
14 2.747063675 0.767218548
15 0.846157901 2.747063675
16 -7.000189559 0.846157901
17 -10.223573298 -7.000189559
18 3.568811079 -10.223573298
19 1.889827923 3.568811079
20 -3.160778277 1.889827923
21 -8.214741359 -3.160778277
22 -2.651431401 -8.214741359
23 2.642929034 -2.651431401
24 2.401070561 2.642929034
25 -1.096953382 2.401070561
26 0.987676651 -1.096953382
27 0.408817594 0.987676651
28 1.452703240 0.408817594
29 1.624019659 1.452703240
30 2.980959492 1.624019659
31 2.971908535 2.980959492
32 2.255756892 2.971908535
33 3.415534754 2.255756892
34 5.073117540 3.415534754
35 5.183969340 5.073117540
36 -1.315515723 5.183969340
37 4.056322944 -1.315515723
38 2.656366747 4.056322944
39 -1.799451687 2.656366747
40 -44.604185478 -1.799451687
41 -0.006256454 -44.604185478
42 2.603430145 -0.006256454
43 1.621313041 2.603430145
44 -0.424252744 1.621313041
45 3.084440473 -0.424252744
46 2.694561713 3.084440473
47 -1.016768532 2.694561713
48 4.200982954 -1.016768532
49 2.502655781 4.200982954
50 6.121608960 2.502655781
51 5.050632264 6.121608960
52 -2.913010106 5.050632264
53 2.468036715 -2.913010106
54 0.977164574 2.468036715
55 0.078157955 0.977164574
56 -0.162239397 0.078157955
57 1.204124213 -0.162239397
58 2.948393385 1.204124213
59 -3.028526107 2.948393385
60 4.252831257 -3.028526107
61 3.612262084 4.252831257
62 -0.108276315 3.612262084
63 -3.156983359 -0.108276315
64 -0.871454863 -3.156983359
65 3.350683379 -0.871454863
66 -2.340491993 3.350683379
67 3.236847516 -2.340491993
68 -2.760168421 3.236847516
69 -2.515603330 -2.760168421
70 7.883326387 -2.515603330
71 -8.836620174 7.883326387
72 4.048517484 -8.836620174
73 2.209975484 4.048517484
74 1.882022463 2.209975484
75 3.977975429 1.882022463
76 -0.173186117 3.977975429
77 0.774746563 -0.173186117
78 1.645420028 0.774746563
79 4.149291848 1.645420028
80 2.199679030 4.149291848
81 0.373110228 2.199679030
82 -1.751394909 0.373110228
83 -4.538088606 -1.751394909
84 2.028362611 -4.538088606
85 NA 2.028362611
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 2.869673050 2.161918705
[2,] 2.368288831 2.869673050
[3,] -2.515168689 2.368288831
[4,] 1.138838197 -2.515168689
[5,] -13.485371021 1.138838197
[6,] -11.198973243 -13.485371021
[7,] 4.280575966 -11.198973243
[8,] 2.930137668 4.280575966
[9,] -0.528822026 2.930137668
[10,] -4.848158731 -0.528822026
[11,] 0.874929089 -4.848158731
[12,] 1.511487720 0.874929089
[13,] 0.767218548 1.511487720
[14,] 2.747063675 0.767218548
[15,] 0.846157901 2.747063675
[16,] -7.000189559 0.846157901
[17,] -10.223573298 -7.000189559
[18,] 3.568811079 -10.223573298
[19,] 1.889827923 3.568811079
[20,] -3.160778277 1.889827923
[21,] -8.214741359 -3.160778277
[22,] -2.651431401 -8.214741359
[23,] 2.642929034 -2.651431401
[24,] 2.401070561 2.642929034
[25,] -1.096953382 2.401070561
[26,] 0.987676651 -1.096953382
[27,] 0.408817594 0.987676651
[28,] 1.452703240 0.408817594
[29,] 1.624019659 1.452703240
[30,] 2.980959492 1.624019659
[31,] 2.971908535 2.980959492
[32,] 2.255756892 2.971908535
[33,] 3.415534754 2.255756892
[34,] 5.073117540 3.415534754
[35,] 5.183969340 5.073117540
[36,] -1.315515723 5.183969340
[37,] 4.056322944 -1.315515723
[38,] 2.656366747 4.056322944
[39,] -1.799451687 2.656366747
[40,] -44.604185478 -1.799451687
[41,] -0.006256454 -44.604185478
[42,] 2.603430145 -0.006256454
[43,] 1.621313041 2.603430145
[44,] -0.424252744 1.621313041
[45,] 3.084440473 -0.424252744
[46,] 2.694561713 3.084440473
[47,] -1.016768532 2.694561713
[48,] 4.200982954 -1.016768532
[49,] 2.502655781 4.200982954
[50,] 6.121608960 2.502655781
[51,] 5.050632264 6.121608960
[52,] -2.913010106 5.050632264
[53,] 2.468036715 -2.913010106
[54,] 0.977164574 2.468036715
[55,] 0.078157955 0.977164574
[56,] -0.162239397 0.078157955
[57,] 1.204124213 -0.162239397
[58,] 2.948393385 1.204124213
[59,] -3.028526107 2.948393385
[60,] 4.252831257 -3.028526107
[61,] 3.612262084 4.252831257
[62,] -0.108276315 3.612262084
[63,] -3.156983359 -0.108276315
[64,] -0.871454863 -3.156983359
[65,] 3.350683379 -0.871454863
[66,] -2.340491993 3.350683379
[67,] 3.236847516 -2.340491993
[68,] -2.760168421 3.236847516
[69,] -2.515603330 -2.760168421
[70,] 7.883326387 -2.515603330
[71,] -8.836620174 7.883326387
[72,] 4.048517484 -8.836620174
[73,] 2.209975484 4.048517484
[74,] 1.882022463 2.209975484
[75,] 3.977975429 1.882022463
[76,] -0.173186117 3.977975429
[77,] 0.774746563 -0.173186117
[78,] 1.645420028 0.774746563
[79,] 4.149291848 1.645420028
[80,] 2.199679030 4.149291848
[81,] 0.373110228 2.199679030
[82,] -1.751394909 0.373110228
[83,] -4.538088606 -1.751394909
[84,] 2.028362611 -4.538088606
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 2.869673050 2.161918705
2 2.368288831 2.869673050
3 -2.515168689 2.368288831
4 1.138838197 -2.515168689
5 -13.485371021 1.138838197
6 -11.198973243 -13.485371021
7 4.280575966 -11.198973243
8 2.930137668 4.280575966
9 -0.528822026 2.930137668
10 -4.848158731 -0.528822026
11 0.874929089 -4.848158731
12 1.511487720 0.874929089
13 0.767218548 1.511487720
14 2.747063675 0.767218548
15 0.846157901 2.747063675
16 -7.000189559 0.846157901
17 -10.223573298 -7.000189559
18 3.568811079 -10.223573298
19 1.889827923 3.568811079
20 -3.160778277 1.889827923
21 -8.214741359 -3.160778277
22 -2.651431401 -8.214741359
23 2.642929034 -2.651431401
24 2.401070561 2.642929034
25 -1.096953382 2.401070561
26 0.987676651 -1.096953382
27 0.408817594 0.987676651
28 1.452703240 0.408817594
29 1.624019659 1.452703240
30 2.980959492 1.624019659
31 2.971908535 2.980959492
32 2.255756892 2.971908535
33 3.415534754 2.255756892
34 5.073117540 3.415534754
35 5.183969340 5.073117540
36 -1.315515723 5.183969340
37 4.056322944 -1.315515723
38 2.656366747 4.056322944
39 -1.799451687 2.656366747
40 -44.604185478 -1.799451687
41 -0.006256454 -44.604185478
42 2.603430145 -0.006256454
43 1.621313041 2.603430145
44 -0.424252744 1.621313041
45 3.084440473 -0.424252744
46 2.694561713 3.084440473
47 -1.016768532 2.694561713
48 4.200982954 -1.016768532
49 2.502655781 4.200982954
50 6.121608960 2.502655781
51 5.050632264 6.121608960
52 -2.913010106 5.050632264
53 2.468036715 -2.913010106
54 0.977164574 2.468036715
55 0.078157955 0.977164574
56 -0.162239397 0.078157955
57 1.204124213 -0.162239397
58 2.948393385 1.204124213
59 -3.028526107 2.948393385
60 4.252831257 -3.028526107
61 3.612262084 4.252831257
62 -0.108276315 3.612262084
63 -3.156983359 -0.108276315
64 -0.871454863 -3.156983359
65 3.350683379 -0.871454863
66 -2.340491993 3.350683379
67 3.236847516 -2.340491993
68 -2.760168421 3.236847516
69 -2.515603330 -2.760168421
70 7.883326387 -2.515603330
71 -8.836620174 7.883326387
72 4.048517484 -8.836620174
73 2.209975484 4.048517484
74 1.882022463 2.209975484
75 3.977975429 1.882022463
76 -0.173186117 3.977975429
77 0.774746563 -0.173186117
78 1.645420028 0.774746563
79 4.149291848 1.645420028
80 2.199679030 4.149291848
81 0.373110228 2.199679030
82 -1.751394909 0.373110228
83 -4.538088606 -1.751394909
84 2.028362611 -4.538088606
> plot(z,main=paste('Residual Lag plot, lowess, and regression line'), ylab='values of Residuals', xlab='lagged values of Residuals')
> lines(lowess(z))
> abline(lm(z))
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/7qinw1324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/81vvl1324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/wessaorg/rcomp/tmp/9i4sm1324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
> plot(mylm, las = 1, sub='Residual Diagnostics')
> par(opar)
> dev.off()
null device
1
> if (n > n25) {
+ postscript(file="/var/wessaorg/rcomp/tmp/10771x1324132041.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
+ plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint')
+ grid()
+ dev.off()
+ }
null device
1
>
> #Note: the /var/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/wessaorg/rcomp/createtable")
>
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Estimated Regression Equation', 1, TRUE)
> a<-table.row.end(a)
> myeq <- colnames(x)[1]
> myeq <- paste(myeq, '[t] = ', sep='')
> for (i in 1:k){
+ if (mysum$coefficients[i,1] > 0) myeq <- paste(myeq, '+', '')
+ myeq <- paste(myeq, mysum$coefficients[i,1], sep=' ')
+ if (rownames(mysum$coefficients)[i] != '(Intercept)') {
+ myeq <- paste(myeq, rownames(mysum$coefficients)[i], sep='')
+ if (rownames(mysum$coefficients)[i] != 't') myeq <- paste(myeq, '[t]', sep='')
+ }
+ }
> myeq <- paste(myeq, ' + e[t]')
> a<-table.row.start(a)
> a<-table.element(a, myeq)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/11mggv1324132041.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,hyperlink('http://www.xycoon.com/ols1.htm','Multiple Linear Regression - Ordinary Least Squares',''), 6, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Variable',header=TRUE)
> a<-table.element(a,'Parameter',header=TRUE)
> a<-table.element(a,'S.D.',header=TRUE)
> a<-table.element(a,'T-STAT
H0: parameter = 0',header=TRUE)
> a<-table.element(a,'2-tail p-value',header=TRUE)
> a<-table.element(a,'1-tail p-value',header=TRUE)
> a<-table.row.end(a)
> for (i in 1:k){
+ a<-table.row.start(a)
+ a<-table.element(a,rownames(mysum$coefficients)[i],header=TRUE)
+ a<-table.element(a,mysum$coefficients[i,1])
+ a<-table.element(a, round(mysum$coefficients[i,2],6))
+ a<-table.element(a, round(mysum$coefficients[i,3],4))
+ a<-table.element(a, round(mysum$coefficients[i,4],6))
+ a<-table.element(a, round(mysum$coefficients[i,4]/2,6))
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/12xo6m1324132041.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Regression Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple R',1,TRUE)
> a<-table.element(a, sqrt(mysum$r.squared))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'R-squared',1,TRUE)
> a<-table.element(a, mysum$r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Adjusted R-squared',1,TRUE)
> a<-table.element(a, mysum$adj.r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (value)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[3])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'p-value',1,TRUE)
> a<-table.element(a, 1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3]))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Residual Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Residual Standard Deviation',1,TRUE)
> a<-table.element(a, mysum$sigma)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Sum Squared Residuals',1,TRUE)
> a<-table.element(a, sum(myerror*myerror))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/13s9r91324132041.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Actuals, Interpolation, and Residuals', 4, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Time or Index', 1, TRUE)
> a<-table.element(a, 'Actuals', 1, TRUE)
> a<-table.element(a, 'Interpolation
Forecast', 1, TRUE)
> a<-table.element(a, 'Residuals
Prediction Error', 1, TRUE)
> a<-table.row.end(a)
> for (i in 1:n) {
+ a<-table.row.start(a)
+ a<-table.element(a,i, 1, TRUE)
+ a<-table.element(a,x[i])
+ a<-table.element(a,x[i]-mysum$resid[i])
+ a<-table.element(a,mysum$resid[i])
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/1403eq1324132041.tab")
> if (n > n25) {
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'p-values',header=TRUE)
+ a<-table.element(a,'Alternative Hypothesis',3,header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'breakpoint index',header=TRUE)
+ a<-table.element(a,'greater',header=TRUE)
+ a<-table.element(a,'2-sided',header=TRUE)
+ a<-table.element(a,'less',header=TRUE)
+ a<-table.row.end(a)
+ for (mypoint in kp3:nmkm3) {
+ a<-table.row.start(a)
+ a<-table.element(a,mypoint,header=TRUE)
+ a<-table.element(a,gqarr[mypoint-kp3+1,1])
+ a<-table.element(a,gqarr[mypoint-kp3+1,2])
+ a<-table.element(a,gqarr[mypoint-kp3+1,3])
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/15m1mn1324132041.tab")
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Meta Analysis of Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'Description',header=TRUE)
+ a<-table.element(a,'# significant tests',header=TRUE)
+ a<-table.element(a,'% significant tests',header=TRUE)
+ a<-table.element(a,'OK/NOK',header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'1% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant1)
+ a<-table.element(a,numsignificant1/numgqtests)
+ if (numsignificant1/numgqtests < 0.01) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'5% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant5)
+ a<-table.element(a,numsignificant5/numgqtests)
+ if (numsignificant5/numgqtests < 0.05) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'10% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant10)
+ a<-table.element(a,numsignificant10/numgqtests)
+ if (numsignificant10/numgqtests < 0.1) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/16xh9g1324132041.tab")
+ }
>
> try(system("convert tmp/1jm961324132041.ps tmp/1jm961324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/243fe1324132041.ps tmp/243fe1324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/3pefi1324132041.ps tmp/3pefi1324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/49h8a1324132041.ps tmp/49h8a1324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/5o3341324132041.ps tmp/5o3341324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/6hqc41324132041.ps tmp/6hqc41324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/7qinw1324132041.ps tmp/7qinw1324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/81vvl1324132041.ps tmp/81vvl1324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/9i4sm1324132041.ps tmp/9i4sm1324132041.png",intern=TRUE))
character(0)
> try(system("convert tmp/10771x1324132041.ps tmp/10771x1324132041.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.348 0.687 4.159