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(527070 + ,250643 + ,276427 + ,45390 + ,114971 + ,208512 + ,0 + ,509846 + ,243422 + ,266424 + ,39215 + ,105531 + ,205708 + ,0 + ,514258 + ,247105 + ,267153 + ,38433 + ,104919 + ,206890 + ,0 + ,516922 + ,248541 + ,268381 + ,37676 + ,104782 + ,207069 + ,1 + ,507561 + ,245039 + ,262522 + ,36055 + ,101281 + ,205305 + ,1 + ,492622 + ,237080 + ,255542 + ,32986 + ,94545 + ,201504 + ,1 + ,490243 + ,237085 + ,253158 + ,30953 + ,93248 + ,200517 + ,1 + ,469357 + ,225554 + ,243803 + ,23558 + ,84031 + ,195771 + ,1 + ,477580 + ,226839 + ,250741 + ,22487 + ,87486 + ,195259 + ,1 + ,528379 + ,247934 + ,280445 + ,43528 + ,115867 + ,197579 + ,1 + ,533590 + ,248333 + ,285257 + ,47913 + ,120327 + ,196985 + ,1 + ,517945 + ,246969 + ,270976 + ,48621 + ,117008 + ,194382 + ,1 + ,506174 + ,245098 + ,261076 + ,42169 + ,108811 + ,191580 + ,1 + ,501866 + ,246263 + ,255603 + ,38444 + ,104519 + ,190765 + ,1 + ,516141 + ,255765 + ,260376 + ,38692 + ,106758 + ,191480 + ,1 + ,528222 + ,264319 + ,263903 + ,38124 + ,109337 + ,192277 + ,1 + ,532638 + ,268347 + ,264291 + ,37886 + ,109078 + ,191632 + ,1 + ,536322 + ,273046 + ,263276 + ,37310 + ,108293 + ,190757 + ,1 + ,536535 + ,273963 + ,262572 + ,34689 + ,106534 + ,190995 + ,1 + ,523597 + ,267430 + ,256167 + ,26450 + ,99197 + ,189081 + ,1 + ,536214 + ,271993 + ,264221 + ,25565 + ,103493 + ,190028 + ,1 + ,586570 + ,292710 + ,293860 + ,46562 + ,130676 + ,196146 + ,1 + ,596594 + ,295881 + ,300713 + ,52653 + ,137448 + ,197070 + ,1 + ,580523 + ,293299 + ,287224 + ,54807 + ,134704 + ,194893 + ,1 + ,564478 + ,288576 + ,275902 + ,47534 + ,123725 + ,193246 + ,1 + ,557560 + ,286445 + ,271115 + ,43565 + ,118277 + ,192484 + ,1 + ,575093 + ,297584 + ,277509 + ,44051 + ,121225 + ,194924 + ,1 + ,580112 + ,300431 + ,279681 + ,42622 + ,120528 + ,197394 + ,1 + ,574761 + ,298522 + ,276239 + ,41761 + ,118240 + ,196598 + ,1 + ,563250 + ,292213 + ,271037 + ,39086 + ,112514 + ,194409 + ,1 + ,551531 + ,285383 + ,266148 + ,35438 + ,107304 + ,193431 + ,1 + ,537034 + ,277537 + ,259497 + ,27356 + ,100001 + ,191942 + ,1 + ,544686 + ,277891 + ,266795 + ,26149 + ,102082 + ,193323 + ,1 + ,600991 + ,302686 + ,298305 + ,47034 + ,130455 + ,199654 + ,1 + ,604378 + ,300653 + ,303725 + ,53091 + ,135574 + ,198422 + ,1 + ,586111 + ,296369 + ,289742 + ,55718 + ,132540 + ,198219 + ,1 + ,563668 + ,287224 + ,276444 + ,47637 + ,119920 + ,197157 + ,1 + ,548604 + ,279998 + ,268606 + ,43237 + ,112454 + ,195115 + ,1 + ,551174 + ,283495 + ,267679 + ,40597 + ,109415 + ,197296 + ,1 + ,555654 + ,285775 + ,269879 + ,39884 + ,109843 + ,198178 + ,0 + ,547970 + ,282329 + ,265641 + ,38504 + ,106365 + ,197787 + ,0 + ,540324 + ,277799 + ,262525 + ,36393 + ,102304 + ,197622 + ,0 + ,530577 + ,271980 + ,258597 + ,33740 + ,97968 + ,196683 + ,0 + ,520579 + ,266730 + ,253849 + ,26131 + ,92462 + ,194590 + ,0 + ,518654 + ,262433 + ,256221 + ,23885 + ,92286 + ,194316 + ,0 + ,572273 + ,285378 + ,286895 + ,43899 + ,120092 + ,199598 + ,0 + ,581302 + ,286692 + ,294610 + ,49871 + ,126656 + ,199055 + ,0 + ,563280 + ,282917 + ,280363 + ,52292 + ,124144 + ,197482 + ,0 + ,547612 + ,277686 + ,269926 + ,45493 + ,114045 + ,196440 + ,0 + ,538712 + ,274371 + ,264341 + ,41124 + ,108120 + ,195338 + ,0 + ,540735 + ,277466 + ,263269 + ,39385 + ,105698 + ,195589 + ,0 + ,561649 + ,290604 + ,271045 + ,41472 + ,111203 + ,198936 + ,0 + ,558685 + ,290770 + ,267915 + ,41688 + ,110030 + ,198262 + ,0 + ,545732 + ,283654 + ,262078 + ,38711 + ,104009 + ,197275 + ,0 + ,536352 + ,278601 + ,257751 + ,36840 + ,99772 + ,196007 + ,0 + ,527676 + ,274405 + ,253271 + ,35141 + ,96301 + ,194447 + ,0 + ,530455 + ,272817 + ,257638 + ,37443 + ,97680 + ,193951 + ,0 + ,581744 + ,294292 + ,287452 + ,51905 + ,121563 + ,198396 + ,0 + ,598714 + ,300562 + ,298152 + ,60016 + ,134210 + ,199486 + ,0 + ,583775 + ,298982 + ,284793 + ,58611 + ,133111 + ,198688 + ,0 + ,571477 + ,296917 + ,274560 + ,52097 + ,124527 + ,196729 + ,0) + ,dim=c(7 + ,61) + ,dimnames=list(c('Totaal' + ,'mannen' + ,'vrouwen' + ,'Beroepsinschakeling' + ,'jongerdan25jaar' + ,'meerdan2jaarinactief' + ,'Dummy') + ,1:61)) > y <- array(NA,dim=c(7,61),dimnames=list(c('Totaal','mannen','vrouwen','Beroepsinschakeling','jongerdan25jaar','meerdan2jaarinactief','Dummy'),1:61)) > 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' > par3 <- 'No Linear Trend' > par2 <- 'Do not include Seasonal Dummies' > par1 <- '1' > #'GNU S' R Code compiled by R2WASP v. 1.0.44 () > #Author: Prof. Dr. P. Wessa > #To cite this work: AUTHOR(S), (YEAR), YOUR SOFTWARE TITLE (vNUMBER) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_YOURPAGE.wasp/ > #Source of accompanying publication: Office for Research, Development, and Education > #Technical description: Write here your technical program description (don't use hard returns!) > library(lattice) > library(lmtest) Loading required package: zoo Attaching package: 'zoo' The following object(s) are masked from 'package:base': as.Date, 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 Totaal mannen vrouwen Beroepsinschakeling jongerdan25jaar 1 527070 250643 276427 45390 114971 2 509846 243422 266424 39215 105531 3 514258 247105 267153 38433 104919 4 516922 248541 268381 37676 104782 5 507561 245039 262522 36055 101281 6 492622 237080 255542 32986 94545 7 490243 237085 253158 30953 93248 8 469357 225554 243803 23558 84031 9 477580 226839 250741 22487 87486 10 528379 247934 280445 43528 115867 11 533590 248333 285257 47913 120327 12 517945 246969 270976 48621 117008 13 506174 245098 261076 42169 108811 14 501866 246263 255603 38444 104519 15 516141 255765 260376 38692 106758 16 528222 264319 263903 38124 109337 17 532638 268347 264291 37886 109078 18 536322 273046 263276 37310 108293 19 536535 273963 262572 34689 106534 20 523597 267430 256167 26450 99197 21 536214 271993 264221 25565 103493 22 586570 292710 293860 46562 130676 23 596594 295881 300713 52653 137448 24 580523 293299 287224 54807 134704 25 564478 288576 275902 47534 123725 26 557560 286445 271115 43565 118277 27 575093 297584 277509 44051 121225 28 580112 300431 279681 42622 120528 29 574761 298522 276239 41761 118240 30 563250 292213 271037 39086 112514 31 551531 285383 266148 35438 107304 32 537034 277537 259497 27356 100001 33 544686 277891 266795 26149 102082 34 600991 302686 298305 47034 130455 35 604378 300653 303725 53091 135574 36 586111 296369 289742 55718 132540 37 563668 287224 276444 47637 119920 38 548604 279998 268606 43237 112454 39 551174 283495 267679 40597 109415 40 555654 285775 269879 39884 109843 41 547970 282329 265641 38504 106365 42 540324 277799 262525 36393 102304 43 530577 271980 258597 33740 97968 44 520579 266730 253849 26131 92462 45 518654 262433 256221 23885 92286 46 572273 285378 286895 43899 120092 47 581302 286692 294610 49871 126656 48 563280 282917 280363 52292 124144 49 547612 277686 269926 45493 114045 50 538712 274371 264341 41124 108120 51 540735 277466 263269 39385 105698 52 561649 290604 271045 41472 111203 53 558685 290770 267915 41688 110030 54 545732 283654 262078 38711 104009 55 536352 278601 257751 36840 99772 56 527676 274405 253271 35141 96301 57 530455 272817 257638 37443 97680 58 581744 294292 287452 51905 121563 59 598714 300562 298152 60016 134210 60 583775 298982 284793 58611 133111 61 571477 296917 274560 52097 124527 meerdan2jaarinactief Dummy 1 208512 0 2 205708 0 3 206890 0 4 207069 1 5 205305 1 6 201504 1 7 200517 1 8 195771 1 9 195259 1 10 197579 1 11 196985 1 12 194382 1 13 191580 1 14 190765 1 15 191480 1 16 192277 1 17 191632 1 18 190757 1 19 190995 1 20 189081 1 21 190028 1 22 196146 1 23 197070 1 24 194893 1 25 193246 1 26 192484 1 27 194924 1 28 197394 1 29 196598 1 30 194409 1 31 193431 1 32 191942 1 33 193323 1 34 199654 1 35 198422 1 36 198219 1 37 197157 1 38 195115 1 39 197296 1 40 198178 0 41 197787 0 42 197622 0 43 196683 0 44 194590 0 45 194316 0 46 199598 0 47 199055 0 48 197482 0 49 196440 0 50 195338 0 51 195589 0 52 198936 0 53 198262 0 54 197275 0 55 196007 0 56 194447 0 57 193951 0 58 198396 0 59 199486 0 60 198688 0 61 196729 0 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) mannen vrouwen -5.962e-11 1.000e+00 1.000e+00 Beroepsinschakeling jongerdan25jaar meerdan2jaarinactief 2.372e-16 -3.513e-16 -7.188e-17 Dummy 1.444e-12 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -6.678e-12 -1.354e-12 -3.081e-13 1.500e-12 5.241e-12 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -5.962e-11 2.045e-11 -2.915e+00 0.00517 ** mannen 1.000e+00 2.719e-17 3.678e+16 < 2e-16 *** vrouwen 1.000e+00 9.926e-17 1.007e+16 < 2e-16 *** Beroepsinschakeling 2.372e-16 1.416e-16 1.675e+00 0.09971 . jongerdan25jaar -3.512e-16 1.794e-16 -1.958e+00 0.05536 . meerdan2jaarinactief -7.188e-17 1.112e-16 -6.470e-01 0.52065 Dummy 1.444e-12 9.351e-13 1.544e+00 0.12832 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 2.449e-12 on 54 degrees of freedom Multiple R-squared: 1, Adjusted R-squared: 1 F-statistic: 1.589e+33 on 6 and 54 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,] 4.019758e-01 8.039516e-01 5.980242e-01 [2,] 3.309909e-06 6.619818e-06 9.999967e-01 [3,] 1.000000e+00 1.344595e-26 6.722973e-27 [4,] 5.404590e-01 9.190819e-01 4.595410e-01 [5,] 1.000000e+00 2.501134e-28 1.250567e-28 [6,] 3.763636e-01 7.527271e-01 6.236364e-01 [7,] 5.289410e-01 9.421181e-01 4.710590e-01 [8,] 1.000000e+00 1.139422e-45 5.697109e-46 [9,] 1.000000e+00 4.562804e-38 2.281402e-38 [10,] 1.017762e-01 2.035524e-01 8.982238e-01 [11,] 1.000000e+00 5.789558e-45 2.894779e-45 [12,] 9.383009e-01 1.233983e-01 6.169914e-02 [13,] 8.710840e-01 2.578320e-01 1.289160e-01 [14,] 8.575866e-01 2.848269e-01 1.424134e-01 [15,] 1.000000e+00 7.154469e-32 3.577235e-32 [16,] 1.000000e+00 1.331811e-26 6.659054e-27 [17,] 9.476467e-01 1.047066e-01 5.235331e-02 [18,] 1.000000e+00 1.535393e-30 7.676967e-31 [19,] 1.000000e+00 8.664181e-33 4.332091e-33 [20,] 1.000000e+00 9.193583e-30 4.596791e-30 [21,] 9.985515e-01 2.896985e-03 1.448492e-03 [22,] 1.906372e-23 3.812745e-23 1.000000e+00 [23,] 1.000000e+00 1.015100e-17 5.075501e-18 [24,] 2.841332e-23 5.682663e-23 1.000000e+00 [25,] 6.958119e-01 6.083762e-01 3.041881e-01 [26,] 1.079034e-01 2.158069e-01 8.920966e-01 [27,] 2.098048e-26 4.196095e-26 1.000000e+00 [28,] 2.980692e-01 5.961384e-01 7.019308e-01 [29,] 1.000000e+00 4.015504e-22 2.007752e-22 [30,] 1.326084e-27 2.652168e-27 1.000000e+00 [31,] 8.046341e-01 3.907318e-01 1.953659e-01 [32,] 1.000000e+00 4.887440e-14 2.443720e-14 [33,] 1.000000e+00 7.753754e-14 3.876877e-14 [34,] 4.800971e-03 9.601941e-03 9.951990e-01 [35,] 3.580671e-35 7.161341e-35 1.000000e+00 [36,] 1.000000e+00 8.280378e-11 4.140189e-11 [37,] 2.275682e-01 4.551364e-01 7.724318e-01 [38,] 2.910346e-01 5.820693e-01 7.089654e-01 [39,] 8.551109e-29 1.710222e-28 1.000000e+00 [40,] 2.846743e-16 5.693486e-16 1.000000e+00 [41,] 8.148950e-01 3.702100e-01 1.851050e-01 [42,] 9.316591e-01 1.366818e-01 6.834089e-02 > postscript(file="/var/fisher/rcomp/tmp/1jzpl1353438037.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/fisher/rcomp/tmp/2oobs1353438037.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/fisher/rcomp/tmp/3s1bj1353438037.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/fisher/rcomp/tmp/47i6e1353438037.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/fisher/rcomp/tmp/51d5d1353438037.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 = 61 Frequency = 1 1 2 3 4 5 -4.066207e-12 -6.678020e-12 4.557405e-12 -1.411847e-12 2.051201e-12 6 7 8 9 10 3.659277e-12 3.312086e-12 -2.167893e-12 5.240785e-12 -3.711853e-12 11 12 13 14 15 -1.818615e-12 -4.073219e-12 4.300940e-12 3.948944e-12 -1.869265e-12 16 17 18 19 20 -3.438743e-13 -5.264969e-13 -1.836707e-12 4.317024e-13 -3.361651e-13 21 22 23 24 25 -5.802759e-14 2.046324e-12 2.087701e-13 2.147844e-12 -1.459355e-12 26 27 28 29 30 2.734801e-13 -2.803857e-12 2.476219e-12 -5.900185e-13 -2.780251e-12 31 32 33 34 35 -2.199732e-12 -1.552296e-14 -6.028713e-13 1.892493e-12 -3.080851e-13 36 37 38 39 40 2.637635e-12 -1.615429e-12 -3.105720e-12 -9.928938e-13 2.317882e-13 41 42 43 44 45 1.352370e-13 -6.633047e-13 7.848654e-14 -6.365548e-13 1.005989e-12 46 47 48 49 50 -7.480978e-13 1.698357e-12 -1.212275e-12 -8.897575e-14 6.556738e-13 51 52 53 54 55 -8.707463e-13 -1.120423e-12 -1.354190e-12 3.219984e-12 1.070859e-12 56 57 58 59 60 -9.602619e-13 -9.888844e-13 1.068978e-12 1.499906e-12 3.663461e-12 61 5.018155e-13 > postscript(file="/var/fisher/rcomp/tmp/6f4vi1353438037.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 = 61 Frequency = 1 lag(myerror, k = 1) myerror 0 -4.066207e-12 NA 1 -6.678020e-12 -4.066207e-12 2 4.557405e-12 -6.678020e-12 3 -1.411847e-12 4.557405e-12 4 2.051201e-12 -1.411847e-12 5 3.659277e-12 2.051201e-12 6 3.312086e-12 3.659277e-12 7 -2.167893e-12 3.312086e-12 8 5.240785e-12 -2.167893e-12 9 -3.711853e-12 5.240785e-12 10 -1.818615e-12 -3.711853e-12 11 -4.073219e-12 -1.818615e-12 12 4.300940e-12 -4.073219e-12 13 3.948944e-12 4.300940e-12 14 -1.869265e-12 3.948944e-12 15 -3.438743e-13 -1.869265e-12 16 -5.264969e-13 -3.438743e-13 17 -1.836707e-12 -5.264969e-13 18 4.317024e-13 -1.836707e-12 19 -3.361651e-13 4.317024e-13 20 -5.802759e-14 -3.361651e-13 21 2.046324e-12 -5.802759e-14 22 2.087701e-13 2.046324e-12 23 2.147844e-12 2.087701e-13 24 -1.459355e-12 2.147844e-12 25 2.734801e-13 -1.459355e-12 26 -2.803857e-12 2.734801e-13 27 2.476219e-12 -2.803857e-12 28 -5.900185e-13 2.476219e-12 29 -2.780251e-12 -5.900185e-13 30 -2.199732e-12 -2.780251e-12 31 -1.552296e-14 -2.199732e-12 32 -6.028713e-13 -1.552296e-14 33 1.892493e-12 -6.028713e-13 34 -3.080851e-13 1.892493e-12 35 2.637635e-12 -3.080851e-13 36 -1.615429e-12 2.637635e-12 37 -3.105720e-12 -1.615429e-12 38 -9.928938e-13 -3.105720e-12 39 2.317882e-13 -9.928938e-13 40 1.352370e-13 2.317882e-13 41 -6.633047e-13 1.352370e-13 42 7.848654e-14 -6.633047e-13 43 -6.365548e-13 7.848654e-14 44 1.005989e-12 -6.365548e-13 45 -7.480978e-13 1.005989e-12 46 1.698357e-12 -7.480978e-13 47 -1.212275e-12 1.698357e-12 48 -8.897575e-14 -1.212275e-12 49 6.556738e-13 -8.897575e-14 50 -8.707463e-13 6.556738e-13 51 -1.120423e-12 -8.707463e-13 52 -1.354190e-12 -1.120423e-12 53 3.219984e-12 -1.354190e-12 54 1.070859e-12 3.219984e-12 55 -9.602619e-13 1.070859e-12 56 -9.888844e-13 -9.602619e-13 57 1.068978e-12 -9.888844e-13 58 1.499906e-12 1.068978e-12 59 3.663461e-12 1.499906e-12 60 5.018155e-13 3.663461e-12 61 NA 5.018155e-13 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] -6.678020e-12 -4.066207e-12 [2,] 4.557405e-12 -6.678020e-12 [3,] -1.411847e-12 4.557405e-12 [4,] 2.051201e-12 -1.411847e-12 [5,] 3.659277e-12 2.051201e-12 [6,] 3.312086e-12 3.659277e-12 [7,] -2.167893e-12 3.312086e-12 [8,] 5.240785e-12 -2.167893e-12 [9,] -3.711853e-12 5.240785e-12 [10,] -1.818615e-12 -3.711853e-12 [11,] -4.073219e-12 -1.818615e-12 [12,] 4.300940e-12 -4.073219e-12 [13,] 3.948944e-12 4.300940e-12 [14,] -1.869265e-12 3.948944e-12 [15,] -3.438743e-13 -1.869265e-12 [16,] -5.264969e-13 -3.438743e-13 [17,] -1.836707e-12 -5.264969e-13 [18,] 4.317024e-13 -1.836707e-12 [19,] -3.361651e-13 4.317024e-13 [20,] -5.802759e-14 -3.361651e-13 [21,] 2.046324e-12 -5.802759e-14 [22,] 2.087701e-13 2.046324e-12 [23,] 2.147844e-12 2.087701e-13 [24,] -1.459355e-12 2.147844e-12 [25,] 2.734801e-13 -1.459355e-12 [26,] -2.803857e-12 2.734801e-13 [27,] 2.476219e-12 -2.803857e-12 [28,] -5.900185e-13 2.476219e-12 [29,] -2.780251e-12 -5.900185e-13 [30,] -2.199732e-12 -2.780251e-12 [31,] -1.552296e-14 -2.199732e-12 [32,] -6.028713e-13 -1.552296e-14 [33,] 1.892493e-12 -6.028713e-13 [34,] -3.080851e-13 1.892493e-12 [35,] 2.637635e-12 -3.080851e-13 [36,] -1.615429e-12 2.637635e-12 [37,] -3.105720e-12 -1.615429e-12 [38,] -9.928938e-13 -3.105720e-12 [39,] 2.317882e-13 -9.928938e-13 [40,] 1.352370e-13 2.317882e-13 [41,] -6.633047e-13 1.352370e-13 [42,] 7.848654e-14 -6.633047e-13 [43,] -6.365548e-13 7.848654e-14 [44,] 1.005989e-12 -6.365548e-13 [45,] -7.480978e-13 1.005989e-12 [46,] 1.698357e-12 -7.480978e-13 [47,] -1.212275e-12 1.698357e-12 [48,] -8.897575e-14 -1.212275e-12 [49,] 6.556738e-13 -8.897575e-14 [50,] -8.707463e-13 6.556738e-13 [51,] -1.120423e-12 -8.707463e-13 [52,] -1.354190e-12 -1.120423e-12 [53,] 3.219984e-12 -1.354190e-12 [54,] 1.070859e-12 3.219984e-12 [55,] -9.602619e-13 1.070859e-12 [56,] -9.888844e-13 -9.602619e-13 [57,] 1.068978e-12 -9.888844e-13 [58,] 1.499906e-12 1.068978e-12 [59,] 3.663461e-12 1.499906e-12 [60,] 5.018155e-13 3.663461e-12 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 -6.678020e-12 -4.066207e-12 2 4.557405e-12 -6.678020e-12 3 -1.411847e-12 4.557405e-12 4 2.051201e-12 -1.411847e-12 5 3.659277e-12 2.051201e-12 6 3.312086e-12 3.659277e-12 7 -2.167893e-12 3.312086e-12 8 5.240785e-12 -2.167893e-12 9 -3.711853e-12 5.240785e-12 10 -1.818615e-12 -3.711853e-12 11 -4.073219e-12 -1.818615e-12 12 4.300940e-12 -4.073219e-12 13 3.948944e-12 4.300940e-12 14 -1.869265e-12 3.948944e-12 15 -3.438743e-13 -1.869265e-12 16 -5.264969e-13 -3.438743e-13 17 -1.836707e-12 -5.264969e-13 18 4.317024e-13 -1.836707e-12 19 -3.361651e-13 4.317024e-13 20 -5.802759e-14 -3.361651e-13 21 2.046324e-12 -5.802759e-14 22 2.087701e-13 2.046324e-12 23 2.147844e-12 2.087701e-13 24 -1.459355e-12 2.147844e-12 25 2.734801e-13 -1.459355e-12 26 -2.803857e-12 2.734801e-13 27 2.476219e-12 -2.803857e-12 28 -5.900185e-13 2.476219e-12 29 -2.780251e-12 -5.900185e-13 30 -2.199732e-12 -2.780251e-12 31 -1.552296e-14 -2.199732e-12 32 -6.028713e-13 -1.552296e-14 33 1.892493e-12 -6.028713e-13 34 -3.080851e-13 1.892493e-12 35 2.637635e-12 -3.080851e-13 36 -1.615429e-12 2.637635e-12 37 -3.105720e-12 -1.615429e-12 38 -9.928938e-13 -3.105720e-12 39 2.317882e-13 -9.928938e-13 40 1.352370e-13 2.317882e-13 41 -6.633047e-13 1.352370e-13 42 7.848654e-14 -6.633047e-13 43 -6.365548e-13 7.848654e-14 44 1.005989e-12 -6.365548e-13 45 -7.480978e-13 1.005989e-12 46 1.698357e-12 -7.480978e-13 47 -1.212275e-12 1.698357e-12 48 -8.897575e-14 -1.212275e-12 49 6.556738e-13 -8.897575e-14 50 -8.707463e-13 6.556738e-13 51 -1.120423e-12 -8.707463e-13 52 -1.354190e-12 -1.120423e-12 53 3.219984e-12 -1.354190e-12 54 1.070859e-12 3.219984e-12 55 -9.602619e-13 1.070859e-12 56 -9.888844e-13 -9.602619e-13 57 1.068978e-12 -9.888844e-13 58 1.499906e-12 1.068978e-12 59 3.663461e-12 1.499906e-12 60 5.018155e-13 3.663461e-12 > 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/fisher/rcomp/tmp/7th6j1353438037.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/fisher/rcomp/tmp/8t0nw1353438037.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/fisher/rcomp/tmp/9ssof1353438037.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/fisher/rcomp/tmp/10ik711353438037.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/fisher/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/fisher/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/fisher/rcomp/tmp/1188k91353438037.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/fisher/rcomp/tmp/12epdi1353438037.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/fisher/rcomp/tmp/13grdn1353438037.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/fisher/rcomp/tmp/14hnwj1353438037.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/fisher/rcomp/tmp/15fagd1353438037.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/fisher/rcomp/tmp/16aed71353438037.tab") + } > > try(system("convert tmp/1jzpl1353438037.ps tmp/1jzpl1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/2oobs1353438037.ps tmp/2oobs1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/3s1bj1353438037.ps tmp/3s1bj1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/47i6e1353438037.ps tmp/47i6e1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/51d5d1353438037.ps tmp/51d5d1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/6f4vi1353438037.ps tmp/6f4vi1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/7th6j1353438037.ps tmp/7th6j1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/8t0nw1353438037.ps tmp/8t0nw1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/9ssof1353438037.ps tmp/9ssof1353438037.png",intern=TRUE)) character(0) > try(system("convert tmp/10ik711353438037.ps tmp/10ik711353438037.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 6.702 1.492 8.184