R version 2.12.0 (2010-10-15)
Copyright (C) 2010 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(989236,1008380,1207763,1368839,1469798,1498721,1761769,1653214,1599104,1421179,1163995,1037735,1015407,1039210,1258049,1469445,1552346,1549144,1785895,1662335,1629440,1467430,1202209,1076982,1039367,1063449,1335135,1491602,1591972,1641248,1898849,1798580,1762444,1622044,1368955,1262973,1195650,1269530,1479279,1607819,1712466,1721766,1949843,1821326,1757802,1590367,1260647,1149235,1016367,1027885,1262159,1520854,1544144,1564709,1821776,1741365,1623386,1498658,1241822,1136029),dim=c(1,60),dimnames=list(c('PassengersBRU'),1:60))
> y <- array(NA,dim=c(1,60),dimnames=list(c('PassengersBRU'),1:60))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par3 = 'Linear Trend'
> par2 = 'Include Monthly Dummies'
> par1 = '1'
> #'GNU S' R Code compiled by R2WASP v. 1.0.44 ()
> #Author: Prof. Dr. P. Wessa
> #To cite this work: AUTHOR(S), (YEAR), YOUR SOFTWARE TITLE (vNUMBER) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_YOURPAGE.wasp/
> #Source of accompanying publication: Office for Research, Development, and Education
> #Technical description: Write here your technical program description (don't use hard returns!)
> library(lattice)
> library(lmtest)
Loading required package: zoo
> 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
PassengersBRU M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 t
1 989236 1 0 0 0 0 0 0 0 0 0 0 1
2 1008380 0 1 0 0 0 0 0 0 0 0 0 2
3 1207763 0 0 1 0 0 0 0 0 0 0 0 3
4 1368839 0 0 0 1 0 0 0 0 0 0 0 4
5 1469798 0 0 0 0 1 0 0 0 0 0 0 5
6 1498721 0 0 0 0 0 1 0 0 0 0 0 6
7 1761769 0 0 0 0 0 0 1 0 0 0 0 7
8 1653214 0 0 0 0 0 0 0 1 0 0 0 8
9 1599104 0 0 0 0 0 0 0 0 1 0 0 9
10 1421179 0 0 0 0 0 0 0 0 0 1 0 10
11 1163995 0 0 0 0 0 0 0 0 0 0 1 11
12 1037735 0 0 0 0 0 0 0 0 0 0 0 12
13 1015407 1 0 0 0 0 0 0 0 0 0 0 13
14 1039210 0 1 0 0 0 0 0 0 0 0 0 14
15 1258049 0 0 1 0 0 0 0 0 0 0 0 15
16 1469445 0 0 0 1 0 0 0 0 0 0 0 16
17 1552346 0 0 0 0 1 0 0 0 0 0 0 17
18 1549144 0 0 0 0 0 1 0 0 0 0 0 18
19 1785895 0 0 0 0 0 0 1 0 0 0 0 19
20 1662335 0 0 0 0 0 0 0 1 0 0 0 20
21 1629440 0 0 0 0 0 0 0 0 1 0 0 21
22 1467430 0 0 0 0 0 0 0 0 0 1 0 22
23 1202209 0 0 0 0 0 0 0 0 0 0 1 23
24 1076982 0 0 0 0 0 0 0 0 0 0 0 24
25 1039367 1 0 0 0 0 0 0 0 0 0 0 25
26 1063449 0 1 0 0 0 0 0 0 0 0 0 26
27 1335135 0 0 1 0 0 0 0 0 0 0 0 27
28 1491602 0 0 0 1 0 0 0 0 0 0 0 28
29 1591972 0 0 0 0 1 0 0 0 0 0 0 29
30 1641248 0 0 0 0 0 1 0 0 0 0 0 30
31 1898849 0 0 0 0 0 0 1 0 0 0 0 31
32 1798580 0 0 0 0 0 0 0 1 0 0 0 32
33 1762444 0 0 0 0 0 0 0 0 1 0 0 33
34 1622044 0 0 0 0 0 0 0 0 0 1 0 34
35 1368955 0 0 0 0 0 0 0 0 0 0 1 35
36 1262973 0 0 0 0 0 0 0 0 0 0 0 36
37 1195650 1 0 0 0 0 0 0 0 0 0 0 37
38 1269530 0 1 0 0 0 0 0 0 0 0 0 38
39 1479279 0 0 1 0 0 0 0 0 0 0 0 39
40 1607819 0 0 0 1 0 0 0 0 0 0 0 40
41 1712466 0 0 0 0 1 0 0 0 0 0 0 41
42 1721766 0 0 0 0 0 1 0 0 0 0 0 42
43 1949843 0 0 0 0 0 0 1 0 0 0 0 43
44 1821326 0 0 0 0 0 0 0 1 0 0 0 44
45 1757802 0 0 0 0 0 0 0 0 1 0 0 45
46 1590367 0 0 0 0 0 0 0 0 0 1 0 46
47 1260647 0 0 0 0 0 0 0 0 0 0 1 47
48 1149235 0 0 0 0 0 0 0 0 0 0 0 48
49 1016367 1 0 0 0 0 0 0 0 0 0 0 49
50 1027885 0 1 0 0 0 0 0 0 0 0 0 50
51 1262159 0 0 1 0 0 0 0 0 0 0 0 51
52 1520854 0 0 0 1 0 0 0 0 0 0 0 52
53 1544144 0 0 0 0 1 0 0 0 0 0 0 53
54 1564709 0 0 0 0 0 1 0 0 0 0 0 54
55 1821776 0 0 0 0 0 0 1 0 0 0 0 55
56 1741365 0 0 0 0 0 0 0 1 0 0 0 56
57 1623386 0 0 0 0 0 0 0 0 1 0 0 57
58 1498658 0 0 0 0 0 0 0 0 0 1 0 58
59 1241822 0 0 0 0 0 0 0 0 0 0 1 59
60 1136029 0 0 0 0 0 0 0 0 0 0 0 60
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) M1 M2 M3 M4 M5
1046424 -55057 -26965 197428 378269 458309
M6 M7 M8 M9 M10 M11
476888 723003 612347 549025 392132 117328
t
2394
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-111250 -41801 -16434 48259 159117
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1046423.8 39446.6 26.528 < 2e-16 ***
M1 -55056.6 47989.0 -1.147 0.257075
M2 -26964.7 47917.3 -0.563 0.576290
M3 197427.9 47852.4 4.126 0.000150 ***
M4 378269.2 47794.2 7.915 3.42e-10 ***
M5 458309.1 47742.7 9.600 1.18e-12 ***
M6 476888.0 47698.1 9.998 3.22e-13 ***
M7 723003.2 47660.4 15.170 < 2e-16 ***
M8 612347.3 47629.4 12.856 < 2e-16 ***
M9 549025.0 47605.4 11.533 2.64e-15 ***
M10 392131.9 47588.2 8.240 1.12e-10 ***
M11 117328.3 47577.9 2.466 0.017367 *
t 2393.5 572.2 4.183 0.000125 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 75220 on 47 degrees of freedom
Multiple R-squared: 0.938, Adjusted R-squared: 0.9221
F-statistic: 59.23 on 12 and 47 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,] 3.069457e-02 6.138913e-02 0.9693054
[2,] 1.021270e-02 2.042539e-02 0.9897873
[3,] 2.492801e-03 4.985602e-03 0.9975072
[4,] 1.002046e-03 2.004092e-03 0.9989980
[5,] 6.559263e-04 1.311853e-03 0.9993441
[6,] 2.050239e-04 4.100478e-04 0.9997950
[7,] 6.401790e-05 1.280358e-04 0.9999360
[8,] 2.121062e-05 4.242124e-05 0.9999788
[9,] 9.915003e-06 1.983001e-05 0.9999901
[10,] 6.545748e-06 1.309150e-05 0.9999935
[11,] 4.703894e-06 9.407788e-06 0.9999953
[12,] 8.764098e-06 1.752820e-05 0.9999912
[13,] 8.767677e-06 1.753535e-05 0.9999912
[14,] 9.392111e-06 1.878422e-05 0.9999906
[15,] 3.981747e-05 7.963495e-05 0.9999602
[16,] 1.665183e-04 3.330367e-04 0.9998335
[17,] 1.286674e-03 2.573347e-03 0.9987133
[18,] 3.797007e-03 7.594014e-03 0.9962030
[19,] 1.992027e-02 3.984055e-02 0.9800797
[20,] 4.484166e-02 8.968331e-02 0.9551583
[21,] 1.015502e-01 2.031004e-01 0.8984498
[22,] 8.705864e-02 1.741173e-01 0.9129414
[23,] 2.513921e-01 5.027841e-01 0.7486079
[24,] 4.617478e-01 9.234955e-01 0.5382522
[25,] 3.487466e-01 6.974931e-01 0.6512534
[26,] 3.873123e-01 7.746246e-01 0.6126877
[27,] 4.479781e-01 8.959562e-01 0.5520219
[28,] 4.390648e-01 8.781297e-01 0.5609352
[29,] 3.015527e-01 6.031054e-01 0.6984473
> postscript(file="/var/www/rcomp/tmp/1b5gh1292683443.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/rcomp/tmp/2b5gh1292683443.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/rcomp/tmp/3mwyk1292683443.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/rcomp/tmp/4mwyk1292683443.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/rcomp/tmp/5mwyk1292683443.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 = 60
Frequency = 1
1 2 3 4 5 6
-4524.733 -15866.133 -43269.333 -65428.133 -46902.533 -38951.933
7 8 9 10 11 12
-24412.733 -24705.333 -17886.533 -41311.933 -26085.933 -37411.133
13 14 15 16 17 18
-7076.067 -13758.467 -21705.667 6455.533 6923.133 -17251.267
19 20 21 22 23 24
-29009.067 -44306.667 -16272.867 -23783.267 -16594.267 -26886.467
25 26 27 28 29 30
-11838.400 -18241.800 26658.000 -109.800 17826.800 46130.400
31 32 33 34 35 36
55222.600 63216.000 88008.800 102108.400 121429.400 130382.200
37 38 39 40 41 42
115722.267 159116.867 142079.667 87384.867 109598.467 97926.067
43 44 45 46 47 48
77494.267 57239.667 54644.467 41709.067 -15600.933 -12078.133
49 50 51 52 53 54
-92283.067 -111250.467 -103762.667 -28302.467 -87445.867 -87853.267
55 56 57 58 59 60
-79295.067 -51443.667 -108493.867 -78722.267 -63148.267 -54006.467
> postscript(file="/var/www/rcomp/tmp/6xnxn1292683443.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 = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 -4524.733 NA
1 -15866.133 -4524.733
2 -43269.333 -15866.133
3 -65428.133 -43269.333
4 -46902.533 -65428.133
5 -38951.933 -46902.533
6 -24412.733 -38951.933
7 -24705.333 -24412.733
8 -17886.533 -24705.333
9 -41311.933 -17886.533
10 -26085.933 -41311.933
11 -37411.133 -26085.933
12 -7076.067 -37411.133
13 -13758.467 -7076.067
14 -21705.667 -13758.467
15 6455.533 -21705.667
16 6923.133 6455.533
17 -17251.267 6923.133
18 -29009.067 -17251.267
19 -44306.667 -29009.067
20 -16272.867 -44306.667
21 -23783.267 -16272.867
22 -16594.267 -23783.267
23 -26886.467 -16594.267
24 -11838.400 -26886.467
25 -18241.800 -11838.400
26 26658.000 -18241.800
27 -109.800 26658.000
28 17826.800 -109.800
29 46130.400 17826.800
30 55222.600 46130.400
31 63216.000 55222.600
32 88008.800 63216.000
33 102108.400 88008.800
34 121429.400 102108.400
35 130382.200 121429.400
36 115722.267 130382.200
37 159116.867 115722.267
38 142079.667 159116.867
39 87384.867 142079.667
40 109598.467 87384.867
41 97926.067 109598.467
42 77494.267 97926.067
43 57239.667 77494.267
44 54644.467 57239.667
45 41709.067 54644.467
46 -15600.933 41709.067
47 -12078.133 -15600.933
48 -92283.067 -12078.133
49 -111250.467 -92283.067
50 -103762.667 -111250.467
51 -28302.467 -103762.667
52 -87445.867 -28302.467
53 -87853.267 -87445.867
54 -79295.067 -87853.267
55 -51443.667 -79295.067
56 -108493.867 -51443.667
57 -78722.267 -108493.867
58 -63148.267 -78722.267
59 -54006.467 -63148.267
60 NA -54006.467
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -15866.133 -4524.733
[2,] -43269.333 -15866.133
[3,] -65428.133 -43269.333
[4,] -46902.533 -65428.133
[5,] -38951.933 -46902.533
[6,] -24412.733 -38951.933
[7,] -24705.333 -24412.733
[8,] -17886.533 -24705.333
[9,] -41311.933 -17886.533
[10,] -26085.933 -41311.933
[11,] -37411.133 -26085.933
[12,] -7076.067 -37411.133
[13,] -13758.467 -7076.067
[14,] -21705.667 -13758.467
[15,] 6455.533 -21705.667
[16,] 6923.133 6455.533
[17,] -17251.267 6923.133
[18,] -29009.067 -17251.267
[19,] -44306.667 -29009.067
[20,] -16272.867 -44306.667
[21,] -23783.267 -16272.867
[22,] -16594.267 -23783.267
[23,] -26886.467 -16594.267
[24,] -11838.400 -26886.467
[25,] -18241.800 -11838.400
[26,] 26658.000 -18241.800
[27,] -109.800 26658.000
[28,] 17826.800 -109.800
[29,] 46130.400 17826.800
[30,] 55222.600 46130.400
[31,] 63216.000 55222.600
[32,] 88008.800 63216.000
[33,] 102108.400 88008.800
[34,] 121429.400 102108.400
[35,] 130382.200 121429.400
[36,] 115722.267 130382.200
[37,] 159116.867 115722.267
[38,] 142079.667 159116.867
[39,] 87384.867 142079.667
[40,] 109598.467 87384.867
[41,] 97926.067 109598.467
[42,] 77494.267 97926.067
[43,] 57239.667 77494.267
[44,] 54644.467 57239.667
[45,] 41709.067 54644.467
[46,] -15600.933 41709.067
[47,] -12078.133 -15600.933
[48,] -92283.067 -12078.133
[49,] -111250.467 -92283.067
[50,] -103762.667 -111250.467
[51,] -28302.467 -103762.667
[52,] -87445.867 -28302.467
[53,] -87853.267 -87445.867
[54,] -79295.067 -87853.267
[55,] -51443.667 -79295.067
[56,] -108493.867 -51443.667
[57,] -78722.267 -108493.867
[58,] -63148.267 -78722.267
[59,] -54006.467 -63148.267
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -15866.133 -4524.733
2 -43269.333 -15866.133
3 -65428.133 -43269.333
4 -46902.533 -65428.133
5 -38951.933 -46902.533
6 -24412.733 -38951.933
7 -24705.333 -24412.733
8 -17886.533 -24705.333
9 -41311.933 -17886.533
10 -26085.933 -41311.933
11 -37411.133 -26085.933
12 -7076.067 -37411.133
13 -13758.467 -7076.067
14 -21705.667 -13758.467
15 6455.533 -21705.667
16 6923.133 6455.533
17 -17251.267 6923.133
18 -29009.067 -17251.267
19 -44306.667 -29009.067
20 -16272.867 -44306.667
21 -23783.267 -16272.867
22 -16594.267 -23783.267
23 -26886.467 -16594.267
24 -11838.400 -26886.467
25 -18241.800 -11838.400
26 26658.000 -18241.800
27 -109.800 26658.000
28 17826.800 -109.800
29 46130.400 17826.800
30 55222.600 46130.400
31 63216.000 55222.600
32 88008.800 63216.000
33 102108.400 88008.800
34 121429.400 102108.400
35 130382.200 121429.400
36 115722.267 130382.200
37 159116.867 115722.267
38 142079.667 159116.867
39 87384.867 142079.667
40 109598.467 87384.867
41 97926.067 109598.467
42 77494.267 97926.067
43 57239.667 77494.267
44 54644.467 57239.667
45 41709.067 54644.467
46 -15600.933 41709.067
47 -12078.133 -15600.933
48 -92283.067 -12078.133
49 -111250.467 -92283.067
50 -103762.667 -111250.467
51 -28302.467 -103762.667
52 -87445.867 -28302.467
53 -87853.267 -87445.867
54 -79295.067 -87853.267
55 -51443.667 -79295.067
56 -108493.867 -51443.667
57 -78722.267 -108493.867
58 -63148.267 -78722.267
59 -54006.467 -63148.267
> 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/rcomp/tmp/77feq1292683443.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/rcomp/tmp/87feq1292683443.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/rcomp/tmp/97feq1292683443.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/rcomp/tmp/100ovt1292683443.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/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/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/rcomp/tmp/1146cz1292683443.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/rcomp/tmp/12ca0b1292683443.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/rcomp/tmp/13lhqe1292683443.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/rcomp/tmp/14pzpk1292683443.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/rcomp/tmp/15ai5p1292683443.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/rcomp/tmp/16vimv1292683443.tab")
+ }
>
> try(system("convert tmp/1b5gh1292683443.ps tmp/1b5gh1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/2b5gh1292683443.ps tmp/2b5gh1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/3mwyk1292683443.ps tmp/3mwyk1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/4mwyk1292683443.ps tmp/4mwyk1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/5mwyk1292683443.ps tmp/5mwyk1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/6xnxn1292683443.ps tmp/6xnxn1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/77feq1292683443.ps tmp/77feq1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/87feq1292683443.ps tmp/87feq1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/97feq1292683443.ps tmp/97feq1292683443.png",intern=TRUE))
character(0)
> try(system("convert tmp/100ovt1292683443.ps tmp/100ovt1292683443.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
2.920 1.760 4.724