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(210907
+ ,79
+ ,94
+ ,112285
+ ,-1
+ ,179321
+ ,108
+ ,103
+ ,101193
+ ,3
+ ,149061
+ ,43
+ ,93
+ ,116174
+ ,0
+ ,237213
+ ,78
+ ,123
+ ,66198
+ ,3
+ ,173326
+ ,86
+ ,148
+ ,71701
+ ,4
+ ,133131
+ ,44
+ ,90
+ ,57793
+ ,0
+ ,258873
+ ,104
+ ,124
+ ,80444
+ ,0
+ ,324799
+ ,158
+ ,168
+ ,97668
+ ,7
+ ,230964
+ ,102
+ ,115
+ ,133824
+ ,1
+ ,236785
+ ,77
+ ,71
+ ,101481
+ ,0
+ ,344297
+ ,80
+ ,108
+ ,67654
+ ,1
+ ,174724
+ ,123
+ ,120
+ ,69112
+ ,4
+ ,174415
+ ,73
+ ,114
+ ,82753
+ ,1
+ ,223632
+ ,105
+ ,120
+ ,72654
+ ,5
+ ,294424
+ ,107
+ ,124
+ ,101494
+ ,13
+ ,325107
+ ,84
+ ,126
+ ,79215
+ ,4
+ ,106408
+ ,33
+ ,37
+ ,31081
+ ,0
+ ,96560
+ ,42
+ ,38
+ ,22996
+ ,0
+ ,265769
+ ,96
+ ,120
+ ,83122
+ ,6
+ ,269651
+ ,106
+ ,93
+ ,70106
+ ,0
+ ,149112
+ ,56
+ ,95
+ ,60578
+ ,1
+ ,152871
+ ,59
+ ,90
+ ,79892
+ ,3
+ ,362301
+ ,76
+ ,110
+ ,100708
+ ,1
+ ,183167
+ ,91
+ ,138
+ ,82875
+ ,0
+ ,277965
+ ,115
+ ,133
+ ,139077
+ ,2
+ ,218946
+ ,76
+ ,96
+ ,80670
+ ,3
+ ,244052
+ ,101
+ ,164
+ ,143558
+ ,4
+ ,341570
+ ,94
+ ,78
+ ,117105
+ ,12
+ ,233328
+ ,92
+ ,102
+ ,120733
+ ,0
+ ,206161
+ ,75
+ ,99
+ ,73107
+ ,3
+ ,311473
+ ,128
+ ,129
+ ,132068
+ ,0
+ ,207176
+ ,56
+ ,114
+ ,87011
+ ,4
+ ,196553
+ ,41
+ ,99
+ ,95260
+ ,-1
+ ,143246
+ ,67
+ ,104
+ ,106671
+ ,2
+ ,182192
+ ,77
+ ,138
+ ,70054
+ ,1
+ ,194979
+ ,66
+ ,151
+ ,74011
+ ,1
+ ,167488
+ ,69
+ ,72
+ ,83737
+ ,0
+ ,143756
+ ,105
+ ,120
+ ,69094
+ ,2
+ ,275541
+ ,116
+ ,115
+ ,93133
+ ,0
+ ,152299
+ ,62
+ ,98
+ ,61370
+ ,2
+ ,193339
+ ,100
+ ,71
+ ,84651
+ ,4
+ ,130585
+ ,67
+ ,107
+ ,95364
+ ,0
+ ,112611
+ ,46
+ ,73
+ ,26706
+ ,0
+ ,148446
+ ,135
+ ,129
+ ,126846
+ ,6
+ ,182079
+ ,124
+ ,118
+ ,102860
+ ,13
+ ,243060
+ ,58
+ ,104
+ ,111813
+ ,4
+ ,162765
+ ,68
+ ,107
+ ,120293
+ ,-1
+ ,85574
+ ,37
+ ,36
+ ,24266
+ ,3
+ ,225060
+ ,93
+ ,139
+ ,109825
+ ,0
+ ,133328
+ ,56
+ ,56
+ ,40909
+ ,2
+ ,100750
+ ,83
+ ,93
+ ,140867
+ ,0
+ ,101523
+ ,59
+ ,87
+ ,61056
+ ,1
+ ,243511
+ ,133
+ ,110
+ ,101338
+ ,1
+ ,152474
+ ,106
+ ,83
+ ,65567
+ ,0
+ ,132487
+ ,71
+ ,98
+ ,40735
+ ,31
+ ,317394
+ ,116
+ ,82
+ ,91413
+ ,2
+ ,244749
+ ,98
+ ,115
+ ,76643
+ ,5
+ ,184510
+ ,64
+ ,140
+ ,110681
+ ,1
+ ,128423
+ ,32
+ ,120
+ ,92696
+ ,1
+ ,97839
+ ,25
+ ,66
+ ,94785
+ ,2
+ ,172494
+ ,46
+ ,139
+ ,86687
+ ,13
+ ,229242
+ ,63
+ ,119
+ ,91721
+ ,5
+ ,351619
+ ,95
+ ,141
+ ,115168
+ ,3
+ ,324598
+ ,113
+ ,133
+ ,135777
+ ,1
+ ,195838
+ ,111
+ ,98
+ ,102372
+ ,1
+ ,254488
+ ,120
+ ,117
+ ,103772
+ ,4
+ ,199476
+ ,87
+ ,105
+ ,135400
+ ,2
+ ,92499
+ ,25
+ ,55
+ ,21399
+ ,0
+ ,224330
+ ,131
+ ,132
+ ,130115
+ ,4
+ ,181633
+ ,47
+ ,73
+ ,64466
+ ,0
+ ,271856
+ ,109
+ ,86
+ ,54990
+ ,0
+ ,95227
+ ,37
+ ,48
+ ,34777
+ ,0
+ ,98146
+ ,15
+ ,48
+ ,27114
+ ,7
+ ,118612
+ ,54
+ ,43
+ ,30080
+ ,3
+ ,65475
+ ,16
+ ,46
+ ,69008
+ ,4
+ ,108446
+ ,22
+ ,65
+ ,46300
+ ,1
+ ,121848
+ ,37
+ ,52
+ ,30594
+ ,0
+ ,76302
+ ,29
+ ,68
+ ,30976
+ ,2
+ ,98104
+ ,55
+ ,47
+ ,25568
+ ,0
+ ,30989
+ ,5
+ ,41
+ ,4154
+ ,0
+ ,31774
+ ,0
+ ,47
+ ,4143
+ ,0
+ ,150580
+ ,27
+ ,71
+ ,45588
+ ,2
+ ,54157
+ ,37
+ ,30
+ ,18625
+ ,1
+ ,59382
+ ,29
+ ,24
+ ,26263
+ ,0
+ ,84105
+ ,17
+ ,63
+ ,20055
+ ,0)
+ ,dim=c(5
+ ,85)
+ ,dimnames=list(c('time_in_rfc'
+ ,'blogged_computations'
+ ,'feedback_messages_p120'
+ ,'totsize'
+ ,'difference_hyperlinks-blogs')
+ ,1:85))
> y <- array(NA,dim=c(5,85),dimnames=list(c('time_in_rfc','blogged_computations','feedback_messages_p120','totsize','difference_hyperlinks-blogs'),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 = '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
> 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
time_in_rfc blogged_computations feedback_messages_p120 totsize
1 210907 79 94 112285
2 179321 108 103 101193
3 149061 43 93 116174
4 237213 78 123 66198
5 173326 86 148 71701
6 133131 44 90 57793
7 258873 104 124 80444
8 324799 158 168 97668
9 230964 102 115 133824
10 236785 77 71 101481
11 344297 80 108 67654
12 174724 123 120 69112
13 174415 73 114 82753
14 223632 105 120 72654
15 294424 107 124 101494
16 325107 84 126 79215
17 106408 33 37 31081
18 96560 42 38 22996
19 265769 96 120 83122
20 269651 106 93 70106
21 149112 56 95 60578
22 152871 59 90 79892
23 362301 76 110 100708
24 183167 91 138 82875
25 277965 115 133 139077
26 218946 76 96 80670
27 244052 101 164 143558
28 341570 94 78 117105
29 233328 92 102 120733
30 206161 75 99 73107
31 311473 128 129 132068
32 207176 56 114 87011
33 196553 41 99 95260
34 143246 67 104 106671
35 182192 77 138 70054
36 194979 66 151 74011
37 167488 69 72 83737
38 143756 105 120 69094
39 275541 116 115 93133
40 152299 62 98 61370
41 193339 100 71 84651
42 130585 67 107 95364
43 112611 46 73 26706
44 148446 135 129 126846
45 182079 124 118 102860
46 243060 58 104 111813
47 162765 68 107 120293
48 85574 37 36 24266
49 225060 93 139 109825
50 133328 56 56 40909
51 100750 83 93 140867
52 101523 59 87 61056
53 243511 133 110 101338
54 152474 106 83 65567
55 132487 71 98 40735
56 317394 116 82 91413
57 244749 98 115 76643
58 184510 64 140 110681
59 128423 32 120 92696
60 97839 25 66 94785
61 172494 46 139 86687
62 229242 63 119 91721
63 351619 95 141 115168
64 324598 113 133 135777
65 195838 111 98 102372
66 254488 120 117 103772
67 199476 87 105 135400
68 92499 25 55 21399
69 224330 131 132 130115
70 181633 47 73 64466
71 271856 109 86 54990
72 95227 37 48 34777
73 98146 15 48 27114
74 118612 54 43 30080
75 65475 16 46 69008
76 108446 22 65 46300
77 121848 37 52 30594
78 76302 29 68 30976
79 98104 55 47 25568
80 30989 5 41 4154
81 31774 0 47 4143
82 150580 27 71 45588
83 54157 37 30 18625
84 59382 29 24 26263
85 84105 17 63 20055
difference_hyperlinks-blogs
1 -1
2 3
3 0
4 3
5 4
6 0
7 0
8 7
9 1
10 0
11 1
12 4
13 1
14 5
15 13
16 4
17 0
18 0
19 6
20 0
21 1
22 3
23 1
24 0
25 2
26 3
27 4
28 12
29 0
30 3
31 0
32 4
33 -1
34 2
35 1
36 1
37 0
38 2
39 0
40 2
41 4
42 0
43 0
44 6
45 13
46 4
47 -1
48 3
49 0
50 2
51 0
52 1
53 1
54 0
55 31
56 2
57 5
58 1
59 1
60 2
61 13
62 5
63 3
64 1
65 1
66 4
67 2
68 0
69 4
70 0
71 0
72 0
73 7
74 3
75 4
76 1
77 0
78 2
79 0
80 0
81 0
82 2
83 1
84 0
85 0
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) blogged_computations
27100.7497 1057.1452
feedback_messages_p120 totsize
515.5426 0.3977
`difference_hyperlinks-blogs`
-93.5784
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-137758 -27388 -8570 27874 158191
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 27100.7497 17564.2879 1.543 0.1268
blogged_computations 1057.1452 240.7064 4.392 3.41e-05 ***
feedback_messages_p120 515.5426 268.5699 1.920 0.0585 .
totsize 0.3977 0.2460 1.617 0.1099
`difference_hyperlinks-blogs` -93.5784 1371.1381 -0.068 0.9458
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 52660 on 80 degrees of freedom
Multiple R-squared: 0.5859, Adjusted R-squared: 0.5652
F-statistic: 28.3 on 4 and 80 DF, p-value: 1.178e-14
> 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.4004091 0.8008181769 0.5995909116
[2,] 0.2403042 0.4806084111 0.7596957944
[3,] 0.2639649 0.5279297016 0.7360351492
[4,] 0.7609402 0.4781195951 0.2390597975
[5,] 0.8629266 0.2741468393 0.1370734197
[6,] 0.8195227 0.3609545710 0.1804772855
[7,] 0.7508784 0.4982431292 0.2491215646
[8,] 0.7979929 0.4040142449 0.2020071224
[9,] 0.9071211 0.1857577979 0.0928788989
[10,] 0.8761193 0.2477614614 0.1238807307
[11,] 0.8343759 0.3312481821 0.1656240910
[12,] 0.7955076 0.4089847223 0.2044923612
[13,] 0.7952432 0.4095135084 0.2047567542
[14,] 0.7470604 0.5058792703 0.2529396351
[15,] 0.6970935 0.6058130492 0.3029065246
[16,] 0.9528004 0.0943991283 0.0471995641
[17,] 0.9493429 0.1013142342 0.0506571171
[18,] 0.9284837 0.1430326899 0.0715163449
[19,] 0.9064870 0.1870259722 0.0935129861
[20,] 0.8884306 0.2231388220 0.1115694110
[21,] 0.9574133 0.0851734906 0.0425867453
[22,] 0.9412090 0.1175820210 0.0587910105
[23,] 0.9216206 0.1567587646 0.0783793823
[24,] 0.9088430 0.1823139086 0.0911569543
[25,] 0.8859502 0.2280996961 0.1140498480
[26,] 0.8695774 0.2608452251 0.1304226126
[27,] 0.8792481 0.2415038323 0.1207519161
[28,] 0.8495229 0.3009541573 0.1504770787
[29,] 0.8101986 0.3796028237 0.1898014119
[30,] 0.7708569 0.4582862983 0.2291431491
[31,] 0.8528475 0.2943049248 0.1471524624
[32,] 0.8288652 0.3422695771 0.1711347885
[33,] 0.7921865 0.4156269694 0.2078134847
[34,] 0.7644256 0.4711488988 0.2355744494
[35,] 0.7819932 0.4360135291 0.2180067646
[36,] 0.7383843 0.5232313345 0.2616156673
[37,] 0.9478925 0.1042149555 0.0521074778
[38,] 0.9686342 0.0627315721 0.0313657861
[39,] 0.9770730 0.0458540063 0.0229270031
[40,] 0.9709190 0.0581619151 0.0290809576
[41,] 0.9601482 0.0797036504 0.0398518252
[42,] 0.9479250 0.1041499567 0.0520749783
[43,] 0.9278515 0.1442970334 0.0721485167
[44,] 0.9769436 0.0461128340 0.0230564170
[45,] 0.9817883 0.0364234735 0.0182117368
[46,] 0.9771030 0.0457939409 0.0228969705
[47,] 0.9862138 0.0275723764 0.0137861882
[48,] 0.9832287 0.0335425592 0.0167712796
[49,] 0.9940115 0.0119769541 0.0059884770
[50,] 0.9901771 0.0196458160 0.0098229080
[51,] 0.9898595 0.0202809173 0.0101404587
[52,] 0.9927487 0.0145026162 0.0072513081
[53,] 0.9881812 0.0236375918 0.0118187959
[54,] 0.9864079 0.0271841956 0.0135920978
[55,] 0.9783065 0.0433869845 0.0216934923
[56,] 0.9923486 0.0153027441 0.0076513721
[57,] 0.9955370 0.0089260787 0.0044630393
[58,] 0.9936005 0.0127989151 0.0063994576
[59,] 0.9877935 0.0244130816 0.0122065408
[60,] 0.9780719 0.0438561940 0.0219280970
[61,] 0.9612082 0.0775836222 0.0387918111
[62,] 0.9995501 0.0008998836 0.0004499418
[63,] 0.9991720 0.0016560386 0.0008280193
[64,] 0.9981808 0.0036384983 0.0018192491
[65,] 0.9950221 0.0099557619 0.0049778809
[66,] 0.9939530 0.0120939654 0.0060469827
[67,] 0.9943715 0.0112570689 0.0056285345
[68,] 0.9886029 0.0227942775 0.0113971387
[69,] 0.9873136 0.0253727539 0.0126863770
[70,] 0.9601586 0.0796827897 0.0398413948
> postscript(file="/var/wessaorg/rcomp/tmp/1nn9r1324326746.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/2a0y61324326746.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/3okv51324326746.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/4olzx1324326746.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/5acbk1324326746.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 6
7082.8704 -55014.7626 -17643.3876 38197.8084 -49129.8144 -9866.5207
7 8 9 10 11 12
25910.3019 5871.8668 -16379.4954 51322.8341 150134.4683 -71381.3799
13 14 15 16 17 18
-21445.4523 -4759.7966 51135.0667 113119.1916 12985.8609 -3676.6864
19 20 21 22 23 24
42822.0975 54667.1253 -10162.9571 -14491.4559 158190.8066 -44237.1856
25 26 27 28 29 30
5603.3762 30209.4192 -31086.2934 129437.0407 8370.5575 19942.6458
31 32 33 34 35 36
30030.8931 27874.3984 37093.3014 -50534.4620 -25219.8252 -9079.9313
37 38 39 40 41 42
-2975.9830 -83500.7647 29486.1889 -15086.8462 -9370.1141 -60432.5960
43 44 45 46 47 48
-11373.6799 -137757.9298 -76631.3866 56936.0932 -39317.2670 -8570.2030
49 50 51 52 53 54
-15691.6923 2074.8967 -118060.2887 -56989.1460 -21107.0124 -55549.3454
55 56 57 58 59 60
-33493.1023 89223.2752 24748.5583 -26346.8663 -31141.9605 -27223.8315
61 62 63 64 65 66
-8153.6600 38183.1398 105877.8237 55569.4565 -39747.5151 -683.1593
67 68 69 70 71 72
-27388.0679 2104.6635 -60679.1736 41574.5012 63320.9300 -9564.5414
73 74 75 76 77 78
20314.1794 575.3761 -29324.3287 6258.4387 16657.8141 -28644.4666
79 80 81 82 83 84
-21538.3096 -24186.7175 -21204.8729 40390.1886 -34837.7495 -21193.4463
85
-1422.0244
> postscript(file="/var/wessaorg/rcomp/tmp/6c6f71324326746.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 7082.8704 NA
1 -55014.7626 7082.8704
2 -17643.3876 -55014.7626
3 38197.8084 -17643.3876
4 -49129.8144 38197.8084
5 -9866.5207 -49129.8144
6 25910.3019 -9866.5207
7 5871.8668 25910.3019
8 -16379.4954 5871.8668
9 51322.8341 -16379.4954
10 150134.4683 51322.8341
11 -71381.3799 150134.4683
12 -21445.4523 -71381.3799
13 -4759.7966 -21445.4523
14 51135.0667 -4759.7966
15 113119.1916 51135.0667
16 12985.8609 113119.1916
17 -3676.6864 12985.8609
18 42822.0975 -3676.6864
19 54667.1253 42822.0975
20 -10162.9571 54667.1253
21 -14491.4559 -10162.9571
22 158190.8066 -14491.4559
23 -44237.1856 158190.8066
24 5603.3762 -44237.1856
25 30209.4192 5603.3762
26 -31086.2934 30209.4192
27 129437.0407 -31086.2934
28 8370.5575 129437.0407
29 19942.6458 8370.5575
30 30030.8931 19942.6458
31 27874.3984 30030.8931
32 37093.3014 27874.3984
33 -50534.4620 37093.3014
34 -25219.8252 -50534.4620
35 -9079.9313 -25219.8252
36 -2975.9830 -9079.9313
37 -83500.7647 -2975.9830
38 29486.1889 -83500.7647
39 -15086.8462 29486.1889
40 -9370.1141 -15086.8462
41 -60432.5960 -9370.1141
42 -11373.6799 -60432.5960
43 -137757.9298 -11373.6799
44 -76631.3866 -137757.9298
45 56936.0932 -76631.3866
46 -39317.2670 56936.0932
47 -8570.2030 -39317.2670
48 -15691.6923 -8570.2030
49 2074.8967 -15691.6923
50 -118060.2887 2074.8967
51 -56989.1460 -118060.2887
52 -21107.0124 -56989.1460
53 -55549.3454 -21107.0124
54 -33493.1023 -55549.3454
55 89223.2752 -33493.1023
56 24748.5583 89223.2752
57 -26346.8663 24748.5583
58 -31141.9605 -26346.8663
59 -27223.8315 -31141.9605
60 -8153.6600 -27223.8315
61 38183.1398 -8153.6600
62 105877.8237 38183.1398
63 55569.4565 105877.8237
64 -39747.5151 55569.4565
65 -683.1593 -39747.5151
66 -27388.0679 -683.1593
67 2104.6635 -27388.0679
68 -60679.1736 2104.6635
69 41574.5012 -60679.1736
70 63320.9300 41574.5012
71 -9564.5414 63320.9300
72 20314.1794 -9564.5414
73 575.3761 20314.1794
74 -29324.3287 575.3761
75 6258.4387 -29324.3287
76 16657.8141 6258.4387
77 -28644.4666 16657.8141
78 -21538.3096 -28644.4666
79 -24186.7175 -21538.3096
80 -21204.8729 -24186.7175
81 40390.1886 -21204.8729
82 -34837.7495 40390.1886
83 -21193.4463 -34837.7495
84 -1422.0244 -21193.4463
85 NA -1422.0244
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -55014.7626 7082.8704
[2,] -17643.3876 -55014.7626
[3,] 38197.8084 -17643.3876
[4,] -49129.8144 38197.8084
[5,] -9866.5207 -49129.8144
[6,] 25910.3019 -9866.5207
[7,] 5871.8668 25910.3019
[8,] -16379.4954 5871.8668
[9,] 51322.8341 -16379.4954
[10,] 150134.4683 51322.8341
[11,] -71381.3799 150134.4683
[12,] -21445.4523 -71381.3799
[13,] -4759.7966 -21445.4523
[14,] 51135.0667 -4759.7966
[15,] 113119.1916 51135.0667
[16,] 12985.8609 113119.1916
[17,] -3676.6864 12985.8609
[18,] 42822.0975 -3676.6864
[19,] 54667.1253 42822.0975
[20,] -10162.9571 54667.1253
[21,] -14491.4559 -10162.9571
[22,] 158190.8066 -14491.4559
[23,] -44237.1856 158190.8066
[24,] 5603.3762 -44237.1856
[25,] 30209.4192 5603.3762
[26,] -31086.2934 30209.4192
[27,] 129437.0407 -31086.2934
[28,] 8370.5575 129437.0407
[29,] 19942.6458 8370.5575
[30,] 30030.8931 19942.6458
[31,] 27874.3984 30030.8931
[32,] 37093.3014 27874.3984
[33,] -50534.4620 37093.3014
[34,] -25219.8252 -50534.4620
[35,] -9079.9313 -25219.8252
[36,] -2975.9830 -9079.9313
[37,] -83500.7647 -2975.9830
[38,] 29486.1889 -83500.7647
[39,] -15086.8462 29486.1889
[40,] -9370.1141 -15086.8462
[41,] -60432.5960 -9370.1141
[42,] -11373.6799 -60432.5960
[43,] -137757.9298 -11373.6799
[44,] -76631.3866 -137757.9298
[45,] 56936.0932 -76631.3866
[46,] -39317.2670 56936.0932
[47,] -8570.2030 -39317.2670
[48,] -15691.6923 -8570.2030
[49,] 2074.8967 -15691.6923
[50,] -118060.2887 2074.8967
[51,] -56989.1460 -118060.2887
[52,] -21107.0124 -56989.1460
[53,] -55549.3454 -21107.0124
[54,] -33493.1023 -55549.3454
[55,] 89223.2752 -33493.1023
[56,] 24748.5583 89223.2752
[57,] -26346.8663 24748.5583
[58,] -31141.9605 -26346.8663
[59,] -27223.8315 -31141.9605
[60,] -8153.6600 -27223.8315
[61,] 38183.1398 -8153.6600
[62,] 105877.8237 38183.1398
[63,] 55569.4565 105877.8237
[64,] -39747.5151 55569.4565
[65,] -683.1593 -39747.5151
[66,] -27388.0679 -683.1593
[67,] 2104.6635 -27388.0679
[68,] -60679.1736 2104.6635
[69,] 41574.5012 -60679.1736
[70,] 63320.9300 41574.5012
[71,] -9564.5414 63320.9300
[72,] 20314.1794 -9564.5414
[73,] 575.3761 20314.1794
[74,] -29324.3287 575.3761
[75,] 6258.4387 -29324.3287
[76,] 16657.8141 6258.4387
[77,] -28644.4666 16657.8141
[78,] -21538.3096 -28644.4666
[79,] -24186.7175 -21538.3096
[80,] -21204.8729 -24186.7175
[81,] 40390.1886 -21204.8729
[82,] -34837.7495 40390.1886
[83,] -21193.4463 -34837.7495
[84,] -1422.0244 -21193.4463
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -55014.7626 7082.8704
2 -17643.3876 -55014.7626
3 38197.8084 -17643.3876
4 -49129.8144 38197.8084
5 -9866.5207 -49129.8144
6 25910.3019 -9866.5207
7 5871.8668 25910.3019
8 -16379.4954 5871.8668
9 51322.8341 -16379.4954
10 150134.4683 51322.8341
11 -71381.3799 150134.4683
12 -21445.4523 -71381.3799
13 -4759.7966 -21445.4523
14 51135.0667 -4759.7966
15 113119.1916 51135.0667
16 12985.8609 113119.1916
17 -3676.6864 12985.8609
18 42822.0975 -3676.6864
19 54667.1253 42822.0975
20 -10162.9571 54667.1253
21 -14491.4559 -10162.9571
22 158190.8066 -14491.4559
23 -44237.1856 158190.8066
24 5603.3762 -44237.1856
25 30209.4192 5603.3762
26 -31086.2934 30209.4192
27 129437.0407 -31086.2934
28 8370.5575 129437.0407
29 19942.6458 8370.5575
30 30030.8931 19942.6458
31 27874.3984 30030.8931
32 37093.3014 27874.3984
33 -50534.4620 37093.3014
34 -25219.8252 -50534.4620
35 -9079.9313 -25219.8252
36 -2975.9830 -9079.9313
37 -83500.7647 -2975.9830
38 29486.1889 -83500.7647
39 -15086.8462 29486.1889
40 -9370.1141 -15086.8462
41 -60432.5960 -9370.1141
42 -11373.6799 -60432.5960
43 -137757.9298 -11373.6799
44 -76631.3866 -137757.9298
45 56936.0932 -76631.3866
46 -39317.2670 56936.0932
47 -8570.2030 -39317.2670
48 -15691.6923 -8570.2030
49 2074.8967 -15691.6923
50 -118060.2887 2074.8967
51 -56989.1460 -118060.2887
52 -21107.0124 -56989.1460
53 -55549.3454 -21107.0124
54 -33493.1023 -55549.3454
55 89223.2752 -33493.1023
56 24748.5583 89223.2752
57 -26346.8663 24748.5583
58 -31141.9605 -26346.8663
59 -27223.8315 -31141.9605
60 -8153.6600 -27223.8315
61 38183.1398 -8153.6600
62 105877.8237 38183.1398
63 55569.4565 105877.8237
64 -39747.5151 55569.4565
65 -683.1593 -39747.5151
66 -27388.0679 -683.1593
67 2104.6635 -27388.0679
68 -60679.1736 2104.6635
69 41574.5012 -60679.1736
70 63320.9300 41574.5012
71 -9564.5414 63320.9300
72 20314.1794 -9564.5414
73 575.3761 20314.1794
74 -29324.3287 575.3761
75 6258.4387 -29324.3287
76 16657.8141 6258.4387
77 -28644.4666 16657.8141
78 -21538.3096 -28644.4666
79 -24186.7175 -21538.3096
80 -21204.8729 -24186.7175
81 40390.1886 -21204.8729
82 -34837.7495 40390.1886
83 -21193.4463 -34837.7495
84 -1422.0244 -21193.4463
> 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/7vw231324326746.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/8x7pl1324326746.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/9vzd61324326746.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/10zbco1324326746.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/11w8uk1324326746.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/12tuh41324326746.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/136toi1324326746.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/14z3qv1324326746.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/15ztyx1324326746.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/161ojx1324326746.tab")
+ }
>
> try(system("convert tmp/1nn9r1324326746.ps tmp/1nn9r1324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/2a0y61324326746.ps tmp/2a0y61324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/3okv51324326746.ps tmp/3okv51324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/4olzx1324326746.ps tmp/4olzx1324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/5acbk1324326746.ps tmp/5acbk1324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/6c6f71324326746.ps tmp/6c6f71324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/7vw231324326746.ps tmp/7vw231324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/8x7pl1324326746.ps tmp/8x7pl1324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/9vzd61324326746.ps tmp/9vzd61324326746.png",intern=TRUE))
character(0)
> try(system("convert tmp/10zbco1324326746.ps tmp/10zbco1324326746.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
3.705 0.721 4.445