R version 2.8.0 (2008-10-20)
Copyright (C) 2008 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
Natural language support but running in an English locale
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> x <- array(list(621,0,587,0,655,0,517,0,646,0,657,0,382,0,345,0,625,0,654,0,606,0,510,0,614,0,647,0,580,0,614,0,636,0,388,0,356,0,639,0,753,0,611,0,639,0,630,0,586,0,695,0,552,0,619,0,681,0,421,0,307,0,754,0,690,0,644,0,643,0,608,0,651,0,691,0,627,0,634,0,731,0,475,0,337,0,803,0,722,0,590,0,724,0,627,0,696,0,825,0,677,0,656,0,785,0,412,0,352,0,839,0,729,0,696,0,641,0,695,0,638,0,762,0,635,0,721,0,854,0,418,0,367,0,824,0,687,0,601,0,676,0,740,0,691,0,683,0,594,0,729,0,731,0,386,0,331,0,706,0,715,0,657,0,653,0,642,0,643,0,718,0,654,0,632,0,731,0,392,0,344,0,792,0,852,0,649,0,629,0,685,0,617,0,715,0,715,0,629,0,916,0,531,1,357,1,917,1,828,1,708,1,858,1,775,1,785,1,1006,1,789,1,734,1,906,1,532,1,387,1,991,1,841,1,892,1,782,1,813,1,793,1,978,1,775,1,797,1,946,1,594,1,438,1,1022,1,868,1,795,1),dim=c(2,130),dimnames=list(c('Y','X1'),1:130))
> y <- array(NA,dim=c(2,130),dimnames=list(c('Y','X1'),1:130))
> 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 X1
1 621 0
2 587 0
3 655 0
4 517 0
5 646 0
6 657 0
7 382 0
8 345 0
9 625 0
10 654 0
11 606 0
12 510 0
13 614 0
14 647 0
15 580 0
16 614 0
17 636 0
18 388 0
19 356 0
20 639 0
21 753 0
22 611 0
23 639 0
24 630 0
25 586 0
26 695 0
27 552 0
28 619 0
29 681 0
30 421 0
31 307 0
32 754 0
33 690 0
34 644 0
35 643 0
36 608 0
37 651 0
38 691 0
39 627 0
40 634 0
41 731 0
42 475 0
43 337 0
44 803 0
45 722 0
46 590 0
47 724 0
48 627 0
49 696 0
50 825 0
51 677 0
52 656 0
53 785 0
54 412 0
55 352 0
56 839 0
57 729 0
58 696 0
59 641 0
60 695 0
61 638 0
62 762 0
63 635 0
64 721 0
65 854 0
66 418 0
67 367 0
68 824 0
69 687 0
70 601 0
71 676 0
72 740 0
73 691 0
74 683 0
75 594 0
76 729 0
77 731 0
78 386 0
79 331 0
80 706 0
81 715 0
82 657 0
83 653 0
84 642 0
85 643 0
86 718 0
87 654 0
88 632 0
89 731 0
90 392 0
91 344 0
92 792 0
93 852 0
94 649 0
95 629 0
96 685 0
97 617 0
98 715 0
99 715 0
100 629 0
101 916 0
102 531 1
103 357 1
104 917 1
105 828 1
106 708 1
107 858 1
108 775 1
109 785 1
110 1006 1
111 789 1
112 734 1
113 906 1
114 532 1
115 387 1
116 991 1
117 841 1
118 892 1
119 782 1
120 813 1
121 793 1
122 978 1
123 775 1
124 797 1
125 946 1
126 594 1
127 438 1
128 1022 1
129 868 1
130 795 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X1
628.5 145.2
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-416.72 -26.27 17.98 86.48 287.48
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 628.52 14.22 44.207 < 2e-16 ***
X1 145.20 30.10 4.823 3.94e-06 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 142.9 on 128 degrees of freedom
Multiple R-squared: 0.1538, Adjusted R-squared: 0.1472
F-statistic: 23.27 on 1 and 128 DF, p-value: 3.938e-06
> 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.101967260 0.20393452 0.8980327
[2,] 0.048190409 0.09638082 0.9518096
[3,] 0.285090923 0.57018185 0.7149091
[4,] 0.497218641 0.99443728 0.5027814
[5,] 0.404765429 0.80953086 0.5952346
[6,] 0.340355819 0.68071164 0.6596442
[7,] 0.252599572 0.50519914 0.7474004
[8,] 0.193556460 0.38711292 0.8064435
[9,] 0.138503420 0.27700684 0.8614966
[10,] 0.105325037 0.21065007 0.8946750
[11,] 0.068778569 0.13755714 0.9312214
[12,] 0.045070157 0.09014031 0.9549298
[13,] 0.030547082 0.06109416 0.9694529
[14,] 0.061083264 0.12216653 0.9389167
[15,] 0.124589447 0.24917889 0.8754106
[16,] 0.099434363 0.19886873 0.9005656
[17,] 0.132943886 0.26588777 0.8670561
[18,] 0.099467070 0.19893414 0.9005329
[19,] 0.076455716 0.15291143 0.9235443
[20,] 0.056600989 0.11320198 0.9433990
[21,] 0.039376149 0.07875230 0.9606239
[22,] 0.035462698 0.07092540 0.9645373
[23,] 0.025063389 0.05012678 0.9749366
[24,] 0.017138963 0.03427793 0.9828610
[25,] 0.014008980 0.02801796 0.9859910
[26,] 0.020223224 0.04044645 0.9797768
[27,] 0.075591534 0.15118307 0.9244085
[28,] 0.090489619 0.18097924 0.9095104
[29,] 0.080322894 0.16064579 0.9196771
[30,] 0.063108603 0.12621721 0.9368914
[31,] 0.048805622 0.09761124 0.9511944
[32,] 0.035876980 0.07175396 0.9641230
[33,] 0.027379930 0.05475986 0.9726201
[34,] 0.023134451 0.04626890 0.9768655
[35,] 0.016589004 0.03317801 0.9834110
[36,] 0.011815753 0.02363151 0.9881842
[37,] 0.011683659 0.02336732 0.9883163
[38,] 0.011853903 0.02370781 0.9881461
[39,] 0.034924156 0.06984831 0.9650758
[40,] 0.050369304 0.10073861 0.9496307
[41,] 0.046816093 0.09363219 0.9531839
[42,] 0.035480695 0.07096139 0.9645193
[43,] 0.032777624 0.06555525 0.9672224
[44,] 0.024349449 0.04869890 0.9756506
[45,] 0.020119216 0.04023843 0.9798808
[46,] 0.031301876 0.06260375 0.9686981
[47,] 0.024553057 0.04910611 0.9754469
[48,] 0.018392850 0.03678570 0.9816072
[49,] 0.021408400 0.04281680 0.9785916
[50,] 0.032195583 0.06439117 0.9678044
[51,] 0.069349928 0.13869986 0.9306501
[52,] 0.097986398 0.19597280 0.9020136
[53,] 0.089002368 0.17800474 0.9109976
[54,] 0.074915563 0.14983113 0.9250844
[55,] 0.058838439 0.11767688 0.9411616
[56,] 0.048487963 0.09697593 0.9515120
[57,] 0.037149926 0.07429985 0.9628501
[58,] 0.036478848 0.07295770 0.9635212
[59,] 0.027502966 0.05500593 0.9724970
[60,] 0.023320836 0.04664167 0.9766792
[61,] 0.036480802 0.07296160 0.9635192
[62,] 0.050826118 0.10165224 0.9491739
[63,] 0.092729806 0.18545961 0.9072702
[64,] 0.110894993 0.22178999 0.8891050
[65,] 0.091956186 0.18391237 0.9080438
[66,] 0.073875641 0.14775128 0.9261244
[67,] 0.059189570 0.11837914 0.9408104
[68,] 0.052883503 0.10576701 0.9471165
[69,] 0.042435103 0.08487021 0.9575649
[70,] 0.033312302 0.06662460 0.9666877
[71,] 0.025528490 0.05105698 0.9744715
[72,] 0.021586969 0.04317394 0.9784130
[73,] 0.018286135 0.03657227 0.9817139
[74,] 0.032716776 0.06543355 0.9672832
[75,] 0.082913574 0.16582715 0.9170864
[76,] 0.068464171 0.13692834 0.9315358
[77,] 0.056822564 0.11364513 0.9431774
[78,] 0.043859940 0.08771988 0.9561401
[79,] 0.033323427 0.06664685 0.9666766
[80,] 0.024919587 0.04983917 0.9750804
[81,] 0.018359546 0.03671909 0.9816405
[82,] 0.014453632 0.02890726 0.9855464
[83,] 0.010357402 0.02071480 0.9896426
[84,] 0.007317778 0.01463556 0.9926822
[85,] 0.005752930 0.01150586 0.9942471
[86,] 0.012124886 0.02424977 0.9878751
[87,] 0.041872815 0.08374563 0.9581272
[88,] 0.038353210 0.07670642 0.9616468
[89,] 0.044937205 0.08987441 0.9550628
[90,] 0.034120767 0.06824153 0.9658792
[91,] 0.026258801 0.05251760 0.9737412
[92,] 0.019405077 0.03881015 0.9805949
[93,] 0.015510181 0.03102036 0.9844898
[94,] 0.011444805 0.02288961 0.9885552
[95,] 0.008400660 0.01680132 0.9915993
[96,] 0.008847203 0.01769441 0.9911528
[97,] 0.009750629 0.01950126 0.9902494
[98,] 0.012787460 0.02557492 0.9872125
[99,] 0.070195426 0.14039085 0.9298046
[100,] 0.101520110 0.20304022 0.8984799
[101,] 0.087466789 0.17493358 0.9125332
[102,] 0.070205323 0.14041065 0.9297947
[103,] 0.058914295 0.11782859 0.9410857
[104,] 0.042910342 0.08582068 0.9570897
[105,] 0.030371280 0.06074256 0.9696287
[106,] 0.046283404 0.09256681 0.9537166
[107,] 0.031948052 0.06389610 0.9680519
[108,] 0.021984840 0.04396968 0.9780152
[109,] 0.018749189 0.03749838 0.9812508
[110,] 0.033678099 0.06735620 0.9663219
[111,] 0.249661765 0.49932353 0.7503382
[112,] 0.287547499 0.57509500 0.7124525
[113,] 0.223954366 0.44790873 0.7760456
[114,] 0.184600311 0.36920062 0.8153997
[115,] 0.129977719 0.25995544 0.8700223
[116,] 0.086430337 0.17286067 0.9135697
[117,] 0.053454507 0.10690901 0.9465455
[118,] 0.060687260 0.12137452 0.9393127
[119,] 0.033348922 0.06669784 0.9666511
[120,] 0.016351099 0.03270220 0.9836489
[121,] 0.015174248 0.03034850 0.9848258
> postscript(file="/var/www/html/freestat/rcomp/tmp/1tyq41293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index')
> points(x[,1]-mysum$resid)
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/24qq71293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/34qq71293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/44qq71293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals')
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/5xzps1293562311.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 = 130
Frequency = 1
1 2 3 4 5 6
-7.5247525 -41.5247525 26.4752475 -111.5247525 17.4752475 28.4752475
7 8 9 10 11 12
-246.5247525 -283.5247525 -3.5247525 25.4752475 -22.5247525 -118.5247525
13 14 15 16 17 18
-14.5247525 18.4752475 -48.5247525 -14.5247525 7.4752475 -240.5247525
19 20 21 22 23 24
-272.5247525 10.4752475 124.4752475 -17.5247525 10.4752475 1.4752475
25 26 27 28 29 30
-42.5247525 66.4752475 -76.5247525 -9.5247525 52.4752475 -207.5247525
31 32 33 34 35 36
-321.5247525 125.4752475 61.4752475 15.4752475 14.4752475 -20.5247525
37 38 39 40 41 42
22.4752475 62.4752475 -1.5247525 5.4752475 102.4752475 -153.5247525
43 44 45 46 47 48
-291.5247525 174.4752475 93.4752475 -38.5247525 95.4752475 -1.5247525
49 50 51 52 53 54
67.4752475 196.4752475 48.4752475 27.4752475 156.4752475 -216.5247525
55 56 57 58 59 60
-276.5247525 210.4752475 100.4752475 67.4752475 12.4752475 66.4752475
61 62 63 64 65 66
9.4752475 133.4752475 6.4752475 92.4752475 225.4752475 -210.5247525
67 68 69 70 71 72
-261.5247525 195.4752475 58.4752475 -27.5247525 47.4752475 111.4752475
73 74 75 76 77 78
62.4752475 54.4752475 -34.5247525 100.4752475 102.4752475 -242.5247525
79 80 81 82 83 84
-297.5247525 77.4752475 86.4752475 28.4752475 24.4752475 13.4752475
85 86 87 88 89 90
14.4752475 89.4752475 25.4752475 3.4752475 102.4752475 -236.5247525
91 92 93 94 95 96
-284.5247525 163.4752475 223.4752475 20.4752475 0.4752475 56.4752475
97 98 99 100 101 102
-11.5247525 86.4752475 86.4752475 0.4752475 287.4752475 -242.7241379
103 104 105 106 107 108
-416.7241379 143.2758621 54.2758621 -65.7241379 84.2758621 1.2758621
109 110 111 112 113 114
11.2758621 232.2758621 15.2758621 -39.7241379 132.2758621 -241.7241379
115 116 117 118 119 120
-386.7241379 217.2758621 67.2758621 118.2758621 8.2758621 39.2758621
121 122 123 124 125 126
19.2758621 204.2758621 1.2758621 23.2758621 172.2758621 -179.7241379
127 128 129 130
-335.7241379 248.2758621 94.2758621 21.2758621
> postscript(file="/var/www/html/freestat/rcomp/tmp/6xzps1293562311.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 = 130
Frequency = 1
lag(myerror, k = 1) myerror
0 -7.5247525 NA
1 -41.5247525 -7.5247525
2 26.4752475 -41.5247525
3 -111.5247525 26.4752475
4 17.4752475 -111.5247525
5 28.4752475 17.4752475
6 -246.5247525 28.4752475
7 -283.5247525 -246.5247525
8 -3.5247525 -283.5247525
9 25.4752475 -3.5247525
10 -22.5247525 25.4752475
11 -118.5247525 -22.5247525
12 -14.5247525 -118.5247525
13 18.4752475 -14.5247525
14 -48.5247525 18.4752475
15 -14.5247525 -48.5247525
16 7.4752475 -14.5247525
17 -240.5247525 7.4752475
18 -272.5247525 -240.5247525
19 10.4752475 -272.5247525
20 124.4752475 10.4752475
21 -17.5247525 124.4752475
22 10.4752475 -17.5247525
23 1.4752475 10.4752475
24 -42.5247525 1.4752475
25 66.4752475 -42.5247525
26 -76.5247525 66.4752475
27 -9.5247525 -76.5247525
28 52.4752475 -9.5247525
29 -207.5247525 52.4752475
30 -321.5247525 -207.5247525
31 125.4752475 -321.5247525
32 61.4752475 125.4752475
33 15.4752475 61.4752475
34 14.4752475 15.4752475
35 -20.5247525 14.4752475
36 22.4752475 -20.5247525
37 62.4752475 22.4752475
38 -1.5247525 62.4752475
39 5.4752475 -1.5247525
40 102.4752475 5.4752475
41 -153.5247525 102.4752475
42 -291.5247525 -153.5247525
43 174.4752475 -291.5247525
44 93.4752475 174.4752475
45 -38.5247525 93.4752475
46 95.4752475 -38.5247525
47 -1.5247525 95.4752475
48 67.4752475 -1.5247525
49 196.4752475 67.4752475
50 48.4752475 196.4752475
51 27.4752475 48.4752475
52 156.4752475 27.4752475
53 -216.5247525 156.4752475
54 -276.5247525 -216.5247525
55 210.4752475 -276.5247525
56 100.4752475 210.4752475
57 67.4752475 100.4752475
58 12.4752475 67.4752475
59 66.4752475 12.4752475
60 9.4752475 66.4752475
61 133.4752475 9.4752475
62 6.4752475 133.4752475
63 92.4752475 6.4752475
64 225.4752475 92.4752475
65 -210.5247525 225.4752475
66 -261.5247525 -210.5247525
67 195.4752475 -261.5247525
68 58.4752475 195.4752475
69 -27.5247525 58.4752475
70 47.4752475 -27.5247525
71 111.4752475 47.4752475
72 62.4752475 111.4752475
73 54.4752475 62.4752475
74 -34.5247525 54.4752475
75 100.4752475 -34.5247525
76 102.4752475 100.4752475
77 -242.5247525 102.4752475
78 -297.5247525 -242.5247525
79 77.4752475 -297.5247525
80 86.4752475 77.4752475
81 28.4752475 86.4752475
82 24.4752475 28.4752475
83 13.4752475 24.4752475
84 14.4752475 13.4752475
85 89.4752475 14.4752475
86 25.4752475 89.4752475
87 3.4752475 25.4752475
88 102.4752475 3.4752475
89 -236.5247525 102.4752475
90 -284.5247525 -236.5247525
91 163.4752475 -284.5247525
92 223.4752475 163.4752475
93 20.4752475 223.4752475
94 0.4752475 20.4752475
95 56.4752475 0.4752475
96 -11.5247525 56.4752475
97 86.4752475 -11.5247525
98 86.4752475 86.4752475
99 0.4752475 86.4752475
100 287.4752475 0.4752475
101 -242.7241379 287.4752475
102 -416.7241379 -242.7241379
103 143.2758621 -416.7241379
104 54.2758621 143.2758621
105 -65.7241379 54.2758621
106 84.2758621 -65.7241379
107 1.2758621 84.2758621
108 11.2758621 1.2758621
109 232.2758621 11.2758621
110 15.2758621 232.2758621
111 -39.7241379 15.2758621
112 132.2758621 -39.7241379
113 -241.7241379 132.2758621
114 -386.7241379 -241.7241379
115 217.2758621 -386.7241379
116 67.2758621 217.2758621
117 118.2758621 67.2758621
118 8.2758621 118.2758621
119 39.2758621 8.2758621
120 19.2758621 39.2758621
121 204.2758621 19.2758621
122 1.2758621 204.2758621
123 23.2758621 1.2758621
124 172.2758621 23.2758621
125 -179.7241379 172.2758621
126 -335.7241379 -179.7241379
127 248.2758621 -335.7241379
128 94.2758621 248.2758621
129 21.2758621 94.2758621
130 NA 21.2758621
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -41.5247525 -7.5247525
[2,] 26.4752475 -41.5247525
[3,] -111.5247525 26.4752475
[4,] 17.4752475 -111.5247525
[5,] 28.4752475 17.4752475
[6,] -246.5247525 28.4752475
[7,] -283.5247525 -246.5247525
[8,] -3.5247525 -283.5247525
[9,] 25.4752475 -3.5247525
[10,] -22.5247525 25.4752475
[11,] -118.5247525 -22.5247525
[12,] -14.5247525 -118.5247525
[13,] 18.4752475 -14.5247525
[14,] -48.5247525 18.4752475
[15,] -14.5247525 -48.5247525
[16,] 7.4752475 -14.5247525
[17,] -240.5247525 7.4752475
[18,] -272.5247525 -240.5247525
[19,] 10.4752475 -272.5247525
[20,] 124.4752475 10.4752475
[21,] -17.5247525 124.4752475
[22,] 10.4752475 -17.5247525
[23,] 1.4752475 10.4752475
[24,] -42.5247525 1.4752475
[25,] 66.4752475 -42.5247525
[26,] -76.5247525 66.4752475
[27,] -9.5247525 -76.5247525
[28,] 52.4752475 -9.5247525
[29,] -207.5247525 52.4752475
[30,] -321.5247525 -207.5247525
[31,] 125.4752475 -321.5247525
[32,] 61.4752475 125.4752475
[33,] 15.4752475 61.4752475
[34,] 14.4752475 15.4752475
[35,] -20.5247525 14.4752475
[36,] 22.4752475 -20.5247525
[37,] 62.4752475 22.4752475
[38,] -1.5247525 62.4752475
[39,] 5.4752475 -1.5247525
[40,] 102.4752475 5.4752475
[41,] -153.5247525 102.4752475
[42,] -291.5247525 -153.5247525
[43,] 174.4752475 -291.5247525
[44,] 93.4752475 174.4752475
[45,] -38.5247525 93.4752475
[46,] 95.4752475 -38.5247525
[47,] -1.5247525 95.4752475
[48,] 67.4752475 -1.5247525
[49,] 196.4752475 67.4752475
[50,] 48.4752475 196.4752475
[51,] 27.4752475 48.4752475
[52,] 156.4752475 27.4752475
[53,] -216.5247525 156.4752475
[54,] -276.5247525 -216.5247525
[55,] 210.4752475 -276.5247525
[56,] 100.4752475 210.4752475
[57,] 67.4752475 100.4752475
[58,] 12.4752475 67.4752475
[59,] 66.4752475 12.4752475
[60,] 9.4752475 66.4752475
[61,] 133.4752475 9.4752475
[62,] 6.4752475 133.4752475
[63,] 92.4752475 6.4752475
[64,] 225.4752475 92.4752475
[65,] -210.5247525 225.4752475
[66,] -261.5247525 -210.5247525
[67,] 195.4752475 -261.5247525
[68,] 58.4752475 195.4752475
[69,] -27.5247525 58.4752475
[70,] 47.4752475 -27.5247525
[71,] 111.4752475 47.4752475
[72,] 62.4752475 111.4752475
[73,] 54.4752475 62.4752475
[74,] -34.5247525 54.4752475
[75,] 100.4752475 -34.5247525
[76,] 102.4752475 100.4752475
[77,] -242.5247525 102.4752475
[78,] -297.5247525 -242.5247525
[79,] 77.4752475 -297.5247525
[80,] 86.4752475 77.4752475
[81,] 28.4752475 86.4752475
[82,] 24.4752475 28.4752475
[83,] 13.4752475 24.4752475
[84,] 14.4752475 13.4752475
[85,] 89.4752475 14.4752475
[86,] 25.4752475 89.4752475
[87,] 3.4752475 25.4752475
[88,] 102.4752475 3.4752475
[89,] -236.5247525 102.4752475
[90,] -284.5247525 -236.5247525
[91,] 163.4752475 -284.5247525
[92,] 223.4752475 163.4752475
[93,] 20.4752475 223.4752475
[94,] 0.4752475 20.4752475
[95,] 56.4752475 0.4752475
[96,] -11.5247525 56.4752475
[97,] 86.4752475 -11.5247525
[98,] 86.4752475 86.4752475
[99,] 0.4752475 86.4752475
[100,] 287.4752475 0.4752475
[101,] -242.7241379 287.4752475
[102,] -416.7241379 -242.7241379
[103,] 143.2758621 -416.7241379
[104,] 54.2758621 143.2758621
[105,] -65.7241379 54.2758621
[106,] 84.2758621 -65.7241379
[107,] 1.2758621 84.2758621
[108,] 11.2758621 1.2758621
[109,] 232.2758621 11.2758621
[110,] 15.2758621 232.2758621
[111,] -39.7241379 15.2758621
[112,] 132.2758621 -39.7241379
[113,] -241.7241379 132.2758621
[114,] -386.7241379 -241.7241379
[115,] 217.2758621 -386.7241379
[116,] 67.2758621 217.2758621
[117,] 118.2758621 67.2758621
[118,] 8.2758621 118.2758621
[119,] 39.2758621 8.2758621
[120,] 19.2758621 39.2758621
[121,] 204.2758621 19.2758621
[122,] 1.2758621 204.2758621
[123,] 23.2758621 1.2758621
[124,] 172.2758621 23.2758621
[125,] -179.7241379 172.2758621
[126,] -335.7241379 -179.7241379
[127,] 248.2758621 -335.7241379
[128,] 94.2758621 248.2758621
[129,] 21.2758621 94.2758621
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -41.5247525 -7.5247525
2 26.4752475 -41.5247525
3 -111.5247525 26.4752475
4 17.4752475 -111.5247525
5 28.4752475 17.4752475
6 -246.5247525 28.4752475
7 -283.5247525 -246.5247525
8 -3.5247525 -283.5247525
9 25.4752475 -3.5247525
10 -22.5247525 25.4752475
11 -118.5247525 -22.5247525
12 -14.5247525 -118.5247525
13 18.4752475 -14.5247525
14 -48.5247525 18.4752475
15 -14.5247525 -48.5247525
16 7.4752475 -14.5247525
17 -240.5247525 7.4752475
18 -272.5247525 -240.5247525
19 10.4752475 -272.5247525
20 124.4752475 10.4752475
21 -17.5247525 124.4752475
22 10.4752475 -17.5247525
23 1.4752475 10.4752475
24 -42.5247525 1.4752475
25 66.4752475 -42.5247525
26 -76.5247525 66.4752475
27 -9.5247525 -76.5247525
28 52.4752475 -9.5247525
29 -207.5247525 52.4752475
30 -321.5247525 -207.5247525
31 125.4752475 -321.5247525
32 61.4752475 125.4752475
33 15.4752475 61.4752475
34 14.4752475 15.4752475
35 -20.5247525 14.4752475
36 22.4752475 -20.5247525
37 62.4752475 22.4752475
38 -1.5247525 62.4752475
39 5.4752475 -1.5247525
40 102.4752475 5.4752475
41 -153.5247525 102.4752475
42 -291.5247525 -153.5247525
43 174.4752475 -291.5247525
44 93.4752475 174.4752475
45 -38.5247525 93.4752475
46 95.4752475 -38.5247525
47 -1.5247525 95.4752475
48 67.4752475 -1.5247525
49 196.4752475 67.4752475
50 48.4752475 196.4752475
51 27.4752475 48.4752475
52 156.4752475 27.4752475
53 -216.5247525 156.4752475
54 -276.5247525 -216.5247525
55 210.4752475 -276.5247525
56 100.4752475 210.4752475
57 67.4752475 100.4752475
58 12.4752475 67.4752475
59 66.4752475 12.4752475
60 9.4752475 66.4752475
61 133.4752475 9.4752475
62 6.4752475 133.4752475
63 92.4752475 6.4752475
64 225.4752475 92.4752475
65 -210.5247525 225.4752475
66 -261.5247525 -210.5247525
67 195.4752475 -261.5247525
68 58.4752475 195.4752475
69 -27.5247525 58.4752475
70 47.4752475 -27.5247525
71 111.4752475 47.4752475
72 62.4752475 111.4752475
73 54.4752475 62.4752475
74 -34.5247525 54.4752475
75 100.4752475 -34.5247525
76 102.4752475 100.4752475
77 -242.5247525 102.4752475
78 -297.5247525 -242.5247525
79 77.4752475 -297.5247525
80 86.4752475 77.4752475
81 28.4752475 86.4752475
82 24.4752475 28.4752475
83 13.4752475 24.4752475
84 14.4752475 13.4752475
85 89.4752475 14.4752475
86 25.4752475 89.4752475
87 3.4752475 25.4752475
88 102.4752475 3.4752475
89 -236.5247525 102.4752475
90 -284.5247525 -236.5247525
91 163.4752475 -284.5247525
92 223.4752475 163.4752475
93 20.4752475 223.4752475
94 0.4752475 20.4752475
95 56.4752475 0.4752475
96 -11.5247525 56.4752475
97 86.4752475 -11.5247525
98 86.4752475 86.4752475
99 0.4752475 86.4752475
100 287.4752475 0.4752475
101 -242.7241379 287.4752475
102 -416.7241379 -242.7241379
103 143.2758621 -416.7241379
104 54.2758621 143.2758621
105 -65.7241379 54.2758621
106 84.2758621 -65.7241379
107 1.2758621 84.2758621
108 11.2758621 1.2758621
109 232.2758621 11.2758621
110 15.2758621 232.2758621
111 -39.7241379 15.2758621
112 132.2758621 -39.7241379
113 -241.7241379 132.2758621
114 -386.7241379 -241.7241379
115 217.2758621 -386.7241379
116 67.2758621 217.2758621
117 118.2758621 67.2758621
118 8.2758621 118.2758621
119 39.2758621 8.2758621
120 19.2758621 39.2758621
121 204.2758621 19.2758621
122 1.2758621 204.2758621
123 23.2758621 1.2758621
124 172.2758621 23.2758621
125 -179.7241379 172.2758621
126 -335.7241379 -179.7241379
127 248.2758621 -335.7241379
128 94.2758621 248.2758621
129 21.2758621 94.2758621
> plot(z,main=paste('Residual Lag plot, lowess, and regression line'), ylab='values of Residuals', xlab='lagged values of Residuals')
> lines(lowess(z))
> abline(lm(z))
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/778od1293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/878od1293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/9ihng1293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
> plot(mylm, las = 1, sub='Residual Diagnostics')
> par(opar)
> dev.off()
null device
1
> if (n > n25) {
+ postscript(file="/var/www/html/freestat/rcomp/tmp/10ihng1293562311.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
+ plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint')
+ grid()
+ dev.off()
+ }
null device
1
>
> #Note: the /var/www/html/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/html/freestat/rcomp/createtable")
>
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Estimated Regression Equation', 1, TRUE)
> a<-table.row.end(a)
> myeq <- colnames(x)[1]
> myeq <- paste(myeq, '[t] = ', sep='')
> for (i in 1:k){
+ if (mysum$coefficients[i,1] > 0) myeq <- paste(myeq, '+', '')
+ myeq <- paste(myeq, mysum$coefficients[i,1], sep=' ')
+ if (rownames(mysum$coefficients)[i] != '(Intercept)') {
+ myeq <- paste(myeq, rownames(mysum$coefficients)[i], sep='')
+ if (rownames(mysum$coefficients)[i] != 't') myeq <- paste(myeq, '[t]', sep='')
+ }
+ }
> myeq <- paste(myeq, ' + e[t]')
> a<-table.row.start(a)
> a<-table.element(a, myeq)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/113iml1293562311.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,hyperlink('http://www.xycoon.com/ols1.htm','Multiple Linear Regression - Ordinary Least Squares',''), 6, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Variable',header=TRUE)
> a<-table.element(a,'Parameter',header=TRUE)
> a<-table.element(a,'S.D.',header=TRUE)
> a<-table.element(a,'T-STAT
H0: parameter = 0',header=TRUE)
> a<-table.element(a,'2-tail p-value',header=TRUE)
> a<-table.element(a,'1-tail p-value',header=TRUE)
> a<-table.row.end(a)
> for (i in 1:k){
+ a<-table.row.start(a)
+ a<-table.element(a,rownames(mysum$coefficients)[i],header=TRUE)
+ a<-table.element(a,mysum$coefficients[i,1])
+ a<-table.element(a, round(mysum$coefficients[i,2],6))
+ a<-table.element(a, round(mysum$coefficients[i,3],4))
+ a<-table.element(a, round(mysum$coefficients[i,4],6))
+ a<-table.element(a, round(mysum$coefficients[i,4]/2,6))
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/127jkr1293562311.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Regression Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple R',1,TRUE)
> a<-table.element(a, sqrt(mysum$r.squared))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'R-squared',1,TRUE)
> a<-table.element(a, mysum$r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Adjusted R-squared',1,TRUE)
> a<-table.element(a, mysum$adj.r.squared)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (value)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE)
> a<-table.element(a, mysum$fstatistic[3])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'p-value',1,TRUE)
> a<-table.element(a, 1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3]))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Residual Statistics', 2, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Residual Standard Deviation',1,TRUE)
> a<-table.element(a, mysum$sigma)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Sum Squared Residuals',1,TRUE)
> a<-table.element(a, sum(myerror*myerror))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/13lsi01293562311.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a, 'Multiple Linear Regression - Actuals, Interpolation, and Residuals', 4, TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Time or Index', 1, TRUE)
> a<-table.element(a, 'Actuals', 1, TRUE)
> a<-table.element(a, 'Interpolation
Forecast', 1, TRUE)
> a<-table.element(a, 'Residuals
Prediction Error', 1, TRUE)
> a<-table.row.end(a)
> for (i in 1:n) {
+ a<-table.row.start(a)
+ a<-table.element(a,i, 1, TRUE)
+ a<-table.element(a,x[i])
+ a<-table.element(a,x[i]-mysum$resid[i])
+ a<-table.element(a,mysum$resid[i])
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/146bh61293562311.tab")
> if (n > n25) {
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'p-values',header=TRUE)
+ a<-table.element(a,'Alternative Hypothesis',3,header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'breakpoint index',header=TRUE)
+ a<-table.element(a,'greater',header=TRUE)
+ a<-table.element(a,'2-sided',header=TRUE)
+ a<-table.element(a,'less',header=TRUE)
+ a<-table.row.end(a)
+ for (mypoint in kp3:nmkm3) {
+ a<-table.row.start(a)
+ a<-table.element(a,mypoint,header=TRUE)
+ a<-table.element(a,gqarr[mypoint-kp3+1,1])
+ a<-table.element(a,gqarr[mypoint-kp3+1,2])
+ a<-table.element(a,gqarr[mypoint-kp3+1,3])
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/www/html/freestat/rcomp/tmp/15atxc1293562311.tab")
+ a<-table.start()
+ a<-table.row.start(a)
+ a<-table.element(a,'Meta Analysis of Goldfeld-Quandt test for Heteroskedasticity',4,TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'Description',header=TRUE)
+ a<-table.element(a,'# significant tests',header=TRUE)
+ a<-table.element(a,'% significant tests',header=TRUE)
+ a<-table.element(a,'OK/NOK',header=TRUE)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'1% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant1)
+ a<-table.element(a,numsignificant1/numgqtests)
+ if (numsignificant1/numgqtests < 0.01) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'5% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant5)
+ a<-table.element(a,numsignificant5/numgqtests)
+ if (numsignificant5/numgqtests < 0.05) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.row.start(a)
+ a<-table.element(a,'10% type I error level',header=TRUE)
+ a<-table.element(a,numsignificant10)
+ a<-table.element(a,numsignificant10/numgqtests)
+ if (numsignificant10/numgqtests < 0.1) dum <- 'OK' else dum <- 'NOK'
+ a<-table.element(a,dum)
+ a<-table.row.end(a)
+ a<-table.end(a)
+ table.save(a,file="/var/www/html/freestat/rcomp/tmp/16old31293562311.tab")
+ }
>
> try(system("convert tmp/1tyq41293562311.ps tmp/1tyq41293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/24qq71293562311.ps tmp/24qq71293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/34qq71293562311.ps tmp/34qq71293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/44qq71293562311.ps tmp/44qq71293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/5xzps1293562311.ps tmp/5xzps1293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/6xzps1293562311.ps tmp/6xzps1293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/778od1293562311.ps tmp/778od1293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/878od1293562311.ps tmp/878od1293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/9ihng1293562311.ps tmp/9ihng1293562311.png",intern=TRUE))
character(0)
> try(system("convert tmp/10ihng1293562311.ps tmp/10ihng1293562311.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
4.860 2.665 5.184