R version 3.3.3 (2017-03-06) -- "Another Canoe"
Copyright (C) 2017 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-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(3.4
+ ,0.4
+ ,0.6
+ ,-.5
+ ,1.126
+ ,1.682
+ ,1.838
+ ,1.2
+ ,-.3
+ ,0.4
+ ,0.6
+ ,1.090
+ ,1.126
+ ,1.682
+ ,3.2
+ ,0.1
+ ,-.3
+ ,0.4
+ ,1.369
+ ,1.090
+ ,1.126
+ ,3.2
+ ,2.2
+ ,0.1
+ ,-.3
+ ,1.666
+ ,1.369
+ ,1.090
+ ,7.6
+ ,3.4
+ ,2.2
+ ,0.1
+ ,2.339
+ ,1.666
+ ,1.369
+ ,10.3
+ ,2.4
+ ,3.4
+ ,2.2
+ ,1.527
+ ,2.339
+ ,1.666
+ ,6.5
+ ,4.4
+ ,2.4
+ ,3.4
+ ,1.573
+ ,1.527
+ ,2.339
+ ,9.4
+ ,6.9
+ ,4.4
+ ,2.4
+ ,0.786
+ ,1.573
+ ,1.527
+ ,3.0
+ ,6.4
+ ,6.9
+ ,4.4
+ ,2.268
+ ,0.786
+ ,1.573
+ ,11.9
+ ,14.5
+ ,6.4
+ ,6.9
+ ,2.540
+ ,2.268
+ ,0.786
+ ,10.0
+ ,14.8
+ ,14.5
+ ,6.4
+ ,3.734
+ ,2.540
+ ,2.268
+ ,14.0
+ ,9.8
+ ,14.8
+ ,14.5
+ ,4.477
+ ,3.734
+ ,2.540
+ ,5.0
+ ,7.3
+ ,9.8
+ ,14.8
+ ,5.446
+ ,4.477
+ ,3.734
+ ,8.7
+ ,3.4
+ ,7.3
+ ,9.8
+ ,3.070
+ ,5.446
+ ,4.477
+ ,3.6
+ ,0.2
+ ,3.4
+ ,7.3
+ ,5.361
+ ,3.070
+ ,5.446
+ ,6.4
+ ,3.5
+ ,0.2
+ ,3.4
+ ,0.802
+ ,5.361
+ ,3.070
+ ,4.3
+ ,2.7
+ ,3.5
+ ,0.2
+ ,2.858
+ ,0.802
+ ,5.361
+ ,3.8
+ ,1.3
+ ,2.7
+ ,3.5
+ ,3.779
+ ,2.858
+ ,0.802
+ ,4.7
+ ,0.5
+ ,1.3
+ ,2.7
+ ,2.869
+ ,3.779
+ ,2.858
+ ,2.4
+ ,1.1
+ ,0.5
+ ,1.3
+ ,5.135
+ ,2.869
+ ,3.779
+ ,4.7
+ ,0.4
+ ,1.1
+ ,0.5
+ ,3.676
+ ,5.135
+ ,2.869
+ ,2.5
+ ,-1.3
+ ,0.4
+ ,1.1
+ ,2.546
+ ,3.676
+ ,5.135
+ ,-11.5
+ ,-17.4
+ ,-1.3
+ ,0.4
+ ,2.297
+ ,2.546
+ ,3.676
+ ,-15.1
+ ,-25.3
+ ,-17.4
+ ,-1.3
+ ,2.252
+ ,2.297
+ ,2.546
+ ,-0.4
+ ,-20.4
+ ,-25.3
+ ,-17.4
+ ,-.020
+ ,2.252
+ ,2.297
+ ,8.3
+ ,-6.0
+ ,-20.4
+ ,-25.3
+ ,3.307
+ ,-.020
+ ,2.252
+ ,10.2
+ ,-2.0
+ ,-6.0
+ ,-20.4
+ ,1.450
+ ,3.307
+ ,-.020
+ ,4.9
+ ,0.1
+ ,-2.0
+ ,-6.0
+ ,0.378
+ ,1.450
+ ,3.307)
+ ,dim=c(7
+ ,28)
+ ,dimnames=list(c('Q'
+ ,'G'
+ ,'G-1'
+ ,'G-2'
+ ,'M'
+ ,'M-1'
+ ,'M-2')
+ ,1:28))
> y <- array(NA,dim=c(7,28),dimnames=list(c('Q','G','G-1','G-2','M','M-1','M-2'),1:28))
> for (i in 1:dim(x)[1])
+ {
+ for (j in 1:dim(x)[2])
+ {
+ y[i,j] <- as.numeric(x[i,j])
+ }
+ }
> par5 = ''
> par4 = ''
> par3 = 'No Linear Trend'
> par2 = 'Do not include Seasonal Dummies'
> par1 = '1'
> par5 <- ''
> par4 <- ''
> par3 <- 'No Linear Trend'
> par2 <- 'Do not include Seasonal Dummies'
> par1 <- '1'
> #'GNU S' R Code compiled by R2WASP v. 1.2.327 (Wed, 08 Jun 2016 16:18:16 +0100)
> #Author: root
> #To cite this work: Wessa P., (2015), Multiple Regression (v1.0.38) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_multipleregression.wasp/
> #Source of accompanying publication: Office for Research, Development, and Education
> #
> library(lattice)
> library(lmtest)
Loading required package: zoo
Attaching package: 'zoo'
The following objects are masked from 'package:base':
as.Date, as.Date.numeric
> n25 <- 25 #minimum number of obs. for Goldfeld-Quandt test
> mywarning <- ''
> par1 <- as.numeric(par1)
> if(is.na(par1)) {
+ par1 <- 1
+ mywarning = 'Warning: you did not specify the column number of the endogenous series! The first column was selected by default.'
+ }
> if (par4=='') par4 <- 0
> par4 <- as.numeric(par4)
> if (par5=='') par5 <- 0
> par5 <- as.numeric(par5)
> x <- na.omit(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'){
+ (n <- n -1)
+ x2 <- array(0, dim=c(n,k), dimnames=list(1:n, paste('(1-B)',colnames(x),sep='')))
+ for (i in 1:n) {
+ for (j in 1:k) {
+ x2[i,j] <- x[i+1,j] - x[i,j]
+ }
+ }
+ x <- x2
+ }
> if (par3 == 'Seasonal Differences (s=12)'){
+ (n <- n - 12)
+ x2 <- array(0, dim=c(n,k), dimnames=list(1:n, paste('(1-B12)',colnames(x),sep='')))
+ for (i in 1:n) {
+ for (j in 1:k) {
+ x2[i,j] <- x[i+12,j] - x[i,j]
+ }
+ }
+ x <- x2
+ }
> if (par3 == 'First and Seasonal Differences (s=12)'){
+ (n <- n -1)
+ x2 <- array(0, dim=c(n,k), dimnames=list(1:n, paste('(1-B)',colnames(x),sep='')))
+ for (i in 1:n) {
+ for (j in 1:k) {
+ x2[i,j] <- x[i+1,j] - x[i,j]
+ }
+ }
+ x <- x2
+ (n <- n - 12)
+ x2 <- array(0, dim=c(n,k), dimnames=list(1:n, paste('(1-B12)',colnames(x),sep='')))
+ for (i in 1:n) {
+ for (j in 1:k) {
+ x2[i,j] <- x[i+12,j] - x[i,j]
+ }
+ }
+ x <- x2
+ }
> if(par4 > 0) {
+ x2 <- array(0, dim=c(n-par4,par4), dimnames=list(1:(n-par4), paste(colnames(x)[par1],'(t-',1:par4,')',sep='')))
+ for (i in 1:(n-par4)) {
+ for (j in 1:par4) {
+ x2[i,j] <- x[i+par4-j,par1]
+ }
+ }
+ x <- cbind(x[(par4+1):n,], x2)
+ n <- n - par4
+ }
> if(par5 > 0) {
+ x2 <- array(0, dim=c(n-par5*12,par5), dimnames=list(1:(n-par5*12), paste(colnames(x)[par1],'(t-',1:par5,'s)',sep='')))
+ for (i in 1:(n-par5*12)) {
+ for (j in 1:par5) {
+ x2[i,j] <- x[i+par5*12-j*12,par1]
+ }
+ }
+ x <- cbind(x[(par5*12+1):n,], x2)
+ n <- n - par5*12
+ }
> 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[n,]))
[1] 7
> if (par3 == 'Linear Trend'){
+ x <- cbind(x, c(1:n))
+ colnames(x)[k+1] <- 't'
+ }
> x
Q G G-1 G-2 M M-1 M-2
1 3.4 0.4 0.6 -0.5 1.126 1.682 1.838
2 1.2 -0.3 0.4 0.6 1.090 1.126 1.682
3 3.2 0.1 -0.3 0.4 1.369 1.090 1.126
4 3.2 2.2 0.1 -0.3 1.666 1.369 1.090
5 7.6 3.4 2.2 0.1 2.339 1.666 1.369
6 10.3 2.4 3.4 2.2 1.527 2.339 1.666
7 6.5 4.4 2.4 3.4 1.573 1.527 2.339
8 9.4 6.9 4.4 2.4 0.786 1.573 1.527
9 3.0 6.4 6.9 4.4 2.268 0.786 1.573
10 11.9 14.5 6.4 6.9 2.540 2.268 0.786
11 10.0 14.8 14.5 6.4 3.734 2.540 2.268
12 14.0 9.8 14.8 14.5 4.477 3.734 2.540
13 5.0 7.3 9.8 14.8 5.446 4.477 3.734
14 8.7 3.4 7.3 9.8 3.070 5.446 4.477
15 3.6 0.2 3.4 7.3 5.361 3.070 5.446
16 6.4 3.5 0.2 3.4 0.802 5.361 3.070
17 4.3 2.7 3.5 0.2 2.858 0.802 5.361
18 3.8 1.3 2.7 3.5 3.779 2.858 0.802
19 4.7 0.5 1.3 2.7 2.869 3.779 2.858
20 2.4 1.1 0.5 1.3 5.135 2.869 3.779
21 4.7 0.4 1.1 0.5 3.676 5.135 2.869
22 2.5 -1.3 0.4 1.1 2.546 3.676 5.135
23 -11.5 -17.4 -1.3 0.4 2.297 2.546 3.676
24 -15.1 -25.3 -17.4 -1.3 2.252 2.297 2.546
25 -0.4 -20.4 -25.3 -17.4 -0.020 2.252 2.297
26 8.3 -6.0 -20.4 -25.3 3.307 -0.020 2.252
27 10.2 -2.0 -6.0 -20.4 1.450 3.307 -0.020
28 4.9 0.1 -2.0 -6.0 0.378 1.450 3.307
> (k <- length(x[n,]))
[1] 7
> head(x)
Q G G-1 G-2 M M-1 M-2
1 3.4 0.4 0.6 -0.5 1.126 1.682 1.838
2 1.2 -0.3 0.4 0.6 1.090 1.126 1.682
3 3.2 0.1 -0.3 0.4 1.369 1.090 1.126
4 3.2 2.2 0.1 -0.3 1.666 1.369 1.090
5 7.6 3.4 2.2 0.1 2.339 1.666 1.369
6 10.3 2.4 3.4 2.2 1.527 2.339 1.666
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) G `G-1` `G-2` M `M-1`
2.23672 0.83790 -0.18007 -0.26449 -0.10729 0.98873
`M-2`
-0.05166
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-2.8900 -1.4763 -0.5414 1.2667 6.9717
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.23672 1.47430 1.517 0.1441
G 0.83790 0.12402 6.756 1.1e-06 ***
`G-1` -0.18007 0.18288 -0.985 0.3360
`G-2` -0.26449 0.12585 -2.102 0.0478 *
M -0.10729 0.42544 -0.252 0.8033
`M-1` 0.98873 0.41241 2.397 0.0259 *
`M-2` -0.05166 0.43598 -0.119 0.9068
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 2.675 on 21 degrees of freedom
Multiple R-squared: 0.8508, Adjusted R-squared: 0.8081
F-statistic: 19.95 on 6 and 21 DF, p-value: 1.135e-07
> 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.3218229 0.64364571 0.67817714
[2,] 0.3173309 0.63466179 0.68266911
[3,] 0.5265439 0.94691215 0.47345607
[4,] 0.5866447 0.82671059 0.41335529
[5,] 0.7613896 0.47722084 0.23861042
[6,] 0.9584528 0.08309445 0.04154722
[7,] 0.9627407 0.07451867 0.03725933
[8,] 0.9406740 0.11865195 0.05932597
[9,] 0.8960902 0.20781968 0.10390984
> postscript(file="/var/wessaorg/rcomp/tmp/12l141495316369.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/2uvlm1495316369.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/3ohlx1495316369.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/4wkbq1495316369.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/5f7l61495316369.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 = 28
Frequency = 1
1 2 3 4 5 6
-0.64335677 -1.46409738 0.05859798 -2.05994649 1.61149632 5.18370906
7 8 9 10 11 12
0.68777690 1.41682536 -2.64555369 -1.43812289 -2.32739171 6.97167398
13 14 15 16 17 18
-1.32358397 2.69698273 1.55984494 -2.89003065 0.27476176 -0.49299209
19 20 21 22 23 24
-0.28840274 -2.41502770 -2.07606154 -1.38060773 -1.36655315 -1.51300446
25 26 27 28
3.18835687 1.21660026 0.04793148 -0.58982468
> postscript(file="/var/wessaorg/rcomp/tmp/6c66x1495316369.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 = 28
Frequency = 1
lag(myerror, k = 1) myerror
0 -0.64335677 NA
1 -1.46409738 -0.64335677
2 0.05859798 -1.46409738
3 -2.05994649 0.05859798
4 1.61149632 -2.05994649
5 5.18370906 1.61149632
6 0.68777690 5.18370906
7 1.41682536 0.68777690
8 -2.64555369 1.41682536
9 -1.43812289 -2.64555369
10 -2.32739171 -1.43812289
11 6.97167398 -2.32739171
12 -1.32358397 6.97167398
13 2.69698273 -1.32358397
14 1.55984494 2.69698273
15 -2.89003065 1.55984494
16 0.27476176 -2.89003065
17 -0.49299209 0.27476176
18 -0.28840274 -0.49299209
19 -2.41502770 -0.28840274
20 -2.07606154 -2.41502770
21 -1.38060773 -2.07606154
22 -1.36655315 -1.38060773
23 -1.51300446 -1.36655315
24 3.18835687 -1.51300446
25 1.21660026 3.18835687
26 0.04793148 1.21660026
27 -0.58982468 0.04793148
28 NA -0.58982468
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -1.46409738 -0.64335677
[2,] 0.05859798 -1.46409738
[3,] -2.05994649 0.05859798
[4,] 1.61149632 -2.05994649
[5,] 5.18370906 1.61149632
[6,] 0.68777690 5.18370906
[7,] 1.41682536 0.68777690
[8,] -2.64555369 1.41682536
[9,] -1.43812289 -2.64555369
[10,] -2.32739171 -1.43812289
[11,] 6.97167398 -2.32739171
[12,] -1.32358397 6.97167398
[13,] 2.69698273 -1.32358397
[14,] 1.55984494 2.69698273
[15,] -2.89003065 1.55984494
[16,] 0.27476176 -2.89003065
[17,] -0.49299209 0.27476176
[18,] -0.28840274 -0.49299209
[19,] -2.41502770 -0.28840274
[20,] -2.07606154 -2.41502770
[21,] -1.38060773 -2.07606154
[22,] -1.36655315 -1.38060773
[23,] -1.51300446 -1.36655315
[24,] 3.18835687 -1.51300446
[25,] 1.21660026 3.18835687
[26,] 0.04793148 1.21660026
[27,] -0.58982468 0.04793148
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -1.46409738 -0.64335677
2 0.05859798 -1.46409738
3 -2.05994649 0.05859798
4 1.61149632 -2.05994649
5 5.18370906 1.61149632
6 0.68777690 5.18370906
7 1.41682536 0.68777690
8 -2.64555369 1.41682536
9 -1.43812289 -2.64555369
10 -2.32739171 -1.43812289
11 6.97167398 -2.32739171
12 -1.32358397 6.97167398
13 2.69698273 -1.32358397
14 1.55984494 2.69698273
15 -2.89003065 1.55984494
16 0.27476176 -2.89003065
17 -0.49299209 0.27476176
18 -0.28840274 -0.49299209
19 -2.41502770 -0.28840274
20 -2.07606154 -2.41502770
21 -1.38060773 -2.07606154
22 -1.36655315 -1.38060773
23 -1.51300446 -1.36655315
24 3.18835687 -1.51300446
25 1.21660026 3.18835687
26 0.04793148 1.21660026
27 -0.58982468 0.04793148
> 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/73gjh1495316369.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/8eznl1495316369.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/9c22p1495316369.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/10dy3p1495316369.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, signif(mysum$coefficients[i,1],6), 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.row.start(a)
> a<-table.element(a, mywarning)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/11e99c1495316369.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,formatC(signif(mysum$coefficients[i,1],5),format='g',flag='+'))
+ a<-table.element(a,formatC(signif(mysum$coefficients[i,2],5),format='g',flag=' '))
+ a<-table.element(a,formatC(signif(mysum$coefficients[i,3],4),format='e',flag='+'))
+ a<-table.element(a,formatC(signif(mysum$coefficients[i,4],4),format='g',flag=' '))
+ a<-table.element(a,formatC(signif(mysum$coefficients[i,4]/2,4),format='g',flag=' '))
+ a<-table.row.end(a)
+ }
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/1220bf1495316369.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,formatC(signif(sqrt(mysum$r.squared),6),format='g',flag=' '))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'R-squared',1,TRUE)
> a<-table.element(a,formatC(signif(mysum$r.squared,6),format='g',flag=' '))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Adjusted R-squared',1,TRUE)
> a<-table.element(a,formatC(signif(mysum$adj.r.squared,6),format='g',flag=' '))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (value)',1,TRUE)
> a<-table.element(a,formatC(signif(mysum$fstatistic[1],6),format='g',flag=' '))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF numerator)',1,TRUE)
> a<-table.element(a, signif(mysum$fstatistic[2],6))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'F-TEST (DF denominator)',1,TRUE)
> a<-table.element(a, signif(mysum$fstatistic[3],6))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'p-value',1,TRUE)
> a<-table.element(a,formatC(signif(1-pf(mysum$fstatistic[1],mysum$fstatistic[2],mysum$fstatistic[3]),6),format='g',flag=' '))
> 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,formatC(signif(mysum$sigma,6),format='g',flag=' '))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a, 'Sum Squared Residuals',1,TRUE)
> a<-table.element(a,formatC(signif(sum(myerror*myerror),6),format='g',flag=' '))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/wessaorg/rcomp/tmp/137tyc1495316369.tab")
> if(n < 200) {
+ 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,formatC(signif(x[i],6),format='g',flag=' '))
+ a<-table.element(a,formatC(signif(x[i]-mysum$resid[i],6),format='g',flag=' '))
+ a<-table.element(a,formatC(signif(mysum$resid[i],6),format='g',flag=' '))
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/14k40j1495316369.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,formatC(signif(gqarr[mypoint-kp3+1,1],6),format='g',flag=' '))
+ a<-table.element(a,formatC(signif(gqarr[mypoint-kp3+1,2],6),format='g',flag=' '))
+ a<-table.element(a,formatC(signif(gqarr[mypoint-kp3+1,3],6),format='g',flag=' '))
+ a<-table.row.end(a)
+ }
+ a<-table.end(a)
+ table.save(a,file="/var/wessaorg/rcomp/tmp/157zs81495316369.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,signif(numsignificant1,6))
+ a<-table.element(a,formatC(signif(numsignificant1/numgqtests,6),format='g',flag=' '))
+ 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,signif(numsignificant5,6))
+ a<-table.element(a,signif(numsignificant5/numgqtests,6))
+ 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,signif(numsignificant10,6))
+ a<-table.element(a,signif(numsignificant10/numgqtests,6))
+ 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/164kwi1495316369.tab")
+ }
+ }
>
> try(system("convert tmp/12l141495316369.ps tmp/12l141495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/2uvlm1495316369.ps tmp/2uvlm1495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/3ohlx1495316369.ps tmp/3ohlx1495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/4wkbq1495316369.ps tmp/4wkbq1495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/5f7l61495316369.ps tmp/5f7l61495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/6c66x1495316369.ps tmp/6c66x1495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/73gjh1495316369.ps tmp/73gjh1495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/8eznl1495316369.ps tmp/8eznl1495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/9c22p1495316369.ps tmp/9c22p1495316369.png",intern=TRUE))
character(0)
> try(system("convert tmp/10dy3p1495316369.ps tmp/10dy3p1495316369.png",intern=TRUE))
character(0)
>
>
> proc.time()
user system elapsed
5.109 0.635 5.874