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.
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(12300.00,0.00,12092.80,0.00,12380.80,0.00,12196.90,0.00,9455.00,0.00,13168.00,0.00,13427.90,0.00,11980.50,0.00,11884.80,0.00,11691.70,0.00,12233.80,0.00,14341.40,0.00,13130.70,0.00,12421.10,0.00,14285.80,0.00,12864.60,0.00,11160.20,0.00,14316.20,0.00,14388.70,0.00,14013.90,0.00,13419.00,0.00,12769.60,0.00,13315.50,0.00,15332.90,0.00,14243.00,0.00,13824.40,0.00,14962.90,0.00,13202.90,0.00,12199.00,0.00,15508.90,0.00,14199.80,0.00,15169.60,0.00,14058.00,0.00,13786.20,0.00,14147.90,0.00,16541.70,0.00,13587.50,0.00,15582.40,0.00,15802.80,0.00,14130.50,0.00,12923.20,0.00,15612.20,1.00,16033.70,1.00,16036.60,1.00,14037.80,1.00,15330.60,1.00,15038.30,1.00,17401.80,1.00,14992.50,1.00,16043.70,1.00,16929.60,1.00,15921.30,1.00,14417.20,1.00,15961.00,1.00,17851.90,1.00,16483.90,1.00,14215.50,1.00,17429.70,1.00,17839.50,1.00,17629.20,1.00),dim=c(2,60),dimnames=list(c('x','y'),1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('x','y'),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 = '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
x y
1 12300.0 0
2 12092.8 0
3 12380.8 0
4 12196.9 0
5 9455.0 0
6 13168.0 0
7 13427.9 0
8 11980.5 0
9 11884.8 0
10 11691.7 0
11 12233.8 0
12 14341.4 0
13 13130.7 0
14 12421.1 0
15 14285.8 0
16 12864.6 0
17 11160.2 0
18 14316.2 0
19 14388.7 0
20 14013.9 0
21 13419.0 0
22 12769.6 0
23 13315.5 0
24 15332.9 0
25 14243.0 0
26 13824.4 0
27 14962.9 0
28 13202.9 0
29 12199.0 0
30 15508.9 0
31 14199.8 0
32 15169.6 0
33 14058.0 0
34 13786.2 0
35 14147.9 0
36 16541.7 0
37 13587.5 0
38 15582.4 0
39 15802.8 0
40 14130.5 0
41 12923.2 0
42 15612.2 1
43 16033.7 1
44 16036.6 1
45 14037.8 1
46 15330.6 1
47 15038.3 1
48 17401.8 1
49 14992.5 1
50 16043.7 1
51 16929.6 1
52 15921.3 1
53 14417.2 1
54 15961.0 1
55 17851.9 1
56 16483.9 1
57 14215.5 1
58 17429.7 1
59 17839.5 1
60 17629.2 1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) y
13474 2589
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-4019.21 -1057.57 -38.04 848.03 3067.49
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 13474.2 210.8 63.913 < 2e-16 ***
y 2589.3 374.6 6.911 4.13e-09 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1350 on 58 degrees of freedom
Multiple R-squared: 0.4516, Adjusted R-squared: 0.4422
F-statistic: 47.77 on 1 and 58 DF, p-value: 4.129e-09
> 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.7358103 0.5283794 0.2641897
[2,] 0.7182314 0.5635372 0.2817686
[3,] 0.7063221 0.5873558 0.2936779
[4,] 0.6119197 0.7761605 0.3880803
[5,] 0.5274907 0.9450185 0.4725093
[6,] 0.4689176 0.9378352 0.5310824
[7,] 0.3931699 0.7863399 0.6068301
[8,] 0.5787939 0.8424122 0.4212061
[9,] 0.5264234 0.9471532 0.4735766
[10,] 0.4647840 0.9295679 0.5352160
[11,] 0.5485933 0.9028133 0.4514067
[12,] 0.4866766 0.9733532 0.5133234
[13,] 0.6202146 0.7595707 0.3797854
[14,] 0.6746846 0.6506308 0.3253154
[15,] 0.7109468 0.5781065 0.2890532
[16,] 0.6961935 0.6076131 0.3038065
[17,] 0.6505447 0.6989106 0.3494553
[18,] 0.6199364 0.7601271 0.3800636
[19,] 0.5767099 0.8465801 0.4232901
[20,] 0.7043979 0.5912042 0.2956021
[21,] 0.6806926 0.6386148 0.3193074
[22,] 0.6352066 0.7295868 0.3647934
[23,] 0.6625395 0.6749209 0.3374605
[24,] 0.6217663 0.7564675 0.3782337
[25,] 0.6900776 0.6198448 0.3099224
[26,] 0.7580732 0.4838535 0.2419268
[27,] 0.7184875 0.5630249 0.2815125
[28,] 0.7292732 0.5414535 0.2707268
[29,] 0.6797662 0.6404676 0.3202338
[30,] 0.6312762 0.7374476 0.3687238
[31,] 0.5791839 0.8416321 0.4208161
[32,] 0.7599527 0.4800945 0.2400473
[33,] 0.7165443 0.5669113 0.2834557
[34,] 0.7385584 0.5228832 0.2614416
[35,] 0.8150415 0.3699169 0.1849585
[36,] 0.7774838 0.4450325 0.2225162
[37,] 0.7126224 0.5747553 0.2873776
[38,] 0.6425135 0.7149730 0.3574865
[39,] 0.5597272 0.8805455 0.4402728
[40,] 0.4721996 0.9443992 0.5278004
[41,] 0.5754391 0.8491218 0.4245609
[42,] 0.5172111 0.9655778 0.4827889
[43,] 0.4929555 0.9859111 0.5070445
[44,] 0.4758557 0.9517113 0.5241443
[45,] 0.4612169 0.9224338 0.5387831
[46,] 0.3652412 0.7304825 0.6347588
[47,] 0.2856883 0.5713767 0.7143117
[48,] 0.2032830 0.4065660 0.7967170
[49,] 0.3101484 0.6202968 0.6898516
[50,] 0.2307296 0.4614593 0.7692704
[51,] 0.1829473 0.3658945 0.8170527
> postscript(file="/var/www/html/freestat/rcomp/tmp/1unei1227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(x[,1], type='l', main='Actuals and Interpolation', ylab='value of Actuals and Interpolation (dots)', xlab='time or index')
> points(x[,1]-mysum$resid)
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/248xh1227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> plot(mysum$resid, type='b', pch=19, main='Residuals', ylab='value of Residuals', xlab='time or index')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/31th01227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> hist(mysum$resid, main='Residual Histogram', xlab='values of Residuals')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/471p91227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> densityplot(~mysum$resid,col='black',main='Residual Density Plot', xlab='values of Residuals')
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/5oyk71227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> qqnorm(mysum$resid, main='Residual Normal Q-Q Plot')
> qqline(mysum$resid)
> grid()
> dev.off()
null device
1
> (myerror <- as.ts(mysum$resid))
Time Series:
Start = 1
End = 60
Frequency = 1
1 2 3 4 5 6
-1174.20732 -1381.40732 -1093.40732 -1277.30732 -4019.20732 -306.20732
7 8 9 10 11 12
-46.30732 -1493.70732 -1589.40732 -1782.50732 -1240.40732 867.19268
13 14 15 16 17 18
-343.50732 -1053.10732 811.59268 -609.60732 -2314.00732 841.99268
19 20 21 22 23 24
914.49268 539.69268 -55.20732 -704.60732 -158.70732 1858.69268
25 26 27 28 29 30
768.79268 350.19268 1488.69268 -271.30732 -1275.20732 2034.69268
31 32 33 34 35 36
725.59268 1695.39268 583.79268 311.99268 673.69268 3067.49268
37 38 39 40 41 42
113.29268 2108.19268 2328.59268 656.29268 -551.00732 -451.27368
43 44 45 46 47 48
-29.77368 -26.87368 -2025.67368 -732.87368 -1025.17368 1338.32632
49 50 51 52 53 54
-1070.97368 -19.77368 866.12632 -142.17368 -1646.27368 -102.47368
55 56 57 58 59 60
1788.42632 420.42632 -1847.97368 1366.22632 1776.02632 1565.72632
> postscript(file="/var/www/html/freestat/rcomp/tmp/6wp2r1227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> dum <- cbind(lag(myerror,k=1),myerror)
> dum
Time Series:
Start = 0
End = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 -1174.20732 NA
1 -1381.40732 -1174.20732
2 -1093.40732 -1381.40732
3 -1277.30732 -1093.40732
4 -4019.20732 -1277.30732
5 -306.20732 -4019.20732
6 -46.30732 -306.20732
7 -1493.70732 -46.30732
8 -1589.40732 -1493.70732
9 -1782.50732 -1589.40732
10 -1240.40732 -1782.50732
11 867.19268 -1240.40732
12 -343.50732 867.19268
13 -1053.10732 -343.50732
14 811.59268 -1053.10732
15 -609.60732 811.59268
16 -2314.00732 -609.60732
17 841.99268 -2314.00732
18 914.49268 841.99268
19 539.69268 914.49268
20 -55.20732 539.69268
21 -704.60732 -55.20732
22 -158.70732 -704.60732
23 1858.69268 -158.70732
24 768.79268 1858.69268
25 350.19268 768.79268
26 1488.69268 350.19268
27 -271.30732 1488.69268
28 -1275.20732 -271.30732
29 2034.69268 -1275.20732
30 725.59268 2034.69268
31 1695.39268 725.59268
32 583.79268 1695.39268
33 311.99268 583.79268
34 673.69268 311.99268
35 3067.49268 673.69268
36 113.29268 3067.49268
37 2108.19268 113.29268
38 2328.59268 2108.19268
39 656.29268 2328.59268
40 -551.00732 656.29268
41 -451.27368 -551.00732
42 -29.77368 -451.27368
43 -26.87368 -29.77368
44 -2025.67368 -26.87368
45 -732.87368 -2025.67368
46 -1025.17368 -732.87368
47 1338.32632 -1025.17368
48 -1070.97368 1338.32632
49 -19.77368 -1070.97368
50 866.12632 -19.77368
51 -142.17368 866.12632
52 -1646.27368 -142.17368
53 -102.47368 -1646.27368
54 1788.42632 -102.47368
55 420.42632 1788.42632
56 -1847.97368 420.42632
57 1366.22632 -1847.97368
58 1776.02632 1366.22632
59 1565.72632 1776.02632
60 NA 1565.72632
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -1381.40732 -1174.20732
[2,] -1093.40732 -1381.40732
[3,] -1277.30732 -1093.40732
[4,] -4019.20732 -1277.30732
[5,] -306.20732 -4019.20732
[6,] -46.30732 -306.20732
[7,] -1493.70732 -46.30732
[8,] -1589.40732 -1493.70732
[9,] -1782.50732 -1589.40732
[10,] -1240.40732 -1782.50732
[11,] 867.19268 -1240.40732
[12,] -343.50732 867.19268
[13,] -1053.10732 -343.50732
[14,] 811.59268 -1053.10732
[15,] -609.60732 811.59268
[16,] -2314.00732 -609.60732
[17,] 841.99268 -2314.00732
[18,] 914.49268 841.99268
[19,] 539.69268 914.49268
[20,] -55.20732 539.69268
[21,] -704.60732 -55.20732
[22,] -158.70732 -704.60732
[23,] 1858.69268 -158.70732
[24,] 768.79268 1858.69268
[25,] 350.19268 768.79268
[26,] 1488.69268 350.19268
[27,] -271.30732 1488.69268
[28,] -1275.20732 -271.30732
[29,] 2034.69268 -1275.20732
[30,] 725.59268 2034.69268
[31,] 1695.39268 725.59268
[32,] 583.79268 1695.39268
[33,] 311.99268 583.79268
[34,] 673.69268 311.99268
[35,] 3067.49268 673.69268
[36,] 113.29268 3067.49268
[37,] 2108.19268 113.29268
[38,] 2328.59268 2108.19268
[39,] 656.29268 2328.59268
[40,] -551.00732 656.29268
[41,] -451.27368 -551.00732
[42,] -29.77368 -451.27368
[43,] -26.87368 -29.77368
[44,] -2025.67368 -26.87368
[45,] -732.87368 -2025.67368
[46,] -1025.17368 -732.87368
[47,] 1338.32632 -1025.17368
[48,] -1070.97368 1338.32632
[49,] -19.77368 -1070.97368
[50,] 866.12632 -19.77368
[51,] -142.17368 866.12632
[52,] -1646.27368 -142.17368
[53,] -102.47368 -1646.27368
[54,] 1788.42632 -102.47368
[55,] 420.42632 1788.42632
[56,] -1847.97368 420.42632
[57,] 1366.22632 -1847.97368
[58,] 1776.02632 1366.22632
[59,] 1565.72632 1776.02632
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -1381.40732 -1174.20732
2 -1093.40732 -1381.40732
3 -1277.30732 -1093.40732
4 -4019.20732 -1277.30732
5 -306.20732 -4019.20732
6 -46.30732 -306.20732
7 -1493.70732 -46.30732
8 -1589.40732 -1493.70732
9 -1782.50732 -1589.40732
10 -1240.40732 -1782.50732
11 867.19268 -1240.40732
12 -343.50732 867.19268
13 -1053.10732 -343.50732
14 811.59268 -1053.10732
15 -609.60732 811.59268
16 -2314.00732 -609.60732
17 841.99268 -2314.00732
18 914.49268 841.99268
19 539.69268 914.49268
20 -55.20732 539.69268
21 -704.60732 -55.20732
22 -158.70732 -704.60732
23 1858.69268 -158.70732
24 768.79268 1858.69268
25 350.19268 768.79268
26 1488.69268 350.19268
27 -271.30732 1488.69268
28 -1275.20732 -271.30732
29 2034.69268 -1275.20732
30 725.59268 2034.69268
31 1695.39268 725.59268
32 583.79268 1695.39268
33 311.99268 583.79268
34 673.69268 311.99268
35 3067.49268 673.69268
36 113.29268 3067.49268
37 2108.19268 113.29268
38 2328.59268 2108.19268
39 656.29268 2328.59268
40 -551.00732 656.29268
41 -451.27368 -551.00732
42 -29.77368 -451.27368
43 -26.87368 -29.77368
44 -2025.67368 -26.87368
45 -732.87368 -2025.67368
46 -1025.17368 -732.87368
47 1338.32632 -1025.17368
48 -1070.97368 1338.32632
49 -19.77368 -1070.97368
50 866.12632 -19.77368
51 -142.17368 866.12632
52 -1646.27368 -142.17368
53 -102.47368 -1646.27368
54 1788.42632 -102.47368
55 420.42632 1788.42632
56 -1847.97368 420.42632
57 1366.22632 -1847.97368
58 1776.02632 1366.22632
59 1565.72632 1776.02632
> 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/70zuh1227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> acf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/83ouj1227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> pacf(mysum$resid, lag.max=length(mysum$resid)/2, main='Residual Partial Autocorrelation Function')
> grid()
> dev.off()
null device
1
> postscript(file="/var/www/html/freestat/rcomp/tmp/9gv5v1227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
> plot(mylm, las = 1, sub='Residual Diagnostics')
> par(opar)
> dev.off()
null device
1
> if (n > n25) {
+ postscript(file="/var/www/html/freestat/rcomp/tmp/104kht1227550834.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
+ plot(kp3:nmkm3,gqarr[,2], main='Goldfeld-Quandt test',ylab='2-sided p-value',xlab='breakpoint')
+ grid()
+ dev.off()
+ }
null device
1
>
> #Note: the /var/www/html/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/118bue1227550834.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/12087f1227550834.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/13qorb1227550834.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/14urxl1227550834.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/15kp0i1227550834.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/16t9ib1227550834.tab")
+ }
>
> system("convert tmp/1unei1227550834.ps tmp/1unei1227550834.png")
> system("convert tmp/248xh1227550834.ps tmp/248xh1227550834.png")
> system("convert tmp/31th01227550834.ps tmp/31th01227550834.png")
> system("convert tmp/471p91227550834.ps tmp/471p91227550834.png")
> system("convert tmp/5oyk71227550834.ps tmp/5oyk71227550834.png")
> system("convert tmp/6wp2r1227550834.ps tmp/6wp2r1227550834.png")
> system("convert tmp/70zuh1227550834.ps tmp/70zuh1227550834.png")
> system("convert tmp/83ouj1227550834.ps tmp/83ouj1227550834.png")
> system("convert tmp/9gv5v1227550834.ps tmp/9gv5v1227550834.png")
> system("convert tmp/104kht1227550834.ps tmp/104kht1227550834.png")
>
>
> proc.time()
user system elapsed
3.709 2.531 4.064