R version 3.3.2 (2016-10-31) -- "Sincere Pumpkin Patch"
Copyright (C) 2016 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.
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.
> source('/home/pw/wessanet/cretab')
>
>
> RC.capture <- function (expression, collapse = NULL) {
+ resultConn <- textConnection('RC.resultText', open = 'w', local=TRUE)
+ sink(resultConn)
+ on.exit(function() {
+ sink()
+ close(resultConn)
+ })
+ expression
+ on.exit(NULL)
+ sink()
+ close(resultConn)
+ return(paste(c(RC.resultText, ''), collapse = collapse, sep = ''))
+ }
> RC.texteval <- function (sourceText, collapse = NULL, echo = TRUE) {
+ sourceConn <- textConnection(sourceText, open = 'r')
+ on.exit(close(sourceConn))
+ result <- RC.capture(source(file = sourceConn, local = FALSE, echo = echo, print.eval = TRUE), collapse = collapse)
+ on.exit(NULL)
+ close(sourceConn)
+ res <- ''
+ for(i in 1:length(result)) {
+ if (result[i]!='') res <- paste(res,result[i],'
+ ',sep='')
+ }
+ return(res)
+ }
>
>
> myrfcuid = 'r0597336'
>
> x <- array(list(6,1,0,0,0,3.2,3.2,1,10.24,7,0,1,0,1,3.3,0,0,10.89,2,0,1,1,1,3,3,1,9,11,0,1,0,1,3.5,0,0,12.25,13,0,1,0,0,3.7,3.7,1,13.69,3,1,0,0,0,2.7,0,0,7.29,17,0,1,1,1,3.6,3.6,1,12.96,10,0,1,0,1,3.5,0,0,12.25,4,1,0,0,0,3.8,3.8,1,14.44,12,0,1,0,0,3.4,0,0,11.56,7,0,0,0,1,3.7,3.7,1,13.69,11,0,1,0,0,3.5,0,0,12.25,3,0,0,1,0,2.8,2.8,1,7.84,5,1,0,1,0,3.8,0,0,14.44,1,0,1,0,0,4.3,4.3,1,18.49,12,0,0,0,1,3.3,0,0,10.89,18,0,0,0,0,3.6,3.6,1,12.96,8,1,0,1,0,3.6,0,0,12.96,6,1,1,0,0,3.3,3.3,1,10.89,1,0,0,0,0,2.8,0,0,7.84),dim=c(9,20),dimnames=list(c('Score','X1','X2','X3','X4','X5','Inter','Geslacht','X6'),1:20))
> y <- array(NA,dim=c(9,20),dimnames=list(c('Score','X1','X2','X3','X4','X5','Inter','Geslacht','X6'),1:20))
> 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 (Mon, 12 Sep 2016 10:14:58 +0200)
> #Author: root
> #To cite this work: Wessa P., (2016), Multiple Regression (v1.0.40) 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
> library(car)
> library(MASS)
> 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] 9
> if (par3 == 'Linear Trend'){
+ x <- cbind(x, c(1:n))
+ colnames(x)[k+1] <- 't'
+ }
> print(x)
Score X1 X2 X3 X4 X5 Inter Geslacht X6
1 6 1 0 0 0 3.2 3.2 1 10.24
2 7 0 1 0 1 3.3 0.0 0 10.89
3 2 0 1 1 1 3.0 3.0 1 9.00
4 11 0 1 0 1 3.5 0.0 0 12.25
5 13 0 1 0 0 3.7 3.7 1 13.69
6 3 1 0 0 0 2.7 0.0 0 7.29
7 17 0 1 1 1 3.6 3.6 1 12.96
8 10 0 1 0 1 3.5 0.0 0 12.25
9 4 1 0 0 0 3.8 3.8 1 14.44
10 12 0 1 0 0 3.4 0.0 0 11.56
11 7 0 0 0 1 3.7 3.7 1 13.69
12 11 0 1 0 0 3.5 0.0 0 12.25
13 3 0 0 1 0 2.8 2.8 1 7.84
14 5 1 0 1 0 3.8 0.0 0 14.44
15 1 0 1 0 0 4.3 4.3 1 18.49
16 12 0 0 0 1 3.3 0.0 0 10.89
17 18 0 0 0 0 3.6 3.6 1 12.96
18 8 1 0 1 0 3.6 0.0 0 12.96
19 6 1 1 0 0 3.3 3.3 1 10.89
20 1 0 0 0 0 2.8 0.0 0 7.84
> (k <- length(x[n,]))
[1] 9
> head(x)
Score X1 X2 X3 X4 X5 Inter Geslacht X6
1 6 1 0 0 0 3.2 3.2 1 10.24
2 7 0 1 0 1 3.3 0.0 0 10.89
3 2 0 1 1 1 3.0 3.0 1 9.00
4 11 0 1 0 1 3.5 0.0 0 12.25
5 13 0 1 0 0 3.7 3.7 1 13.69
6 3 1 0 0 0 2.7 0.0 0 7.29
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X1 X2 X3 X4 X5
-257.7942 -4.8654 -0.2172 3.2287 -2.9267 163.2222
Inter Geslacht X6
13.1261 -44.4619 -24.6204
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-5.2041 -1.4614 -0.1287 1.5080 4.4825
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -257.7942 57.7528 -4.464 0.000957 ***
X1 -4.8654 2.1552 -2.258 0.045288 *
X2 -0.2172 1.8490 -0.117 0.908594
X3 3.2287 2.6060 1.239 0.241156
X4 -2.9267 2.1233 -1.378 0.195464
X5 163.2222 36.7108 4.446 0.000985 ***
Inter 13.1261 7.3277 1.791 0.100767
Geslacht -44.4619 24.9491 -1.782 0.102322
X6 -24.6204 5.8039 -4.242 0.001384 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 3.502 on 11 degrees of freedom
Multiple R-squared: 0.7181, Adjusted R-squared: 0.5131
F-statistic: 3.502 on 8 and 11 DF, p-value: 0.02904
> 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
+ }
> postscript(file="/home/pw/wessanet/rcomp/tmp/1rs3a1485164468.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="/home/pw/wessanet/rcomp/tmp/2wcai1485164468.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="/home/pw/wessanet/rcomp/tmp/3vvrn1485164468.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> sresid <- studres(mylm)
> hist(sresid, freq=FALSE, main='Distribution of Studentized Residuals')
> xfit<-seq(min(sresid),max(sresid),length=40)
> yfit<-dnorm(xfit)
> lines(xfit, yfit)
> grid()
> dev.off()
null device
1
> postscript(file="/home/pw/wessanet/rcomp/tmp/4obeq1485164468.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="/home/pw/wessanet/rcomp/tmp/5sb541485164468.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556)
> qqPlot(mylm, main='QQ Plot')
> grid()
> dev.off()
null device
1
> (myerror <- as.ts(mysum$resid))
Time Series:
Start = 1
End = 20
Frequency = 1
1 2 3 4 5 6
0.91980226 -2.57912298 -3.28994446 2.26015366 0.03781393 4.44226886
7 8 9 10 11 12
3.39780336 1.26015366 -3.48355908 -0.33241665 -3.25268906 -0.66657454
13 14 15 16 17 18
1.27610744 -0.29512213 0.40667972 2.20364582 4.48253160 -1.08884421
19 20
-0.49454570 -5.20414151
> postscript(file="/home/pw/wessanet/rcomp/tmp/6nzip1485164468.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 = 20
Frequency = 1
lag(myerror, k = 1) myerror
0 0.91980226 NA
1 -2.57912298 0.91980226
2 -3.28994446 -2.57912298
3 2.26015366 -3.28994446
4 0.03781393 2.26015366
5 4.44226886 0.03781393
6 3.39780336 4.44226886
7 1.26015366 3.39780336
8 -3.48355908 1.26015366
9 -0.33241665 -3.48355908
10 -3.25268906 -0.33241665
11 -0.66657454 -3.25268906
12 1.27610744 -0.66657454
13 -0.29512213 1.27610744
14 0.40667972 -0.29512213
15 2.20364582 0.40667972
16 4.48253160 2.20364582
17 -1.08884421 4.48253160
18 -0.49454570 -1.08884421
19 -5.20414151 -0.49454570
20 NA -5.20414151
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -2.57912298 0.91980226
[2,] -3.28994446 -2.57912298
[3,] 2.26015366 -3.28994446
[4,] 0.03781393 2.26015366
[5,] 4.44226886 0.03781393
[6,] 3.39780336 4.44226886
[7,] 1.26015366 3.39780336
[8,] -3.48355908 1.26015366
[9,] -0.33241665 -3.48355908
[10,] -3.25268906 -0.33241665
[11,] -0.66657454 -3.25268906
[12,] 1.27610744 -0.66657454
[13,] -0.29512213 1.27610744
[14,] 0.40667972 -0.29512213
[15,] 2.20364582 0.40667972
[16,] 4.48253160 2.20364582
[17,] -1.08884421 4.48253160
[18,] -0.49454570 -1.08884421
[19,] -5.20414151 -0.49454570
> z <- as.data.frame(dum1)
> print(z)
lag(myerror, k = 1) myerror
1 -2.57912298 0.91980226
2 -3.28994446 -2.57912298
3 2.26015366 -3.28994446
4 0.03781393 2.26015366
5 4.44226886 0.03781393
6 3.39780336 4.44226886
7 1.26015366 3.39780336
8 -3.48355908 1.26015366
9 -0.33241665 -3.48355908
10 -3.25268906 -0.33241665
11 -0.66657454 -3.25268906
12 1.27610744 -0.66657454
13 -0.29512213 1.27610744
14 0.40667972 -0.29512213
15 2.20364582 0.40667972
16 4.48253160 2.20364582
17 -1.08884421 4.48253160
18 -0.49454570 -1.08884421
19 -5.20414151 -0.49454570
> 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="/home/pw/wessanet/rcomp/tmp/77o5q1485164468.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="/home/pw/wessanet/rcomp/tmp/8741a1485164468.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="/home/pw/wessanet/rcomp/tmp/9atcs1485164468.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="/home/pw/wessanet/rcomp/tmp/106tr91485164468.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()
+ }
>
> 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="/home/pw/wessanet/rcomp/tmp/11oukl1485164468.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,'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="/home/pw/wessanet/rcomp/tmp/12mgw71485164468.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="/home/pw/wessanet/rcomp/tmp/13geah1485164468.tab")
> myr <- as.numeric(mysum$resid)
> myr
[1] 0.91980226 -2.57912298 -3.28994446 2.26015366 0.03781393 4.44226886
[7] 3.39780336 1.26015366 -3.48355908 -0.33241665 -3.25268906 -0.66657454
[13] 1.27610744 -0.29512213 0.40667972 2.20364582 4.48253160 -1.08884421
[19] -0.49454570 -5.20414151
> 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="/home/pw/wessanet/rcomp/tmp/144nzd1485164468.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="/home/pw/wessanet/rcomp/tmp/15tvt41485164468.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="/home/pw/wessanet/rcomp/tmp/16y3vr1485164468.tab")
+ }
+ }
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,'Ramsey RESET F-Test for powers (2 and 3) of fitted values',1,TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> reset_test_fitted <- resettest(mylm,power=2:3,type='fitted')
> a<-table.element(a,paste('
',RC.texteval('reset_test_fitted'),'',sep=''))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Ramsey RESET F-Test for powers (2 and 3) of regressors',1,TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> reset_test_regressors <- resettest(mylm,power=2:3,type='regressor')
Warning message:
In pf(reset, df1, df2, lower.tail = FALSE) : NaNs produced
> a<-table.element(a,paste('',RC.texteval('reset_test_regressors'),'',sep=''))
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Ramsey RESET F-Test for powers (2 and 3) of principal components',1,TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> reset_test_principal_components <- resettest(mylm,power=2:3,type='princomp')
> a<-table.element(a,paste('',RC.texteval('reset_test_principal_components'),'',sep=''))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/home/pw/wessanet/rcomp/tmp/17gtp21485164468.tab")
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,'Variance Inflation Factors (Multicollinearity)',1,TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> vif <- vif(mylm)
> a<-table.element(a,paste('',RC.texteval('vif'),'',sep=''))
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/home/pw/wessanet/rcomp/tmp/18d5h71485164468.tab")
>
> try(system("convert /home/pw/wessanet/rcomp/tmp/1rs3a1485164468.ps /home/pw/wessanet/rcomp/tmp/1rs3a1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/2wcai1485164468.ps /home/pw/wessanet/rcomp/tmp/2wcai1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/3vvrn1485164468.ps /home/pw/wessanet/rcomp/tmp/3vvrn1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/4obeq1485164468.ps /home/pw/wessanet/rcomp/tmp/4obeq1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/5sb541485164468.ps /home/pw/wessanet/rcomp/tmp/5sb541485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/6nzip1485164468.ps /home/pw/wessanet/rcomp/tmp/6nzip1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/77o5q1485164468.ps /home/pw/wessanet/rcomp/tmp/77o5q1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/8741a1485164468.ps /home/pw/wessanet/rcomp/tmp/8741a1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/9atcs1485164468.ps /home/pw/wessanet/rcomp/tmp/9atcs1485164468.png",intern=TRUE))
character(0)
> try(system("convert /home/pw/wessanet/rcomp/tmp/106tr91485164468.ps /home/pw/wessanet/rcomp/tmp/106tr91485164468.png",intern=TRUE))
convert: unable to open image `/home/pw/wessanet/rcomp/tmp/106tr91485164468.ps': No such file or directory @ error/blob.c/OpenBlob/2712.
convert: no images defined `/home/pw/wessanet/rcomp/tmp/106tr91485164468.png' @ error/convert.c/ConvertImageCommand/3210.
character(0)
attr(,"status")
[1] 1
Warning message:
running command 'convert /home/pw/wessanet/rcomp/tmp/106tr91485164468.ps /home/pw/wessanet/rcomp/tmp/106tr91485164468.png' had status 1
>
> proc.time()
user system elapsed
5.492 0.484 6.098