R version 2.15.2 (2012-10-26) -- "Trick or Treat" Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i686-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(1 + ,1 + ,1 + ,6 + ,6 + ,100 + ,100 + ,6 + ,6 + ,9 + ,9 + ,1 + ,2 + ,2 + ,9 + ,9 + ,99 + ,99 + ,2 + ,2 + ,8 + ,8 + ,1 + ,3 + ,3 + ,7 + ,7 + ,108 + ,108 + ,4 + ,4 + ,3 + ,3 + ,1 + ,4 + ,4 + ,8 + ,8 + ,103 + ,103 + ,0 + ,0 + ,4 + ,4 + ,1 + ,5 + ,5 + ,1 + ,1 + ,99 + ,99 + ,8 + ,8 + ,7 + ,7 + ,1 + ,6 + ,6 + ,9 + ,9 + ,115 + ,115 + ,0 + ,0 + ,7 + ,7 + ,1 + ,7 + ,7 + ,9 + ,9 + ,90 + ,90 + ,8 + ,8 + ,1 + ,1 + ,1 + ,8 + ,8 + ,7 + ,7 + ,95 + ,95 + ,9 + ,9 + ,9 + ,9 + ,1 + ,9 + ,9 + ,2 + ,2 + ,114 + ,114 + ,4 + ,4 + ,4 + ,4 + ,1 + ,10 + ,10 + ,9 + ,9 + ,108 + ,108 + ,2 + ,2 + ,9 + ,9 + ,1 + ,11 + ,11 + ,8 + ,8 + ,112 + ,112 + ,6 + ,6 + ,3 + ,3 + ,1 + ,12 + ,12 + ,3 + ,3 + ,109 + ,109 + ,1 + ,1 + ,3 + ,3 + ,1 + ,13 + ,13 + ,0 + ,0 + ,105 + ,105 + ,0 + ,0 + ,3 + ,3 + ,1 + ,14 + ,14 + ,7 + ,7 + ,105 + ,105 + ,0 + ,0 + ,2 + ,2 + ,1 + ,15 + ,15 + ,5 + ,5 + ,118 + ,118 + ,5 + ,5 + ,8 + ,8 + ,1 + ,16 + ,16 + ,7 + ,7 + ,103 + ,103 + ,7 + ,7 + ,6 + ,6 + ,1 + ,17 + ,17 + ,9 + ,9 + ,112 + ,112 + ,5 + ,5 + ,2 + ,2 + ,1 + ,18 + ,18 + ,6 + ,6 + ,116 + ,116 + ,6 + ,6 + ,6 + ,6 + ,1 + ,19 + ,19 + ,4 + ,4 + ,96 + ,96 + ,6 + ,6 + ,6 + ,6 + ,1 + ,20 + ,20 + ,5 + ,5 + ,101 + ,101 + ,9 + ,9 + ,0 + ,0 + ,1 + ,21 + ,21 + ,8 + ,8 + ,116 + ,116 + ,5 + ,5 + ,4 + ,4 + ,1 + ,22 + ,22 + ,5 + ,5 + ,119 + ,119 + ,3 + ,3 + ,9 + ,9 + ,1 + ,23 + ,23 + ,9 + ,9 + ,115 + ,115 + ,4 + ,4 + ,5 + ,5 + ,1 + ,24 + ,24 + ,0 + ,0 + ,108 + ,108 + ,5 + ,5 + ,2 + ,2 + ,1 + ,25 + ,25 + ,0 + ,0 + ,111 + ,111 + ,5 + ,5 + ,8 + ,8 + ,1 + ,26 + ,26 + ,3 + ,3 + ,108 + ,108 + ,8 + ,8 + ,3 + ,3 + ,1 + ,27 + ,27 + ,8 + ,8 + ,121 + ,121 + ,8 + ,8 + ,9 + ,9 + ,1 + ,28 + ,28 + ,1 + ,1 + ,109 + ,109 + ,6 + ,6 + ,8 + ,8 + ,1 + ,29 + ,29 + ,3 + ,3 + ,112 + ,112 + ,2 + ,2 + ,8 + ,8 + ,1 + ,30 + ,30 + ,2 + ,2 + ,119 + ,119 + ,6 + ,6 + ,8 + ,8 + ,1 + ,31 + ,31 + ,5 + ,5 + ,104 + ,104 + ,1 + ,1 + ,5 + ,5 + ,1 + ,32 + ,32 + ,2 + ,2 + ,105 + ,105 + ,3 + ,3 + ,4 + ,4 + ,1 + ,33 + ,33 + ,5 + ,5 + ,115 + ,115 + ,0 + ,0 + ,4 + ,4 + ,1 + ,34 + ,34 + ,4 + ,4 + ,124 + ,124 + ,1 + ,1 + ,1 + ,1 + ,1 + ,35 + ,35 + ,3 + ,3 + ,116 + ,116 + ,8 + ,8 + ,6 + ,6 + ,1 + ,36 + ,36 + ,0 + ,0 + ,107 + ,107 + ,5 + ,5 + ,2 + ,2 + ,1 + ,37 + ,37 + ,7 + ,7 + ,115 + ,115 + ,6 + ,6 + ,1 + ,1 + ,1 + ,38 + ,38 + ,8 + ,8 + ,116 + ,116 + ,2 + ,2 + ,3 + ,3 + ,1 + ,39 + ,39 + ,8 + ,8 + ,116 + ,116 + ,3 + ,3 + ,8 + ,8 + ,1 + ,40 + ,40 + ,3 + ,3 + ,119 + ,119 + ,0 + ,0 + ,9 + ,9 + ,1 + ,41 + ,41 + ,1 + ,1 + ,111 + ,111 + ,9 + ,9 + ,1 + ,1 + ,1 + ,42 + ,42 + ,9 + ,9 + ,118 + ,118 + ,6 + ,6 + ,7 + ,7 + ,1 + ,43 + ,43 + ,0 + ,0 + ,106 + ,106 + ,9 + ,9 + ,2 + ,2 + ,1 + ,44 + ,44 + ,8 + ,8 + ,103 + ,103 + ,2 + ,2 + ,5 + ,5 + ,1 + ,45 + ,45 + ,8 + ,8 + ,118 + ,118 + ,6 + ,6 + ,0 + ,0 + ,1 + ,46 + ,46 + ,7 + ,7 + ,118 + ,118 + ,7 + ,7 + ,5 + ,5 + ,1 + ,47 + ,47 + ,4 + ,4 + ,102 + ,102 + ,8 + ,8 + ,0 + ,0 + ,1 + ,48 + ,48 + ,3 + ,3 + ,100 + ,100 + ,6 + ,6 + ,1 + ,1 + ,1 + ,49 + ,49 + ,0 + ,0 + ,94 + ,94 + ,9 + ,9 + ,6 + ,6 + ,1 + ,50 + ,50 + ,2 + ,2 + ,94 + ,94 + ,5 + ,5 + ,3 + ,3 + ,1 + ,51 + ,51 + ,1 + ,1 + ,102 + ,102 + ,9 + ,9 + ,9 + ,9 + ,1 + ,52 + ,52 + ,1 + ,1 + ,95 + ,95 + ,3 + ,3 + ,3 + ,3 + ,1 + ,53 + ,53 + ,8 + ,8 + ,92 + ,92 + ,5 + ,5 + ,5 + ,5 + ,1 + ,54 + ,54 + ,7 + ,7 + ,102 + ,102 + ,7 + ,7 + ,8 + ,8 + ,1 + ,55 + ,55 + ,6 + ,6 + ,91 + ,91 + ,5 + ,5 + ,7 + ,7 + ,1 + ,56 + ,56 + ,1 + ,1 + ,89 + ,89 + ,5 + ,5 + ,4 + ,4 + ,1 + ,57 + ,57 + ,5 + ,5 + ,104 + ,104 + ,2 + ,2 + ,8 + ,8 + ,1 + ,58 + ,58 + ,1 + ,1 + ,105 + ,105 + ,2 + ,2 + ,1 + ,1 + ,1 + ,59 + ,59 + ,1 + ,1 + ,99 + ,99 + ,0 + ,0 + ,2 + ,2 + ,1 + ,60 + ,60 + ,7 + ,7 + ,95 + ,95 + ,5 + ,5 + ,0 + ,0 + ,1 + ,61 + ,61 + ,3 + ,3 + ,90 + ,90 + ,5 + ,5 + ,8 + ,8 + ,1 + ,62 + ,62 + ,8 + ,8 + ,96 + ,96 + ,1 + ,1 + ,7 + ,7 + ,1 + ,63 + ,63 + ,5 + ,5 + ,113 + ,113 + ,0 + ,0 + ,5 + ,5 + ,1 + ,64 + ,64 + ,7 + ,7 + ,101 + ,101 + ,9 + ,9 + ,0 + ,0 + ,1 + ,65 + ,65 + ,5 + ,5 + ,101 + ,101 + ,4 + ,4 + ,9 + ,9 + ,1 + ,66 + ,66 + ,7 + ,7 + ,113 + ,113 + ,6 + ,6 + ,8 + ,8 + ,1 + ,67 + ,67 + ,2 + ,2 + ,96 + ,96 + ,6 + ,6 + ,2 + ,2 + ,1 + ,68 + ,68 + ,4 + ,4 + ,97 + ,97 + ,8 + ,8 + ,2 + ,2 + ,1 + ,69 + ,69 + ,0 + ,0 + ,114 + ,114 + ,9 + ,9 + ,9 + ,9 + ,1 + ,70 + ,70 + ,0 + ,0 + ,112 + ,112 + ,5 + ,5 + ,5 + ,5 + ,1 + ,71 + ,71 + ,5 + ,5 + ,108 + ,108 + ,4 + ,4 + ,9 + ,9 + ,1 + ,72 + ,72 + ,3 + ,3 + ,107 + ,107 + ,0 + ,0 + ,0 + ,0 + ,1 + ,73 + ,73 + ,1 + ,1 + ,103 + ,103 + ,5 + ,5 + ,9 + ,9 + ,1 + ,74 + ,74 + ,1 + ,1 + ,107 + ,107 + ,5 + ,5 + ,0 + ,0 + ,1 + ,75 + ,75 + ,3 + ,3 + ,122 + ,122 + ,3 + ,3 + ,9 + ,9) + ,dim=c(11 + ,75) + ,dimnames=list(c('pop' + ,'md' + ,'md-p' + ,'steenkool' + ,'steenkool_p' + ,'aardolie' + ,'aardolie_p' + ,'uranium' + ,'uranium_p' + ,'metaal' + ,'metaal_p ') + ,1:75)) > y <- array(NA,dim=c(11,75),dimnames=list(c('pop','md','md-p','steenkool','steenkool_p','aardolie','aardolie_p','uranium','uranium_p','metaal','metaal_p '),1:75)) > 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 = '8' > library(lattice) > library(lmtest) Loading required package: zoo Attaching package: 'zoo' The following object(s) are masked from 'package:base': as.Date, 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 uranium pop md md-p steenkool steenkool_p aardolie aardolie_p uranium_p 1 6 1 1 1 6 6 100 100 6 2 2 1 2 2 9 9 99 99 2 3 4 1 3 3 7 7 108 108 4 4 0 1 4 4 8 8 103 103 0 5 8 1 5 5 1 1 99 99 8 6 0 1 6 6 9 9 115 115 0 7 8 1 7 7 9 9 90 90 8 8 9 1 8 8 7 7 95 95 9 9 4 1 9 9 2 2 114 114 4 10 2 1 10 10 9 9 108 108 2 11 6 1 11 11 8 8 112 112 6 12 1 1 12 12 3 3 109 109 1 13 0 1 13 13 0 0 105 105 0 14 0 1 14 14 7 7 105 105 0 15 5 1 15 15 5 5 118 118 5 16 7 1 16 16 7 7 103 103 7 17 5 1 17 17 9 9 112 112 5 18 6 1 18 18 6 6 116 116 6 19 6 1 19 19 4 4 96 96 6 20 9 1 20 20 5 5 101 101 9 21 5 1 21 21 8 8 116 116 5 22 3 1 22 22 5 5 119 119 3 23 4 1 23 23 9 9 115 115 4 24 5 1 24 24 0 0 108 108 5 25 5 1 25 25 0 0 111 111 5 26 8 1 26 26 3 3 108 108 8 27 8 1 27 27 8 8 121 121 8 28 6 1 28 28 1 1 109 109 6 29 2 1 29 29 3 3 112 112 2 30 6 1 30 30 2 2 119 119 6 31 1 1 31 31 5 5 104 104 1 32 3 1 32 32 2 2 105 105 3 33 0 1 33 33 5 5 115 115 0 34 1 1 34 34 4 4 124 124 1 35 8 1 35 35 3 3 116 116 8 36 5 1 36 36 0 0 107 107 5 37 6 1 37 37 7 7 115 115 6 38 2 1 38 38 8 8 116 116 2 39 3 1 39 39 8 8 116 116 3 40 0 1 40 40 3 3 119 119 0 41 9 1 41 41 1 1 111 111 9 42 6 1 42 42 9 9 118 118 6 43 9 1 43 43 0 0 106 106 9 44 2 1 44 44 8 8 103 103 2 45 6 1 45 45 8 8 118 118 6 46 7 1 46 46 7 7 118 118 7 47 8 1 47 47 4 4 102 102 8 48 6 1 48 48 3 3 100 100 6 49 9 1 49 49 0 0 94 94 9 50 5 1 50 50 2 2 94 94 5 51 9 1 51 51 1 1 102 102 9 52 3 1 52 52 1 1 95 95 3 53 5 1 53 53 8 8 92 92 5 54 7 1 54 54 7 7 102 102 7 55 5 1 55 55 6 6 91 91 5 56 5 1 56 56 1 1 89 89 5 57 2 1 57 57 5 5 104 104 2 58 2 1 58 58 1 1 105 105 2 59 0 1 59 59 1 1 99 99 0 60 5 1 60 60 7 7 95 95 5 61 5 1 61 61 3 3 90 90 5 62 1 1 62 62 8 8 96 96 1 63 0 1 63 63 5 5 113 113 0 64 9 1 64 64 7 7 101 101 9 65 4 1 65 65 5 5 101 101 4 66 6 1 66 66 7 7 113 113 6 67 6 1 67 67 2 2 96 96 6 68 8 1 68 68 4 4 97 97 8 69 9 1 69 69 0 0 114 114 9 70 5 1 70 70 0 0 112 112 5 71 4 1 71 71 5 5 108 108 4 72 0 1 72 72 3 3 107 107 0 73 5 1 73 73 1 1 103 103 5 74 5 1 74 74 1 1 107 107 5 75 3 1 75 75 3 3 122 122 3 metaal metaal_p\r 1 9 9 2 8 8 3 3 3 4 4 4 5 7 7 6 7 7 7 1 1 8 9 9 9 4 4 10 9 9 11 3 3 12 3 3 13 3 3 14 2 2 15 8 8 16 6 6 17 2 2 18 6 6 19 6 6 20 0 0 21 4 4 22 9 9 23 5 5 24 2 2 25 8 8 26 3 3 27 9 9 28 8 8 29 8 8 30 8 8 31 5 5 32 4 4 33 4 4 34 1 1 35 6 6 36 2 2 37 1 1 38 3 3 39 8 8 40 9 9 41 1 1 42 7 7 43 2 2 44 5 5 45 0 0 46 5 5 47 0 0 48 1 1 49 6 6 50 3 3 51 9 9 52 3 3 53 5 5 54 8 8 55 7 7 56 4 4 57 8 8 58 1 1 59 2 2 60 0 0 61 8 8 62 7 7 63 5 5 64 0 0 65 9 9 66 8 8 67 2 2 68 2 2 69 9 9 70 5 5 71 9 9 72 0 0 73 9 9 74 0 0 75 9 9 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) pop md `md-p` steenkool 0.000e+00 NA 0.000e+00 NA 0.000e+00 steenkool_p aardolie aardolie_p uranium_p metaal NA 0.000e+00 NA 1.000e+00 -8.995e-18 `metaal_p\\r` NA > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -1.533e-15 -6.263e-17 -3.770e-18 1.035e-16 5.277e-16 Coefficients: (5 not defined because of singularities) Estimate Std. Error t value Pr(>|t|) (Intercept) 0.000e+00 3.705e-16 0.000e+00 1.000 pop NA NA NA NA md 0.000e+00 1.361e-18 0.000e+00 1.000 `md-p` NA NA NA NA steenkool 0.000e+00 9.857e-18 0.000e+00 1.000 steenkool_p NA NA NA NA aardolie 0.000e+00 3.278e-18 0.000e+00 1.000 aardolie_p NA NA NA NA uranium_p 1.000e+00 1.016e-17 9.839e+16 <2e-16 *** metaal -8.995e-18 9.232e-18 -9.740e-01 0.333 `metaal_p\\r` NA NA NA NA --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 2.405e-16 on 69 degrees of freedom Multiple R-squared: 1, Adjusted R-squared: 1 F-statistic: 2.016e+33 on 5 and 69 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,] 9.973617e-01 5.276664e-03 2.638332e-03 [2,] 9.842718e-01 3.145635e-02 1.572817e-02 [3,] 1.946816e-02 3.893633e-02 9.805318e-01 [4,] 1.000000e+00 1.207738e-09 6.038692e-10 [5,] 6.875886e-01 6.248228e-01 3.124114e-01 [6,] 9.023573e-01 1.952855e-01 9.764275e-02 [7,] 3.750192e-01 7.500384e-01 6.249808e-01 [8,] 9.549532e-01 9.009359e-02 4.504680e-02 [9,] 9.994662e-01 1.067503e-03 5.337515e-04 [10,] 9.915678e-01 1.686441e-02 8.432204e-03 [11,] 3.947762e-04 7.895523e-04 9.996052e-01 [12,] 4.225601e-03 8.451202e-03 9.957744e-01 [13,] 9.999154e-01 1.691419e-04 8.457094e-05 [14,] 5.064556e-02 1.012911e-01 9.493544e-01 [15,] 2.287459e-03 4.574917e-03 9.977125e-01 [16,] 1.059822e-02 2.119644e-02 9.894018e-01 [17,] 1.394988e-03 2.789976e-03 9.986050e-01 [18,] 9.523919e-01 9.521629e-02 4.760815e-02 [19,] 2.715861e-03 5.431723e-03 9.972841e-01 [20,] 9.999979e-01 4.185260e-06 2.092630e-06 [21,] 5.045533e-01 9.908933e-01 4.954467e-01 [22,] 4.182534e-04 8.365069e-04 9.995817e-01 [23,] 9.999966e-01 6.846839e-06 3.423420e-06 [24,] 6.117213e-03 1.223443e-02 9.938828e-01 [25,] 9.993543e-01 1.291306e-03 6.456531e-04 [26,] 9.910316e-01 1.793683e-02 8.968416e-03 [27,] 1.000000e+00 5.581753e-09 2.790876e-09 [28,] 1.000000e+00 1.232098e-08 6.160490e-09 [29,] 1.523376e-04 3.046751e-04 9.998477e-01 [30,] 9.834464e-01 3.310723e-02 1.655362e-02 [31,] 7.536148e-08 1.507230e-07 9.999999e-01 [32,] 7.249475e-10 1.449895e-09 1.000000e+00 [33,] 5.210347e-01 9.579307e-01 4.789653e-01 [34,] 1.150104e-03 2.300208e-03 9.988499e-01 [35,] 9.341833e-02 1.868367e-01 9.065817e-01 [36,] 6.011833e-01 7.976335e-01 3.988167e-01 [37,] 1.736750e-01 3.473501e-01 8.263250e-01 [38,] 7.209588e-02 1.441918e-01 9.279041e-01 [39,] 9.920166e-01 1.596672e-02 7.983362e-03 [40,] 5.596009e-01 8.807983e-01 4.403991e-01 [41,] 9.999619e-01 7.618270e-05 3.809135e-05 [42,] 1.341788e-12 2.683575e-12 1.000000e+00 [43,] 9.863405e-01 2.731906e-02 1.365953e-02 [44,] 7.975860e-01 4.048280e-01 2.024140e-01 [45,] 8.439590e-04 1.687918e-03 9.991560e-01 [46,] 2.141796e-01 4.283591e-01 7.858204e-01 [47,] 8.083545e-07 1.616709e-06 9.999992e-01 [48,] 9.322236e-01 1.355527e-01 6.777635e-02 > postscript(file="/var/wessaorg/rcomp/tmp/1o2iv1353427664.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/2nbgj1353427664.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/37elf1353427664.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/467me1353427664.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/5deib1353427664.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 = 75 Frequency = 1 1 2 3 4 5 2.974999e-16 5.276569e-16 -3.362190e-16 1.512909e-16 -1.533063e-15 6 7 8 9 10 2.740009e-16 -9.746531e-18 -8.000598e-17 1.068820e-16 1.272721e-16 11 12 13 14 15 -6.035392e-17 -1.747588e-16 1.609248e-16 -1.066397e-16 7.356488e-17 16 17 18 19 20 -3.766125e-18 6.422507e-17 -2.478598e-17 8.889301e-17 -9.051145e-17 21 22 23 24 25 2.071326e-18 -2.540163e-17 -3.996012e-17 1.001829e-16 1.406986e-16 26 27 28 29 30 2.561890e-16 1.239382e-17 1.581200e-16 -5.892811e-17 5.704764e-17 31 32 33 34 35 -4.780707e-17 -2.820709e-17 -2.342015e-16 1.525148e-16 -3.001291e-17 36 37 38 39 40 8.159760e-17 -4.049543e-17 -4.172801e-17 -1.739927e-16 -6.484991e-17 41 42 43 44 45 2.746864e-16 -1.167880e-16 3.983477e-16 -1.258416e-16 -8.636746e-17 46 47 48 49 50 -2.990220e-17 1.936219e-16 -6.428363e-17 -6.098184e-17 3.654287e-17 51 52 53 54 55 2.118099e-16 1.170618e-16 -5.208124e-17 1.949972e-17 -2.082902e-17 56 57 58 59 60 2.841749e-18 2.179716e-16 9.726184e-17 8.796828e-17 -1.006370e-16 61 62 63 64 65 2.599202e-17 8.812870e-17 -2.605869e-16 1.104741e-16 1.997104e-17 66 67 68 69 70 1.254762e-17 -3.745572e-17 -4.720149e-17 1.115423e-16 -7.194523e-17 71 72 73 74 75 -1.519386e-16 -2.599741e-16 6.551183e-18 -5.489818e-17 -1.186992e-16 > postscript(file="/var/wessaorg/rcomp/tmp/6lm0t1353427664.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 = 75 Frequency = 1 lag(myerror, k = 1) myerror 0 2.974999e-16 NA 1 5.276569e-16 2.974999e-16 2 -3.362190e-16 5.276569e-16 3 1.512909e-16 -3.362190e-16 4 -1.533063e-15 1.512909e-16 5 2.740009e-16 -1.533063e-15 6 -9.746531e-18 2.740009e-16 7 -8.000598e-17 -9.746531e-18 8 1.068820e-16 -8.000598e-17 9 1.272721e-16 1.068820e-16 10 -6.035392e-17 1.272721e-16 11 -1.747588e-16 -6.035392e-17 12 1.609248e-16 -1.747588e-16 13 -1.066397e-16 1.609248e-16 14 7.356488e-17 -1.066397e-16 15 -3.766125e-18 7.356488e-17 16 6.422507e-17 -3.766125e-18 17 -2.478598e-17 6.422507e-17 18 8.889301e-17 -2.478598e-17 19 -9.051145e-17 8.889301e-17 20 2.071326e-18 -9.051145e-17 21 -2.540163e-17 2.071326e-18 22 -3.996012e-17 -2.540163e-17 23 1.001829e-16 -3.996012e-17 24 1.406986e-16 1.001829e-16 25 2.561890e-16 1.406986e-16 26 1.239382e-17 2.561890e-16 27 1.581200e-16 1.239382e-17 28 -5.892811e-17 1.581200e-16 29 5.704764e-17 -5.892811e-17 30 -4.780707e-17 5.704764e-17 31 -2.820709e-17 -4.780707e-17 32 -2.342015e-16 -2.820709e-17 33 1.525148e-16 -2.342015e-16 34 -3.001291e-17 1.525148e-16 35 8.159760e-17 -3.001291e-17 36 -4.049543e-17 8.159760e-17 37 -4.172801e-17 -4.049543e-17 38 -1.739927e-16 -4.172801e-17 39 -6.484991e-17 -1.739927e-16 40 2.746864e-16 -6.484991e-17 41 -1.167880e-16 2.746864e-16 42 3.983477e-16 -1.167880e-16 43 -1.258416e-16 3.983477e-16 44 -8.636746e-17 -1.258416e-16 45 -2.990220e-17 -8.636746e-17 46 1.936219e-16 -2.990220e-17 47 -6.428363e-17 1.936219e-16 48 -6.098184e-17 -6.428363e-17 49 3.654287e-17 -6.098184e-17 50 2.118099e-16 3.654287e-17 51 1.170618e-16 2.118099e-16 52 -5.208124e-17 1.170618e-16 53 1.949972e-17 -5.208124e-17 54 -2.082902e-17 1.949972e-17 55 2.841749e-18 -2.082902e-17 56 2.179716e-16 2.841749e-18 57 9.726184e-17 2.179716e-16 58 8.796828e-17 9.726184e-17 59 -1.006370e-16 8.796828e-17 60 2.599202e-17 -1.006370e-16 61 8.812870e-17 2.599202e-17 62 -2.605869e-16 8.812870e-17 63 1.104741e-16 -2.605869e-16 64 1.997104e-17 1.104741e-16 65 1.254762e-17 1.997104e-17 66 -3.745572e-17 1.254762e-17 67 -4.720149e-17 -3.745572e-17 68 1.115423e-16 -4.720149e-17 69 -7.194523e-17 1.115423e-16 70 -1.519386e-16 -7.194523e-17 71 -2.599741e-16 -1.519386e-16 72 6.551183e-18 -2.599741e-16 73 -5.489818e-17 6.551183e-18 74 -1.186992e-16 -5.489818e-17 75 NA -1.186992e-16 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] 5.276569e-16 2.974999e-16 [2,] -3.362190e-16 5.276569e-16 [3,] 1.512909e-16 -3.362190e-16 [4,] -1.533063e-15 1.512909e-16 [5,] 2.740009e-16 -1.533063e-15 [6,] -9.746531e-18 2.740009e-16 [7,] -8.000598e-17 -9.746531e-18 [8,] 1.068820e-16 -8.000598e-17 [9,] 1.272721e-16 1.068820e-16 [10,] -6.035392e-17 1.272721e-16 [11,] -1.747588e-16 -6.035392e-17 [12,] 1.609248e-16 -1.747588e-16 [13,] -1.066397e-16 1.609248e-16 [14,] 7.356488e-17 -1.066397e-16 [15,] -3.766125e-18 7.356488e-17 [16,] 6.422507e-17 -3.766125e-18 [17,] -2.478598e-17 6.422507e-17 [18,] 8.889301e-17 -2.478598e-17 [19,] -9.051145e-17 8.889301e-17 [20,] 2.071326e-18 -9.051145e-17 [21,] -2.540163e-17 2.071326e-18 [22,] -3.996012e-17 -2.540163e-17 [23,] 1.001829e-16 -3.996012e-17 [24,] 1.406986e-16 1.001829e-16 [25,] 2.561890e-16 1.406986e-16 [26,] 1.239382e-17 2.561890e-16 [27,] 1.581200e-16 1.239382e-17 [28,] -5.892811e-17 1.581200e-16 [29,] 5.704764e-17 -5.892811e-17 [30,] -4.780707e-17 5.704764e-17 [31,] -2.820709e-17 -4.780707e-17 [32,] -2.342015e-16 -2.820709e-17 [33,] 1.525148e-16 -2.342015e-16 [34,] -3.001291e-17 1.525148e-16 [35,] 8.159760e-17 -3.001291e-17 [36,] -4.049543e-17 8.159760e-17 [37,] -4.172801e-17 -4.049543e-17 [38,] -1.739927e-16 -4.172801e-17 [39,] -6.484991e-17 -1.739927e-16 [40,] 2.746864e-16 -6.484991e-17 [41,] -1.167880e-16 2.746864e-16 [42,] 3.983477e-16 -1.167880e-16 [43,] -1.258416e-16 3.983477e-16 [44,] -8.636746e-17 -1.258416e-16 [45,] -2.990220e-17 -8.636746e-17 [46,] 1.936219e-16 -2.990220e-17 [47,] -6.428363e-17 1.936219e-16 [48,] -6.098184e-17 -6.428363e-17 [49,] 3.654287e-17 -6.098184e-17 [50,] 2.118099e-16 3.654287e-17 [51,] 1.170618e-16 2.118099e-16 [52,] -5.208124e-17 1.170618e-16 [53,] 1.949972e-17 -5.208124e-17 [54,] -2.082902e-17 1.949972e-17 [55,] 2.841749e-18 -2.082902e-17 [56,] 2.179716e-16 2.841749e-18 [57,] 9.726184e-17 2.179716e-16 [58,] 8.796828e-17 9.726184e-17 [59,] -1.006370e-16 8.796828e-17 [60,] 2.599202e-17 -1.006370e-16 [61,] 8.812870e-17 2.599202e-17 [62,] -2.605869e-16 8.812870e-17 [63,] 1.104741e-16 -2.605869e-16 [64,] 1.997104e-17 1.104741e-16 [65,] 1.254762e-17 1.997104e-17 [66,] -3.745572e-17 1.254762e-17 [67,] -4.720149e-17 -3.745572e-17 [68,] 1.115423e-16 -4.720149e-17 [69,] -7.194523e-17 1.115423e-16 [70,] -1.519386e-16 -7.194523e-17 [71,] -2.599741e-16 -1.519386e-16 [72,] 6.551183e-18 -2.599741e-16 [73,] -5.489818e-17 6.551183e-18 [74,] -1.186992e-16 -5.489818e-17 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 5.276569e-16 2.974999e-16 2 -3.362190e-16 5.276569e-16 3 1.512909e-16 -3.362190e-16 4 -1.533063e-15 1.512909e-16 5 2.740009e-16 -1.533063e-15 6 -9.746531e-18 2.740009e-16 7 -8.000598e-17 -9.746531e-18 8 1.068820e-16 -8.000598e-17 9 1.272721e-16 1.068820e-16 10 -6.035392e-17 1.272721e-16 11 -1.747588e-16 -6.035392e-17 12 1.609248e-16 -1.747588e-16 13 -1.066397e-16 1.609248e-16 14 7.356488e-17 -1.066397e-16 15 -3.766125e-18 7.356488e-17 16 6.422507e-17 -3.766125e-18 17 -2.478598e-17 6.422507e-17 18 8.889301e-17 -2.478598e-17 19 -9.051145e-17 8.889301e-17 20 2.071326e-18 -9.051145e-17 21 -2.540163e-17 2.071326e-18 22 -3.996012e-17 -2.540163e-17 23 1.001829e-16 -3.996012e-17 24 1.406986e-16 1.001829e-16 25 2.561890e-16 1.406986e-16 26 1.239382e-17 2.561890e-16 27 1.581200e-16 1.239382e-17 28 -5.892811e-17 1.581200e-16 29 5.704764e-17 -5.892811e-17 30 -4.780707e-17 5.704764e-17 31 -2.820709e-17 -4.780707e-17 32 -2.342015e-16 -2.820709e-17 33 1.525148e-16 -2.342015e-16 34 -3.001291e-17 1.525148e-16 35 8.159760e-17 -3.001291e-17 36 -4.049543e-17 8.159760e-17 37 -4.172801e-17 -4.049543e-17 38 -1.739927e-16 -4.172801e-17 39 -6.484991e-17 -1.739927e-16 40 2.746864e-16 -6.484991e-17 41 -1.167880e-16 2.746864e-16 42 3.983477e-16 -1.167880e-16 43 -1.258416e-16 3.983477e-16 44 -8.636746e-17 -1.258416e-16 45 -2.990220e-17 -8.636746e-17 46 1.936219e-16 -2.990220e-17 47 -6.428363e-17 1.936219e-16 48 -6.098184e-17 -6.428363e-17 49 3.654287e-17 -6.098184e-17 50 2.118099e-16 3.654287e-17 51 1.170618e-16 2.118099e-16 52 -5.208124e-17 1.170618e-16 53 1.949972e-17 -5.208124e-17 54 -2.082902e-17 1.949972e-17 55 2.841749e-18 -2.082902e-17 56 2.179716e-16 2.841749e-18 57 9.726184e-17 2.179716e-16 58 8.796828e-17 9.726184e-17 59 -1.006370e-16 8.796828e-17 60 2.599202e-17 -1.006370e-16 61 8.812870e-17 2.599202e-17 62 -2.605869e-16 8.812870e-17 63 1.104741e-16 -2.605869e-16 64 1.997104e-17 1.104741e-16 65 1.254762e-17 1.997104e-17 66 -3.745572e-17 1.254762e-17 67 -4.720149e-17 -3.745572e-17 68 1.115423e-16 -4.720149e-17 69 -7.194523e-17 1.115423e-16 70 -1.519386e-16 -7.194523e-17 71 -2.599741e-16 -1.519386e-16 72 6.551183e-18 -2.599741e-16 73 -5.489818e-17 6.551183e-18 74 -1.186992e-16 -5.489818e-17 > 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/7im8u1353427664.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/8epbt1353427664.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/9ks3z1353427664.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/100ppg1353427664.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, 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='') + } + } Error: subscript out of bounds Execution halted