R version 2.9.0 (2009-04-17) Copyright (C) 2009 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(2649.2 + ,31077 + ,0 + ,2579.4 + ,31293 + ,0 + ,2504.6 + ,30236 + ,0 + ,2462.3 + ,30160 + ,0 + ,2467.4 + ,32436 + ,0 + ,2446.7 + ,30695 + ,0 + ,2656.3 + ,27525 + ,0 + ,2626.2 + ,26434 + ,0 + ,2482.6 + ,25739 + ,0 + ,2539.9 + ,25204 + ,0 + ,2502.7 + ,24977 + ,0 + ,2466.9 + ,24320 + ,0 + ,2513.2 + ,22680 + ,1 + ,2443.3 + ,22052 + ,1 + ,2293.4 + ,21467 + ,1 + ,2070.8 + ,21383 + ,1 + ,2029.6 + ,21777 + ,1 + ,2052 + ,21928 + ,1 + ,1864.4 + ,21814 + ,1 + ,1670.1 + ,22937 + ,1 + ,1811 + ,23595 + ,1 + ,1905.4 + ,20830 + ,1 + ,1862.8 + ,19650 + ,1 + ,2014.5 + ,19195 + ,1 + ,2197.8 + ,19644 + ,1 + ,2962.3 + ,18483 + ,0 + ,3047 + ,18079 + ,0 + ,3032.6 + ,19178 + ,0 + ,3504.4 + ,18391 + ,0 + ,3801.1 + ,18441 + ,0 + ,3857.6 + ,18584 + ,0 + ,3674.4 + ,20108 + ,0 + ,3721 + ,20148 + ,0 + ,3844.5 + ,19394 + ,0 + ,4116.7 + ,17745 + ,0 + ,4105.2 + ,17696 + ,0 + ,4435.2 + ,17032 + ,0 + ,4296.5 + ,16438 + ,0 + ,4202.5 + ,15683 + ,0 + ,4562.8 + ,15594 + ,0 + ,4621.4 + ,15713 + ,0 + ,4697 + ,15937 + ,0 + ,4591.3 + ,16171 + ,0 + ,4357 + ,15928 + ,0 + ,4502.6 + ,16348 + ,0 + ,4443.9 + ,15579 + ,0 + ,4290.9 + ,15305 + ,0 + ,4199.8 + ,15648 + ,0 + ,4138.5 + ,14954 + ,0 + ,3970.1 + ,15137 + ,0 + ,3862.3 + ,15839 + ,0 + ,3701.6 + ,16050 + ,0 + ,3570.12 + ,15168 + ,0 + ,3801.06 + ,17064 + ,0 + ,3895.51 + ,16005 + ,0 + ,3917.96 + ,14886 + ,0 + ,3813.06 + ,14931 + ,0 + ,3667.03 + ,14544 + ,0 + ,3494.17 + ,13812 + ,0 + ,3364 + ,13031 + ,0 + ,3295.3 + ,12574 + ,0 + ,3277.0 + ,11964 + ,0 + ,3257.2 + ,11451 + ,0 + ,3161.7 + ,11346 + ,0 + ,3097.3 + ,11353 + ,0 + ,3061.3 + ,10702 + ,0 + ,3119.3 + ,10646 + ,0 + ,3106.22 + ,10556 + ,0 + ,3080.58 + ,10463 + ,0 + ,2981.85 + ,10407 + ,0) + ,dim=c(3 + ,70) + ,dimnames=list(c('Bel20' + ,'Goudprijs' + ,'Crisis') + ,1:70)) > y <- array(NA,dim=c(3,70),dimnames=list(c('Bel20','Goudprijs','Crisis'),1:70)) > 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' > 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 Bel20 Goudprijs Crisis 1 2649.20 31077 0 2 2579.40 31293 0 3 2504.60 30236 0 4 2462.30 30160 0 5 2467.40 32436 0 6 2446.70 30695 0 7 2656.30 27525 0 8 2626.20 26434 0 9 2482.60 25739 0 10 2539.90 25204 0 11 2502.70 24977 0 12 2466.90 24320 0 13 2513.20 22680 1 14 2443.30 22052 1 15 2293.40 21467 1 16 2070.80 21383 1 17 2029.60 21777 1 18 2052.00 21928 1 19 1864.40 21814 1 20 1670.10 22937 1 21 1811.00 23595 1 22 1905.40 20830 1 23 1862.80 19650 1 24 2014.50 19195 1 25 2197.80 19644 1 26 2962.30 18483 0 27 3047.00 18079 0 28 3032.60 19178 0 29 3504.40 18391 0 30 3801.10 18441 0 31 3857.60 18584 0 32 3674.40 20108 0 33 3721.00 20148 0 34 3844.50 19394 0 35 4116.70 17745 0 36 4105.20 17696 0 37 4435.20 17032 0 38 4296.50 16438 0 39 4202.50 15683 0 40 4562.80 15594 0 41 4621.40 15713 0 42 4697.00 15937 0 43 4591.30 16171 0 44 4357.00 15928 0 45 4502.60 16348 0 46 4443.90 15579 0 47 4290.90 15305 0 48 4199.80 15648 0 49 4138.50 14954 0 50 3970.10 15137 0 51 3862.30 15839 0 52 3701.60 16050 0 53 3570.12 15168 0 54 3801.06 17064 0 55 3895.51 16005 0 56 3917.96 14886 0 57 3813.06 14931 0 58 3667.03 14544 0 59 3494.17 13812 0 60 3364.00 13031 0 61 3295.30 12574 0 62 3277.00 11964 0 63 3257.20 11451 0 64 3161.70 11346 0 65 3097.30 11353 0 66 3061.30 10702 0 67 3119.30 10646 0 68 3106.22 10556 0 69 3080.58 10463 0 70 2981.85 10407 0 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) Goudprijs Crisis 4613.2218 -0.0612 -1244.0672 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -994.51 -397.60 -31.53 408.99 1059.06 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 4.613e+03 2.283e+02 20.205 < 2e-16 *** Goudprijs -6.120e-02 1.199e-02 -5.106 2.93e-06 *** Crisis -1.244e+03 1.717e+02 -7.247 5.47e-10 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 542.8 on 67 degrees of freedom Multiple R-squared: 0.6028, Adjusted R-squared: 0.591 F-statistic: 50.84 on 2 and 67 DF, p-value: 3.684e-14 > 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,] 8.343082e-03 1.668616e-02 9.916569e-01 [2,] 1.818705e-03 3.637410e-03 9.981813e-01 [3,] 2.870938e-04 5.741877e-04 9.997129e-01 [4,] 1.308132e-04 2.616263e-04 9.998692e-01 [5,] 2.831672e-05 5.663345e-05 9.999717e-01 [6,] 8.960722e-06 1.792144e-05 9.999910e-01 [7,] 6.170766e-06 1.234153e-05 9.999938e-01 [8,] 1.141329e-06 2.282657e-06 9.999989e-01 [9,] 2.477402e-07 4.954805e-07 9.999998e-01 [10,] 2.255083e-07 4.510167e-07 9.999998e-01 [11,] 2.798179e-06 5.596359e-06 9.999972e-01 [12,] 5.966764e-06 1.193353e-05 9.999940e-01 [13,] 4.415200e-06 8.830400e-06 9.999956e-01 [14,] 1.370926e-05 2.741851e-05 9.999863e-01 [15,] 1.214120e-04 2.428239e-04 9.998786e-01 [16,] 1.252004e-04 2.504008e-04 9.998748e-01 [17,] 6.998534e-05 1.399707e-04 9.999300e-01 [18,] 4.233559e-05 8.467119e-05 9.999577e-01 [19,] 1.673335e-05 3.346669e-05 9.999833e-01 [20,] 8.031001e-06 1.606200e-05 9.999920e-01 [21,] 2.424689e-05 4.849377e-05 9.999758e-01 [22,] 5.179204e-05 1.035841e-04 9.999482e-01 [23,] 1.892060e-04 3.784121e-04 9.998108e-01 [24,] 1.612037e-03 3.224073e-03 9.983880e-01 [25,] 1.354144e-02 2.708287e-02 9.864586e-01 [26,] 4.220069e-02 8.440137e-02 9.577993e-01 [27,] 1.221215e-01 2.442430e-01 8.778785e-01 [28,] 3.584274e-01 7.168548e-01 6.415726e-01 [29,] 6.852266e-01 6.295469e-01 3.147734e-01 [30,] 7.873787e-01 4.252425e-01 2.126213e-01 [31,] 8.558057e-01 2.883886e-01 1.441943e-01 [32,] 8.949470e-01 2.101059e-01 1.050530e-01 [33,] 8.870476e-01 2.259048e-01 1.129524e-01 [34,] 8.588471e-01 2.823059e-01 1.411529e-01 [35,] 9.083411e-01 1.833179e-01 9.165895e-02 [36,] 9.534032e-01 9.319365e-02 4.659683e-02 [37,] 9.856870e-01 2.862593e-02 1.431296e-02 [38,] 9.935630e-01 1.287405e-02 6.437026e-03 [39,] 9.936912e-01 1.261768e-02 6.308838e-03 [40,] 9.965920e-01 6.816096e-03 3.408048e-03 [41,] 9.993214e-01 1.357147e-03 6.785734e-04 [42,] 9.998704e-01 2.592372e-04 1.296186e-04 [43,] 9.999534e-01 9.310237e-05 4.655118e-05 [44,] 9.999982e-01 3.661097e-06 1.830548e-06 [45,] 9.999995e-01 9.917948e-07 4.958974e-07 [46,] 9.999987e-01 2.653119e-06 1.326559e-06 [47,] 9.999979e-01 4.233684e-06 2.116842e-06 [48,] 9.999981e-01 3.823993e-06 1.911996e-06 [49,] 9.999998e-01 3.541120e-07 1.770560e-07 [50,] 9.999993e-01 1.471433e-06 7.357164e-07 [51,] 9.999998e-01 3.162103e-07 1.581052e-07 [52,] 9.999999e-01 2.323608e-07 1.161804e-07 [53,] 9.999997e-01 5.610374e-07 2.805187e-07 [54,] 9.999985e-01 3.007153e-06 1.503576e-06 [55,] 9.999922e-01 1.559608e-05 7.798039e-06 [56,] 9.999663e-01 6.746387e-05 3.373194e-05 [57,] 9.998065e-01 3.870964e-04 1.935482e-04 [58,] 9.995745e-01 8.509873e-04 4.254936e-04 [59,] 9.972368e-01 5.526499e-03 2.763249e-03 > postscript(file="/var/www/html/rcomp/tmp/1kq6z1290542319.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/rcomp/tmp/2kq6z1290542319.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/rcomp/tmp/3vzn21290542319.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/rcomp/tmp/4vzn21290542319.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/rcomp/tmp/5vzn21290542319.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 = 70 Frequency = 1 1 2 3 4 5 6 -62.237868 -118.819561 -258.303591 -305.254477 -160.872687 -288.114689 7 8 9 10 11 12 -272.505582 -369.370271 -555.501397 -530.941185 -582.032647 -658.038330 13 14 15 16 17 18 531.967685 423.636682 237.937101 10.196648 -6.892181 24.748395 19 20 21 22 23 24 -169.827934 -295.404977 -114.238098 -189.044666 -303.855787 -179.999905 25 26 27 28 29 30 30.777038 -519.838597 -459.861726 -407.007470 16.631384 316.391177 31 32 33 34 35 36 381.642186 291.704685 340.752519 418.110837 589.398855 574.900257 37 38 39 40 41 42 864.266203 689.215859 549.012980 903.866548 969.748856 1059.056730 43 44 45 46 47 48 967.676563 718.505967 889.808231 784.048610 614.280943 544.171125 49 50 51 52 53 54 440.401195 283.200038 218.359535 70.571863 -114.882890 232.084471 55 56 57 58 59 60 261.728049 215.699876 113.553690 -56.159110 -273.814484 -451.778454 61 62 63 64 65 66 -548.444965 -604.074443 -655.267922 -757.193488 -821.165117 -897.003625 67 68 69 70 -842.430593 -861.018221 -892.349437 -994.506405 > postscript(file="/var/www/html/rcomp/tmp/66q451290542319.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 = 70 Frequency = 1 lag(myerror, k = 1) myerror 0 -62.237868 NA 1 -118.819561 -62.237868 2 -258.303591 -118.819561 3 -305.254477 -258.303591 4 -160.872687 -305.254477 5 -288.114689 -160.872687 6 -272.505582 -288.114689 7 -369.370271 -272.505582 8 -555.501397 -369.370271 9 -530.941185 -555.501397 10 -582.032647 -530.941185 11 -658.038330 -582.032647 12 531.967685 -658.038330 13 423.636682 531.967685 14 237.937101 423.636682 15 10.196648 237.937101 16 -6.892181 10.196648 17 24.748395 -6.892181 18 -169.827934 24.748395 19 -295.404977 -169.827934 20 -114.238098 -295.404977 21 -189.044666 -114.238098 22 -303.855787 -189.044666 23 -179.999905 -303.855787 24 30.777038 -179.999905 25 -519.838597 30.777038 26 -459.861726 -519.838597 27 -407.007470 -459.861726 28 16.631384 -407.007470 29 316.391177 16.631384 30 381.642186 316.391177 31 291.704685 381.642186 32 340.752519 291.704685 33 418.110837 340.752519 34 589.398855 418.110837 35 574.900257 589.398855 36 864.266203 574.900257 37 689.215859 864.266203 38 549.012980 689.215859 39 903.866548 549.012980 40 969.748856 903.866548 41 1059.056730 969.748856 42 967.676563 1059.056730 43 718.505967 967.676563 44 889.808231 718.505967 45 784.048610 889.808231 46 614.280943 784.048610 47 544.171125 614.280943 48 440.401195 544.171125 49 283.200038 440.401195 50 218.359535 283.200038 51 70.571863 218.359535 52 -114.882890 70.571863 53 232.084471 -114.882890 54 261.728049 232.084471 55 215.699876 261.728049 56 113.553690 215.699876 57 -56.159110 113.553690 58 -273.814484 -56.159110 59 -451.778454 -273.814484 60 -548.444965 -451.778454 61 -604.074443 -548.444965 62 -655.267922 -604.074443 63 -757.193488 -655.267922 64 -821.165117 -757.193488 65 -897.003625 -821.165117 66 -842.430593 -897.003625 67 -861.018221 -842.430593 68 -892.349437 -861.018221 69 -994.506405 -892.349437 70 NA -994.506405 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] -118.819561 -62.237868 [2,] -258.303591 -118.819561 [3,] -305.254477 -258.303591 [4,] -160.872687 -305.254477 [5,] -288.114689 -160.872687 [6,] -272.505582 -288.114689 [7,] -369.370271 -272.505582 [8,] -555.501397 -369.370271 [9,] -530.941185 -555.501397 [10,] -582.032647 -530.941185 [11,] -658.038330 -582.032647 [12,] 531.967685 -658.038330 [13,] 423.636682 531.967685 [14,] 237.937101 423.636682 [15,] 10.196648 237.937101 [16,] -6.892181 10.196648 [17,] 24.748395 -6.892181 [18,] -169.827934 24.748395 [19,] -295.404977 -169.827934 [20,] -114.238098 -295.404977 [21,] -189.044666 -114.238098 [22,] -303.855787 -189.044666 [23,] -179.999905 -303.855787 [24,] 30.777038 -179.999905 [25,] -519.838597 30.777038 [26,] -459.861726 -519.838597 [27,] -407.007470 -459.861726 [28,] 16.631384 -407.007470 [29,] 316.391177 16.631384 [30,] 381.642186 316.391177 [31,] 291.704685 381.642186 [32,] 340.752519 291.704685 [33,] 418.110837 340.752519 [34,] 589.398855 418.110837 [35,] 574.900257 589.398855 [36,] 864.266203 574.900257 [37,] 689.215859 864.266203 [38,] 549.012980 689.215859 [39,] 903.866548 549.012980 [40,] 969.748856 903.866548 [41,] 1059.056730 969.748856 [42,] 967.676563 1059.056730 [43,] 718.505967 967.676563 [44,] 889.808231 718.505967 [45,] 784.048610 889.808231 [46,] 614.280943 784.048610 [47,] 544.171125 614.280943 [48,] 440.401195 544.171125 [49,] 283.200038 440.401195 [50,] 218.359535 283.200038 [51,] 70.571863 218.359535 [52,] -114.882890 70.571863 [53,] 232.084471 -114.882890 [54,] 261.728049 232.084471 [55,] 215.699876 261.728049 [56,] 113.553690 215.699876 [57,] -56.159110 113.553690 [58,] -273.814484 -56.159110 [59,] -451.778454 -273.814484 [60,] -548.444965 -451.778454 [61,] -604.074443 -548.444965 [62,] -655.267922 -604.074443 [63,] -757.193488 -655.267922 [64,] -821.165117 -757.193488 [65,] -897.003625 -821.165117 [66,] -842.430593 -897.003625 [67,] -861.018221 -842.430593 [68,] -892.349437 -861.018221 [69,] -994.506405 -892.349437 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 -118.819561 -62.237868 2 -258.303591 -118.819561 3 -305.254477 -258.303591 4 -160.872687 -305.254477 5 -288.114689 -160.872687 6 -272.505582 -288.114689 7 -369.370271 -272.505582 8 -555.501397 -369.370271 9 -530.941185 -555.501397 10 -582.032647 -530.941185 11 -658.038330 -582.032647 12 531.967685 -658.038330 13 423.636682 531.967685 14 237.937101 423.636682 15 10.196648 237.937101 16 -6.892181 10.196648 17 24.748395 -6.892181 18 -169.827934 24.748395 19 -295.404977 -169.827934 20 -114.238098 -295.404977 21 -189.044666 -114.238098 22 -303.855787 -189.044666 23 -179.999905 -303.855787 24 30.777038 -179.999905 25 -519.838597 30.777038 26 -459.861726 -519.838597 27 -407.007470 -459.861726 28 16.631384 -407.007470 29 316.391177 16.631384 30 381.642186 316.391177 31 291.704685 381.642186 32 340.752519 291.704685 33 418.110837 340.752519 34 589.398855 418.110837 35 574.900257 589.398855 36 864.266203 574.900257 37 689.215859 864.266203 38 549.012980 689.215859 39 903.866548 549.012980 40 969.748856 903.866548 41 1059.056730 969.748856 42 967.676563 1059.056730 43 718.505967 967.676563 44 889.808231 718.505967 45 784.048610 889.808231 46 614.280943 784.048610 47 544.171125 614.280943 48 440.401195 544.171125 49 283.200038 440.401195 50 218.359535 283.200038 51 70.571863 218.359535 52 -114.882890 70.571863 53 232.084471 -114.882890 54 261.728049 232.084471 55 215.699876 261.728049 56 113.553690 215.699876 57 -56.159110 113.553690 58 -273.814484 -56.159110 59 -451.778454 -273.814484 60 -548.444965 -451.778454 61 -604.074443 -548.444965 62 -655.267922 -604.074443 63 -757.193488 -655.267922 64 -821.165117 -757.193488 65 -897.003625 -821.165117 66 -842.430593 -897.003625 67 -861.018221 -842.430593 68 -892.349437 -861.018221 69 -994.506405 -892.349437 > 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/rcomp/tmp/7y0lq1290542319.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/rcomp/tmp/8y0lq1290542319.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/rcomp/tmp/9y0lq1290542319.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/rcomp/tmp/10rrlt1290542319.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/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/html/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/rcomp/tmp/11ds1h1290542319.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/rcomp/tmp/12ys041290542319.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/rcomp/tmp/13nbfg1290542319.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/rcomp/tmp/14xkwj1290542319.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/rcomp/tmp/151lvp1290542319.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/rcomp/tmp/16xuay1290542319.tab") + } > > try(system("convert tmp/1kq6z1290542319.ps tmp/1kq6z1290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/2kq6z1290542319.ps tmp/2kq6z1290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/3vzn21290542319.ps tmp/3vzn21290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/4vzn21290542319.ps tmp/4vzn21290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/5vzn21290542319.ps tmp/5vzn21290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/66q451290542319.ps tmp/66q451290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/7y0lq1290542319.ps tmp/7y0lq1290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/8y0lq1290542319.ps tmp/8y0lq1290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/9y0lq1290542319.ps tmp/9y0lq1290542319.png",intern=TRUE)) character(0) > try(system("convert tmp/10rrlt1290542319.ps tmp/10rrlt1290542319.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 2.505 1.568 5.859