R version 2.8.0 (2008-10-20) Copyright (C) 2008 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. 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. > x <- array(list(695 + ,0 + ,641 + ,696 + ,729 + ,839 + ,627 + ,608 + ,638 + ,0 + ,695 + ,641 + ,696 + ,729 + ,696 + ,651 + ,762 + ,0 + ,638 + ,695 + ,641 + ,696 + ,825 + ,691 + ,635 + ,0 + ,762 + ,638 + ,695 + ,641 + ,677 + ,627 + ,721 + ,0 + ,635 + ,762 + ,638 + ,695 + ,656 + ,634 + ,854 + ,0 + ,721 + ,635 + ,762 + ,638 + ,785 + ,731 + ,418 + ,0 + ,854 + ,721 + ,635 + ,762 + ,412 + ,475 + ,367 + ,0 + ,418 + ,854 + ,721 + ,635 + ,352 + ,337 + ,824 + ,0 + ,367 + ,418 + ,854 + ,721 + ,839 + ,803 + ,687 + ,0 + ,824 + ,367 + ,418 + ,854 + ,729 + ,722 + ,601 + ,0 + ,687 + ,824 + ,367 + ,418 + ,696 + ,590 + ,676 + ,0 + ,601 + ,687 + ,824 + ,367 + ,641 + ,724 + ,740 + ,0 + ,676 + ,601 + ,687 + ,824 + ,695 + ,627 + ,691 + ,0 + ,740 + ,676 + ,601 + ,687 + ,638 + ,696 + ,683 + ,0 + ,691 + ,740 + ,676 + ,601 + ,762 + ,825 + ,594 + ,0 + ,683 + ,691 + ,740 + ,676 + ,635 + ,677 + ,729 + ,0 + ,594 + ,683 + ,691 + ,740 + ,721 + ,656 + ,731 + ,0 + ,729 + ,594 + ,683 + ,691 + ,854 + ,785 + ,386 + ,0 + ,731 + ,729 + ,594 + ,683 + ,418 + ,412 + ,331 + ,0 + ,386 + ,731 + ,729 + ,594 + ,367 + ,352 + ,706 + ,0 + ,331 + ,386 + ,731 + ,729 + ,824 + ,839 + ,715 + ,0 + ,706 + ,331 + ,386 + ,731 + ,687 + ,729 + ,657 + ,0 + ,715 + ,706 + ,331 + ,386 + ,601 + ,696 + ,653 + ,0 + ,657 + ,715 + ,706 + ,331 + ,676 + ,641 + ,642 + ,0 + ,653 + ,657 + ,715 + ,706 + ,740 + ,695 + ,643 + ,0 + ,642 + ,653 + ,657 + ,715 + ,691 + ,638 + ,718 + ,0 + ,643 + ,642 + ,653 + ,657 + ,683 + ,762 + ,654 + ,0 + ,718 + ,643 + ,642 + ,653 + ,594 + ,635 + ,632 + ,0 + ,654 + ,718 + ,643 + ,642 + ,729 + ,721 + ,731 + ,0 + ,632 + ,654 + ,718 + ,643 + ,731 + ,854 + ,392 + ,1 + ,731 + ,632 + ,654 + ,718 + ,386 + ,418 + ,344 + ,1 + ,392 + ,731 + ,632 + ,654 + ,331 + ,367 + ,792 + ,1 + ,344 + ,392 + ,731 + ,632 + ,706 + ,824 + ,852 + ,1 + ,792 + ,344 + ,392 + ,731 + ,715 + ,687 + ,649 + ,1 + ,852 + ,792 + ,344 + ,392 + ,657 + ,601 + ,629 + ,1 + ,649 + ,852 + ,792 + ,344 + ,653 + ,676 + ,685 + ,1 + ,629 + ,649 + ,852 + ,792 + ,642 + ,740 + ,617 + ,1 + ,685 + ,629 + ,649 + ,852 + ,643 + ,691 + ,715 + ,1 + ,617 + ,685 + ,629 + ,649 + ,718 + ,683 + ,715 + ,1 + ,715 + ,617 + ,685 + ,629 + ,654 + ,594 + ,629 + ,1 + ,715 + ,715 + ,617 + ,685 + ,632 + ,729 + ,916 + ,1 + ,629 + ,715 + ,715 + ,617 + ,731 + ,731 + ,531 + ,1 + ,916 + ,629 + ,715 + ,715 + ,392 + ,386 + ,357 + ,1 + ,531 + ,916 + ,629 + ,715 + ,344 + ,331 + ,917 + ,1 + ,357 + ,531 + ,916 + ,629 + ,792 + ,706 + ,828 + ,1 + ,917 + ,357 + ,531 + ,916 + ,852 + ,715 + ,708 + ,1 + ,828 + ,917 + ,357 + ,531 + ,649 + ,657 + ,858 + ,1 + ,708 + ,828 + ,917 + ,357 + ,629 + ,653 + ,775 + ,1 + ,858 + ,708 + ,828 + ,917 + ,685 + ,642 + ,785 + ,1 + ,775 + ,858 + ,708 + ,828 + ,617 + ,643 + ,1006 + ,1 + ,785 + ,775 + ,858 + ,708 + ,715 + ,718 + ,789 + ,1 + ,1006 + ,785 + ,775 + ,858 + ,715 + ,654 + ,734 + ,1 + ,789 + ,1006 + ,785 + ,775 + ,629 + ,632 + ,906 + ,1 + ,734 + ,789 + ,1006 + ,785 + ,916 + ,731 + ,532 + ,1 + ,906 + ,734 + ,789 + ,1006 + ,531 + ,392 + ,387 + ,1 + ,532 + ,906 + ,734 + ,789 + ,357 + ,344 + ,991 + ,1 + ,387 + ,532 + ,906 + ,734 + ,917 + ,792 + ,841 + ,1 + ,991 + ,387 + ,532 + ,906 + ,828 + ,852 + ,892 + ,1 + ,841 + ,991 + ,387 + ,532 + ,708 + ,649 + ,782 + ,1 + ,892 + ,841 + ,991 + ,387 + ,858 + ,629) + ,dim=c(8 + ,60) + ,dimnames=list(c('faillissement' + ,'crisis' + ,'t-1' + ,'t-2' + ,'t-3' + ,'t-4' + ,'t-12' + ,'t-24') + ,1:60)) > y <- array(NA,dim=c(8,60),dimnames=list(c('faillissement','crisis','t-1','t-2','t-3','t-4','t-12','t-24'),1:60)) > 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' > #'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.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 faillissement crisis t-1 t-2 t-3 t-4 t-12 t-24 1 695 0 641 696 729 839 627 608 2 638 0 695 641 696 729 696 651 3 762 0 638 695 641 696 825 691 4 635 0 762 638 695 641 677 627 5 721 0 635 762 638 695 656 634 6 854 0 721 635 762 638 785 731 7 418 0 854 721 635 762 412 475 8 367 0 418 854 721 635 352 337 9 824 0 367 418 854 721 839 803 10 687 0 824 367 418 854 729 722 11 601 0 687 824 367 418 696 590 12 676 0 601 687 824 367 641 724 13 740 0 676 601 687 824 695 627 14 691 0 740 676 601 687 638 696 15 683 0 691 740 676 601 762 825 16 594 0 683 691 740 676 635 677 17 729 0 594 683 691 740 721 656 18 731 0 729 594 683 691 854 785 19 386 0 731 729 594 683 418 412 20 331 0 386 731 729 594 367 352 21 706 0 331 386 731 729 824 839 22 715 0 706 331 386 731 687 729 23 657 0 715 706 331 386 601 696 24 653 0 657 715 706 331 676 641 25 642 0 653 657 715 706 740 695 26 643 0 642 653 657 715 691 638 27 718 0 643 642 653 657 683 762 28 654 0 718 643 642 653 594 635 29 632 0 654 718 643 642 729 721 30 731 0 632 654 718 643 731 854 31 392 1 731 632 654 718 386 418 32 344 1 392 731 632 654 331 367 33 792 1 344 392 731 632 706 824 34 852 1 792 344 392 731 715 687 35 649 1 852 792 344 392 657 601 36 629 1 649 852 792 344 653 676 37 685 1 629 649 852 792 642 740 38 617 1 685 629 649 852 643 691 39 715 1 617 685 629 649 718 683 40 715 1 715 617 685 629 654 594 41 629 1 715 715 617 685 632 729 42 916 1 629 715 715 617 731 731 43 531 1 916 629 715 715 392 386 44 357 1 531 916 629 715 344 331 45 917 1 357 531 916 629 792 706 46 828 1 917 357 531 916 852 715 47 708 1 828 917 357 531 649 657 48 858 1 708 828 917 357 629 653 49 775 1 858 708 828 917 685 642 50 785 1 775 858 708 828 617 643 51 1006 1 785 775 858 708 715 718 52 789 1 1006 785 775 858 715 654 53 734 1 789 1006 785 775 629 632 54 906 1 734 789 1006 785 916 731 55 532 1 906 734 789 1006 531 392 56 387 1 532 906 734 789 357 344 57 991 1 387 532 906 734 917 792 58 841 1 991 387 532 906 828 852 59 892 1 841 991 387 532 708 649 60 782 1 892 841 991 387 858 629 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) crisis `t-1` `t-2` `t-3` `t-4` -1.227e+02 7.846e+01 5.238e-02 4.446e-02 6.345e-02 -7.502e-03 `t-12` `t-24` 6.689e-01 3.456e-01 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -116.658 -48.727 -7.443 50.785 199.203 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -1.227e+02 1.186e+02 -1.034 0.305735 crisis 7.846e+01 2.048e+01 3.832 0.000344 *** `t-1` 5.238e-02 7.240e-02 0.723 0.472633 `t-2` 4.446e-02 8.217e-02 0.541 0.590807 `t-3` 6.345e-02 7.297e-02 0.869 0.388572 `t-4` -7.502e-03 7.504e-02 -0.100 0.920745 `t-12` 6.689e-01 1.432e-01 4.671 2.15e-05 *** `t-24` 3.456e-01 1.540e-01 2.244 0.029146 * --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 73.88 on 52 degrees of freedom Multiple R-squared: 0.81, Adjusted R-squared: 0.7844 F-statistic: 31.66 on 7 and 52 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,] 0.36629683 0.73259366 0.6337032 [2,] 0.23551812 0.47103624 0.7644819 [3,] 0.20249549 0.40499097 0.7975045 [4,] 0.11822235 0.23644471 0.8817776 [5,] 0.30626662 0.61253323 0.6937334 [6,] 0.27030215 0.54060431 0.7296978 [7,] 0.19340368 0.38680736 0.8065963 [8,] 0.23894033 0.47788066 0.7610597 [9,] 0.18082212 0.36164423 0.8191779 [10,] 0.14058047 0.28116094 0.8594195 [11,] 0.13692774 0.27385548 0.8630723 [12,] 0.16341215 0.32682430 0.8365879 [13,] 0.16314844 0.32629689 0.8368516 [14,] 0.11299041 0.22598083 0.8870096 [15,] 0.10634885 0.21269769 0.8936512 [16,] 0.07175346 0.14350692 0.9282465 [17,] 0.04872630 0.09745260 0.9512737 [18,] 0.04593059 0.09186119 0.9540694 [19,] 0.04173304 0.08346608 0.9582670 [20,] 0.02585658 0.05171317 0.9741434 [21,] 0.01633876 0.03267752 0.9836612 [22,] 0.01003576 0.02007152 0.9899642 [23,] 0.00843734 0.01687468 0.9915627 [24,] 0.01629443 0.03258886 0.9837056 [25,] 0.01024478 0.02048956 0.9897552 [26,] 0.01733295 0.03466590 0.9826671 [27,] 0.01568979 0.03137957 0.9843102 [28,] 0.02334075 0.04668150 0.9766592 [29,] 0.01890207 0.03780414 0.9810979 [30,] 0.01198951 0.02397903 0.9880105 [31,] 0.07650263 0.15300526 0.9234974 [32,] 0.13511061 0.27022122 0.8648894 [33,] 0.12981659 0.25963317 0.8701834 [34,] 0.10191041 0.20382083 0.8980896 [35,] 0.08266987 0.16533974 0.9173301 [36,] 0.06395201 0.12790402 0.9360480 [37,] 0.06346432 0.12692865 0.9365357 [38,] 0.05213773 0.10427545 0.9478623 [39,] 0.02453573 0.04907145 0.9754643 > postscript(file="/var/www/html/freestat/rcomp/tmp/19fz61292773340.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/www/html/freestat/rcomp/tmp/2jog91292773340.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/www/html/freestat/rcomp/tmp/3jog91292773340.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/www/html/freestat/rcomp/tmp/4jog91292773340.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/www/html/freestat/rcomp/tmp/5ugfc1292773340.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 = 60 Frequency = 1 1 2 3 4 5 6 83.758830 -33.365883 -5.643042 -19.337167 83.451061 89.496384 7 8 9 10 11 12 -10.361028 36.973105 21.473527 -6.972854 -38.462714 8.237041 13 14 15 16 17 18 81.653576 44.676563 -96.520253 -50.333969 43.008045 -91.500986 19 20 21 22 23 24 -16.508207 -7.913894 -87.761622 55.588808 50.272162 -6.452926 25 26 27 28 29 30 -73.889264 -15.917318 21.839546 57.948148 -84.192660 -33.243010 31 32 33 34 35 36 -68.873711 -48.191978 2.207338 104.447629 -52.598082 -116.657682 37 38 39 40 41 42 -65.789585 -106.240559 -54.821658 12.926260 -104.631665 113.237781 43 44 45 46 47 48 63.721336 -46.304269 91.841603 -35.417215 -11.680455 126.485873 49 50 51 52 53 54 17.156986 76.918266 199.202854 -1.310675 9.097358 -46.491865 55 56 57 58 59 60 -36.978336 -35.206089 52.323488 -59.054592 129.755187 -109.073541 > postscript(file="/var/www/html/freestat/rcomp/tmp/6ugfc1292773340.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 = 60 Frequency = 1 lag(myerror, k = 1) myerror 0 83.758830 NA 1 -33.365883 83.758830 2 -5.643042 -33.365883 3 -19.337167 -5.643042 4 83.451061 -19.337167 5 89.496384 83.451061 6 -10.361028 89.496384 7 36.973105 -10.361028 8 21.473527 36.973105 9 -6.972854 21.473527 10 -38.462714 -6.972854 11 8.237041 -38.462714 12 81.653576 8.237041 13 44.676563 81.653576 14 -96.520253 44.676563 15 -50.333969 -96.520253 16 43.008045 -50.333969 17 -91.500986 43.008045 18 -16.508207 -91.500986 19 -7.913894 -16.508207 20 -87.761622 -7.913894 21 55.588808 -87.761622 22 50.272162 55.588808 23 -6.452926 50.272162 24 -73.889264 -6.452926 25 -15.917318 -73.889264 26 21.839546 -15.917318 27 57.948148 21.839546 28 -84.192660 57.948148 29 -33.243010 -84.192660 30 -68.873711 -33.243010 31 -48.191978 -68.873711 32 2.207338 -48.191978 33 104.447629 2.207338 34 -52.598082 104.447629 35 -116.657682 -52.598082 36 -65.789585 -116.657682 37 -106.240559 -65.789585 38 -54.821658 -106.240559 39 12.926260 -54.821658 40 -104.631665 12.926260 41 113.237781 -104.631665 42 63.721336 113.237781 43 -46.304269 63.721336 44 91.841603 -46.304269 45 -35.417215 91.841603 46 -11.680455 -35.417215 47 126.485873 -11.680455 48 17.156986 126.485873 49 76.918266 17.156986 50 199.202854 76.918266 51 -1.310675 199.202854 52 9.097358 -1.310675 53 -46.491865 9.097358 54 -36.978336 -46.491865 55 -35.206089 -36.978336 56 52.323488 -35.206089 57 -59.054592 52.323488 58 129.755187 -59.054592 59 -109.073541 129.755187 60 NA -109.073541 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] -33.365883 83.758830 [2,] -5.643042 -33.365883 [3,] -19.337167 -5.643042 [4,] 83.451061 -19.337167 [5,] 89.496384 83.451061 [6,] -10.361028 89.496384 [7,] 36.973105 -10.361028 [8,] 21.473527 36.973105 [9,] -6.972854 21.473527 [10,] -38.462714 -6.972854 [11,] 8.237041 -38.462714 [12,] 81.653576 8.237041 [13,] 44.676563 81.653576 [14,] -96.520253 44.676563 [15,] -50.333969 -96.520253 [16,] 43.008045 -50.333969 [17,] -91.500986 43.008045 [18,] -16.508207 -91.500986 [19,] -7.913894 -16.508207 [20,] -87.761622 -7.913894 [21,] 55.588808 -87.761622 [22,] 50.272162 55.588808 [23,] -6.452926 50.272162 [24,] -73.889264 -6.452926 [25,] -15.917318 -73.889264 [26,] 21.839546 -15.917318 [27,] 57.948148 21.839546 [28,] -84.192660 57.948148 [29,] -33.243010 -84.192660 [30,] -68.873711 -33.243010 [31,] -48.191978 -68.873711 [32,] 2.207338 -48.191978 [33,] 104.447629 2.207338 [34,] -52.598082 104.447629 [35,] -116.657682 -52.598082 [36,] -65.789585 -116.657682 [37,] -106.240559 -65.789585 [38,] -54.821658 -106.240559 [39,] 12.926260 -54.821658 [40,] -104.631665 12.926260 [41,] 113.237781 -104.631665 [42,] 63.721336 113.237781 [43,] -46.304269 63.721336 [44,] 91.841603 -46.304269 [45,] -35.417215 91.841603 [46,] -11.680455 -35.417215 [47,] 126.485873 -11.680455 [48,] 17.156986 126.485873 [49,] 76.918266 17.156986 [50,] 199.202854 76.918266 [51,] -1.310675 199.202854 [52,] 9.097358 -1.310675 [53,] -46.491865 9.097358 [54,] -36.978336 -46.491865 [55,] -35.206089 -36.978336 [56,] 52.323488 -35.206089 [57,] -59.054592 52.323488 [58,] 129.755187 -59.054592 [59,] -109.073541 129.755187 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 -33.365883 83.758830 2 -5.643042 -33.365883 3 -19.337167 -5.643042 4 83.451061 -19.337167 5 89.496384 83.451061 6 -10.361028 89.496384 7 36.973105 -10.361028 8 21.473527 36.973105 9 -6.972854 21.473527 10 -38.462714 -6.972854 11 8.237041 -38.462714 12 81.653576 8.237041 13 44.676563 81.653576 14 -96.520253 44.676563 15 -50.333969 -96.520253 16 43.008045 -50.333969 17 -91.500986 43.008045 18 -16.508207 -91.500986 19 -7.913894 -16.508207 20 -87.761622 -7.913894 21 55.588808 -87.761622 22 50.272162 55.588808 23 -6.452926 50.272162 24 -73.889264 -6.452926 25 -15.917318 -73.889264 26 21.839546 -15.917318 27 57.948148 21.839546 28 -84.192660 57.948148 29 -33.243010 -84.192660 30 -68.873711 -33.243010 31 -48.191978 -68.873711 32 2.207338 -48.191978 33 104.447629 2.207338 34 -52.598082 104.447629 35 -116.657682 -52.598082 36 -65.789585 -116.657682 37 -106.240559 -65.789585 38 -54.821658 -106.240559 39 12.926260 -54.821658 40 -104.631665 12.926260 41 113.237781 -104.631665 42 63.721336 113.237781 43 -46.304269 63.721336 44 91.841603 -46.304269 45 -35.417215 91.841603 46 -11.680455 -35.417215 47 126.485873 -11.680455 48 17.156986 126.485873 49 76.918266 17.156986 50 199.202854 76.918266 51 -1.310675 199.202854 52 9.097358 -1.310675 53 -46.491865 9.097358 54 -36.978336 -46.491865 55 -35.206089 -36.978336 56 52.323488 -35.206089 57 -59.054592 52.323488 58 129.755187 -59.054592 59 -109.073541 129.755187 > 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/freestat/rcomp/tmp/75pwf1292773340.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/www/html/freestat/rcomp/tmp/85pwf1292773340.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/www/html/freestat/rcomp/tmp/9fyei1292773340.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/www/html/freestat/rcomp/tmp/10fyei1292773340.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/www/html/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/html/freestat/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/freestat/rcomp/tmp/111hc51292773340.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/freestat/rcomp/tmp/124zbt1292773340.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/freestat/rcomp/tmp/1309921292773340.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/freestat/rcomp/tmp/14mrp81292773340.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/freestat/rcomp/tmp/157ane1292773340.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/freestat/rcomp/tmp/16as4k1292773340.tab") + } > > try(system("convert tmp/19fz61292773340.ps tmp/19fz61292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/2jog91292773340.ps tmp/2jog91292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/3jog91292773340.ps tmp/3jog91292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/4jog91292773340.ps tmp/4jog91292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/5ugfc1292773340.ps tmp/5ugfc1292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/6ugfc1292773340.ps tmp/6ugfc1292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/75pwf1292773340.ps tmp/75pwf1292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/85pwf1292773340.ps tmp/85pwf1292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/9fyei1292773340.ps tmp/9fyei1292773340.png",intern=TRUE)) character(0) > try(system("convert tmp/10fyei1292773340.ps tmp/10fyei1292773340.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 3.968 2.529 4.380