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(474605,0,470390,0,461251,0,454724,0,455626,0,516847,0,525192,0,522975,0,518585,0,509239,0,512238,0,519164,0,517009,0,509933,0,509127,0,500875,0,506971,0,569323,0,579714,0,577992,0,565644,0,547344,0,554788,0,562325,0,560854,0,555332,0,543599,0,536662,0,542722,0,593530,0,610763,0,612613,0,611324,0,594167,0,595454,0,590865,0,589379,0,584428,0,573100,0,567456,0,569028,0,620735,0,628884,0,628232,0,612117,0,595404,0,597141,0,593408,0,590072,0,579799,0,574205,0,572775,0,572942,0,619567,0,625809,0,619916,0,587625,0,565724,0,557274,0,560576,0,548854,0,531673,0,525919,0,511038,0,498662,0,555362,0,564591,0,541667,0,527070,0,509846,0,514258,0,516922,0,507561,0,492622,0,490243,0,469357,0,477580,0,528379,0,533590,0,517945,1,506174,1,501866,1,516441,1,528222,1,532638,1),dim=c(2,85),dimnames=list(c('Werkzoekend','Crisis'),1:85))
>  y <- array(NA,dim=c(2,85),dimnames=list(c('Werkzoekend','Crisis'),1:85))
>  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
   Werkzoekend Crisis
1       474605      0
2       470390      0
3       461251      0
4       454724      0
5       455626      0
6       516847      0
7       525192      0
8       522975      0
9       518585      0
10      509239      0
11      512238      0
12      519164      0
13      517009      0
14      509933      0
15      509127      0
16      500875      0
17      506971      0
18      569323      0
19      579714      0
20      577992      0
21      565644      0
22      547344      0
23      554788      0
24      562325      0
25      560854      0
26      555332      0
27      543599      0
28      536662      0
29      542722      0
30      593530      0
31      610763      0
32      612613      0
33      611324      0
34      594167      0
35      595454      0
36      590865      0
37      589379      0
38      584428      0
39      573100      0
40      567456      0
41      569028      0
42      620735      0
43      628884      0
44      628232      0
45      612117      0
46      595404      0
47      597141      0
48      593408      0
49      590072      0
50      579799      0
51      574205      0
52      572775      0
53      572942      0
54      619567      0
55      625809      0
56      619916      0
57      587625      0
58      565724      0
59      557274      0
60      560576      0
61      548854      0
62      531673      0
63      525919      0
64      511038      0
65      498662      0
66      555362      0
67      564591      0
68      541667      0
69      527070      0
70      509846      0
71      514258      0
72      516922      0
73      507561      0
74      492622      0
75      490243      0
76      469357      0
77      477580      0
78      528379      0
79      533590      0
80      517945      1
81      506174      1
82      501866      1
83      516441      1
84      528222      1
85      532638      1
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept)       Crisis  
     549627       -32413  
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
   Min     1Q Median     3Q    Max 
-94903 -32705   5161  30172  79257 
Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   549627       4944 111.160   <2e-16 ***
Crisis        -32413      18610  -1.742   0.0853 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 
Residual standard error: 43950 on 83 degrees of freedom
Multiple R-squared: 0.03526,	Adjusted R-squared: 0.02363 
F-statistic: 3.033 on 1 and 83 DF,  p-value: 0.08527 
> 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.02325514 0.0465102882 0.9767448559
 [2,] 0.21937802 0.4387560336 0.7806219832
 [3,] 0.34444049 0.6888809704 0.6555595148
 [4,] 0.36688453 0.7337690599 0.6331154700
 [5,] 0.33598494 0.6719698767 0.6640150617
 [6,] 0.27133995 0.5426798958 0.7286600521
 [7,] 0.21991951 0.4398390159 0.7800804921
 [8,] 0.18753582 0.3750716368 0.8124641816
 [9,] 0.15264413 0.3052882601 0.8473558700
[10,] 0.11594214 0.2318842889 0.8840578555
[11,] 0.08687852 0.1737570389 0.9131214805
[12,] 0.06454551 0.1290910142 0.9354544929
[13,] 0.04795771 0.0959154157 0.9520422921
[14,] 0.14714243 0.2942848660 0.8528575670
[15,] 0.31495096 0.6299019130 0.6850490435
[16,] 0.44598986 0.8919797216 0.5540101392
[17,] 0.48416383 0.9683276656 0.5158361672
[18,] 0.45403680 0.9080736056 0.5459631972
[19,] 0.43723705 0.8744741080 0.5627629460
[20,] 0.43509240 0.8701848018 0.5649075991
[21,] 0.42180501 0.8436100154 0.5781949923
[22,] 0.39158674 0.7831734825 0.6084132588
[23,] 0.34453026 0.6890605237 0.6554697381
[24,] 0.29654080 0.5930815934 0.7034592033
[25,] 0.25436453 0.5087290519 0.7456354741
[26,] 0.32767174 0.6553434808 0.6723282596
[27,] 0.47093031 0.9418606222 0.5290696889
[28,] 0.59976417 0.8004716585 0.4002358292
[29,] 0.69484289 0.6103142181 0.3051571090
[30,] 0.71393896 0.5721220783 0.2860610392
[31,] 0.73109920 0.5378016056 0.2689008028
[32,] 0.73262561 0.5347487750 0.2673743875
[33,] 0.72814617 0.5437076539 0.2718538269
[34,] 0.71125117 0.5774976639 0.2887488320
[35,] 0.67330382 0.6533923613 0.3266961806
[36,] 0.62568950 0.7486209903 0.3743104952
[37,] 0.57748233 0.8450353490 0.4225176745
[38,] 0.67019241 0.6596151759 0.3298075879
[39,] 0.78161426 0.4367714875 0.2183857437
[40,] 0.86693229 0.2661354162 0.1330677081
[41,] 0.89851038 0.2029792482 0.1014896241
[42,] 0.90182086 0.1963582792 0.0981791396
[43,] 0.90862977 0.1827404631 0.0913702315
[44,] 0.91179934 0.1764013119 0.0882006559
[45,] 0.91232914 0.1753417285 0.0876708642
[46,] 0.90257345 0.1948530983 0.0974265492
[47,] 0.88710691 0.2257861776 0.1128930888
[48,] 0.86938912 0.2612217571 0.1306108785
[49,] 0.85126894 0.2974621177 0.1487310588
[50,] 0.92822313 0.1435537376 0.0717768688
[51,] 0.98303491 0.0339301813 0.0169650906
[52,] 0.99808867 0.0038226588 0.0019113294
[53,] 0.99936686 0.0012662789 0.0006331394
[54,] 0.99951793 0.0009641497 0.0004820748
[55,] 0.99953812 0.0009237509 0.0004618755
[56,] 0.99967024 0.0006595174 0.0003297587
[57,] 0.99964927 0.0007014555 0.0003507277
[58,] 0.99941128 0.0011774316 0.0005887158
[59,] 0.99895344 0.0020931255 0.0010465627
[60,] 0.99816645 0.0036670982 0.0018335491
[61,] 0.99744155 0.0051168970 0.0025584485
[62,] 0.99815056 0.0036988796 0.0018494398
[63,] 0.99960115 0.0007977100 0.0003988550
[64,] 0.99972592 0.0005481502 0.0002740751
[65,] 0.99962634 0.0007473274 0.0003736637
[66,] 0.99917752 0.0016449638 0.0008224819
[67,] 0.99836723 0.0032655473 0.0016327736
[68,] 0.99713612 0.0057277570 0.0028638785
[69,] 0.99410308 0.0117938433 0.0058969217
[70,] 0.98829605 0.0234078920 0.0117039460
[71,] 0.97848106 0.0430378846 0.0215189423
[72,] 0.98843526 0.0231294800 0.0115647400
[73,] 0.99919328 0.0016134464 0.0008067232
[74,] 0.99663608 0.0067278337 0.0033639169
[75,] 0.98636617 0.0272676545 0.0136338273
[76,] 0.94999121 0.1000175713 0.0500087857
> postscript(file="/var/www/html/rcomp/tmp/1jcr21258750534.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/2uy5j1258750534.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/3vbl41258750534.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/4oysb1258750534.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/5slv71258750534.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 = 85 
Frequency = 1 
          1           2           3           4           5           6 
-75022.2911 -79237.2911 -88376.2911 -94903.2911 -94001.2911 -32780.2911 
          7           8           9          10          11          12 
-24435.2911 -26652.2911 -31042.2911 -40388.2911 -37389.2911 -30463.2911 
         13          14          15          16          17          18 
-32618.2911 -39694.2911 -40500.2911 -48752.2911 -42656.2911  19695.7089 
         19          20          21          22          23          24 
 30086.7089  28364.7089  16016.7089  -2283.2911   5160.7089  12697.7089 
         25          26          27          28          29          30 
 11226.7089   5704.7089  -6028.2911 -12965.2911  -6905.2911  43902.7089 
         31          32          33          34          35          36 
 61135.7089  62985.7089  61696.7089  44539.7089  45826.7089  41237.7089 
         37          38          39          40          41          42 
 39751.7089  34800.7089  23472.7089  17828.7089  19400.7089  71107.7089 
         43          44          45          46          47          48 
 79256.7089  78604.7089  62489.7089  45776.7089  47513.7089  43780.7089 
         49          50          51          52          53          54 
 40444.7089  30171.7089  24577.7089  23147.7089  23314.7089  69939.7089 
         55          56          57          58          59          60 
 76181.7089  70288.7089  37997.7089  16096.7089   7646.7089  10948.7089 
         61          62          63          64          65          66 
  -773.2911 -17954.2911 -23708.2911 -38589.2911 -50965.2911   5734.7089 
         67          68          69          70          71          72 
 14963.7089  -7960.2911 -22557.2911 -39781.2911 -35369.2911 -32705.2911 
         73          74          75          76          77          78 
-42066.2911 -57005.2911 -59384.2911 -80270.2911 -72047.2911 -21248.2911 
         79          80          81          82          83          84 
-16037.2911    730.6667 -11040.3333 -15348.3333   -773.3333  11007.6667 
         85 
 15423.6667 
> postscript(file="/var/www/html/rcomp/tmp/6959i1258750534.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 = 85 
Frequency = 1 
   lag(myerror, k = 1)     myerror
 0         -75022.2911          NA
 1         -79237.2911 -75022.2911
 2         -88376.2911 -79237.2911
 3         -94903.2911 -88376.2911
 4         -94001.2911 -94903.2911
 5         -32780.2911 -94001.2911
 6         -24435.2911 -32780.2911
 7         -26652.2911 -24435.2911
 8         -31042.2911 -26652.2911
 9         -40388.2911 -31042.2911
10         -37389.2911 -40388.2911
11         -30463.2911 -37389.2911
12         -32618.2911 -30463.2911
13         -39694.2911 -32618.2911
14         -40500.2911 -39694.2911
15         -48752.2911 -40500.2911
16         -42656.2911 -48752.2911
17          19695.7089 -42656.2911
18          30086.7089  19695.7089
19          28364.7089  30086.7089
20          16016.7089  28364.7089
21          -2283.2911  16016.7089
22           5160.7089  -2283.2911
23          12697.7089   5160.7089
24          11226.7089  12697.7089
25           5704.7089  11226.7089
26          -6028.2911   5704.7089
27         -12965.2911  -6028.2911
28          -6905.2911 -12965.2911
29          43902.7089  -6905.2911
30          61135.7089  43902.7089
31          62985.7089  61135.7089
32          61696.7089  62985.7089
33          44539.7089  61696.7089
34          45826.7089  44539.7089
35          41237.7089  45826.7089
36          39751.7089  41237.7089
37          34800.7089  39751.7089
38          23472.7089  34800.7089
39          17828.7089  23472.7089
40          19400.7089  17828.7089
41          71107.7089  19400.7089
42          79256.7089  71107.7089
43          78604.7089  79256.7089
44          62489.7089  78604.7089
45          45776.7089  62489.7089
46          47513.7089  45776.7089
47          43780.7089  47513.7089
48          40444.7089  43780.7089
49          30171.7089  40444.7089
50          24577.7089  30171.7089
51          23147.7089  24577.7089
52          23314.7089  23147.7089
53          69939.7089  23314.7089
54          76181.7089  69939.7089
55          70288.7089  76181.7089
56          37997.7089  70288.7089
57          16096.7089  37997.7089
58           7646.7089  16096.7089
59          10948.7089   7646.7089
60           -773.2911  10948.7089
61         -17954.2911   -773.2911
62         -23708.2911 -17954.2911
63         -38589.2911 -23708.2911
64         -50965.2911 -38589.2911
65           5734.7089 -50965.2911
66          14963.7089   5734.7089
67          -7960.2911  14963.7089
68         -22557.2911  -7960.2911
69         -39781.2911 -22557.2911
70         -35369.2911 -39781.2911
71         -32705.2911 -35369.2911
72         -42066.2911 -32705.2911
73         -57005.2911 -42066.2911
74         -59384.2911 -57005.2911
75         -80270.2911 -59384.2911
76         -72047.2911 -80270.2911
77         -21248.2911 -72047.2911
78         -16037.2911 -21248.2911
79            730.6667 -16037.2911
80         -11040.3333    730.6667
81         -15348.3333 -11040.3333
82           -773.3333 -15348.3333
83          11007.6667   -773.3333
84          15423.6667  11007.6667
85                  NA  15423.6667
> dum1 <- dum[2:length(myerror),]
> dum1
      lag(myerror, k = 1)     myerror
 [1,]         -79237.2911 -75022.2911
 [2,]         -88376.2911 -79237.2911
 [3,]         -94903.2911 -88376.2911
 [4,]         -94001.2911 -94903.2911
 [5,]         -32780.2911 -94001.2911
 [6,]         -24435.2911 -32780.2911
 [7,]         -26652.2911 -24435.2911
 [8,]         -31042.2911 -26652.2911
 [9,]         -40388.2911 -31042.2911
[10,]         -37389.2911 -40388.2911
[11,]         -30463.2911 -37389.2911
[12,]         -32618.2911 -30463.2911
[13,]         -39694.2911 -32618.2911
[14,]         -40500.2911 -39694.2911
[15,]         -48752.2911 -40500.2911
[16,]         -42656.2911 -48752.2911
[17,]          19695.7089 -42656.2911
[18,]          30086.7089  19695.7089
[19,]          28364.7089  30086.7089
[20,]          16016.7089  28364.7089
[21,]          -2283.2911  16016.7089
[22,]           5160.7089  -2283.2911
[23,]          12697.7089   5160.7089
[24,]          11226.7089  12697.7089
[25,]           5704.7089  11226.7089
[26,]          -6028.2911   5704.7089
[27,]         -12965.2911  -6028.2911
[28,]          -6905.2911 -12965.2911
[29,]          43902.7089  -6905.2911
[30,]          61135.7089  43902.7089
[31,]          62985.7089  61135.7089
[32,]          61696.7089  62985.7089
[33,]          44539.7089  61696.7089
[34,]          45826.7089  44539.7089
[35,]          41237.7089  45826.7089
[36,]          39751.7089  41237.7089
[37,]          34800.7089  39751.7089
[38,]          23472.7089  34800.7089
[39,]          17828.7089  23472.7089
[40,]          19400.7089  17828.7089
[41,]          71107.7089  19400.7089
[42,]          79256.7089  71107.7089
[43,]          78604.7089  79256.7089
[44,]          62489.7089  78604.7089
[45,]          45776.7089  62489.7089
[46,]          47513.7089  45776.7089
[47,]          43780.7089  47513.7089
[48,]          40444.7089  43780.7089
[49,]          30171.7089  40444.7089
[50,]          24577.7089  30171.7089
[51,]          23147.7089  24577.7089
[52,]          23314.7089  23147.7089
[53,]          69939.7089  23314.7089
[54,]          76181.7089  69939.7089
[55,]          70288.7089  76181.7089
[56,]          37997.7089  70288.7089
[57,]          16096.7089  37997.7089
[58,]           7646.7089  16096.7089
[59,]          10948.7089   7646.7089
[60,]           -773.2911  10948.7089
[61,]         -17954.2911   -773.2911
[62,]         -23708.2911 -17954.2911
[63,]         -38589.2911 -23708.2911
[64,]         -50965.2911 -38589.2911
[65,]           5734.7089 -50965.2911
[66,]          14963.7089   5734.7089
[67,]          -7960.2911  14963.7089
[68,]         -22557.2911  -7960.2911
[69,]         -39781.2911 -22557.2911
[70,]         -35369.2911 -39781.2911
[71,]         -32705.2911 -35369.2911
[72,]         -42066.2911 -32705.2911
[73,]         -57005.2911 -42066.2911
[74,]         -59384.2911 -57005.2911
[75,]         -80270.2911 -59384.2911
[76,]         -72047.2911 -80270.2911
[77,]         -21248.2911 -72047.2911
[78,]         -16037.2911 -21248.2911
[79,]            730.6667 -16037.2911
[80,]         -11040.3333    730.6667
[81,]         -15348.3333 -11040.3333
[82,]           -773.3333 -15348.3333
[83,]          11007.6667   -773.3333
[84,]          15423.6667  11007.6667
> z <- as.data.frame(dum1)
> z
   lag(myerror, k = 1)     myerror
1          -79237.2911 -75022.2911
2          -88376.2911 -79237.2911
3          -94903.2911 -88376.2911
4          -94001.2911 -94903.2911
5          -32780.2911 -94001.2911
6          -24435.2911 -32780.2911
7          -26652.2911 -24435.2911
8          -31042.2911 -26652.2911
9          -40388.2911 -31042.2911
10         -37389.2911 -40388.2911
11         -30463.2911 -37389.2911
12         -32618.2911 -30463.2911
13         -39694.2911 -32618.2911
14         -40500.2911 -39694.2911
15         -48752.2911 -40500.2911
16         -42656.2911 -48752.2911
17          19695.7089 -42656.2911
18          30086.7089  19695.7089
19          28364.7089  30086.7089
20          16016.7089  28364.7089
21          -2283.2911  16016.7089
22           5160.7089  -2283.2911
23          12697.7089   5160.7089
24          11226.7089  12697.7089
25           5704.7089  11226.7089
26          -6028.2911   5704.7089
27         -12965.2911  -6028.2911
28          -6905.2911 -12965.2911
29          43902.7089  -6905.2911
30          61135.7089  43902.7089
31          62985.7089  61135.7089
32          61696.7089  62985.7089
33          44539.7089  61696.7089
34          45826.7089  44539.7089
35          41237.7089  45826.7089
36          39751.7089  41237.7089
37          34800.7089  39751.7089
38          23472.7089  34800.7089
39          17828.7089  23472.7089
40          19400.7089  17828.7089
41          71107.7089  19400.7089
42          79256.7089  71107.7089
43          78604.7089  79256.7089
44          62489.7089  78604.7089
45          45776.7089  62489.7089
46          47513.7089  45776.7089
47          43780.7089  47513.7089
48          40444.7089  43780.7089
49          30171.7089  40444.7089
50          24577.7089  30171.7089
51          23147.7089  24577.7089
52          23314.7089  23147.7089
53          69939.7089  23314.7089
54          76181.7089  69939.7089
55          70288.7089  76181.7089
56          37997.7089  70288.7089
57          16096.7089  37997.7089
58           7646.7089  16096.7089
59          10948.7089   7646.7089
60           -773.2911  10948.7089
61         -17954.2911   -773.2911
62         -23708.2911 -17954.2911
63         -38589.2911 -23708.2911
64         -50965.2911 -38589.2911
65           5734.7089 -50965.2911
66          14963.7089   5734.7089
67          -7960.2911  14963.7089
68         -22557.2911  -7960.2911
69         -39781.2911 -22557.2911
70         -35369.2911 -39781.2911
71         -32705.2911 -35369.2911
72         -42066.2911 -32705.2911
73         -57005.2911 -42066.2911
74         -59384.2911 -57005.2911
75         -80270.2911 -59384.2911
76         -72047.2911 -80270.2911
77         -21248.2911 -72047.2911
78         -16037.2911 -21248.2911
79            730.6667 -16037.2911
80         -11040.3333    730.6667
81         -15348.3333 -11040.3333
82           -773.3333 -15348.3333
83          11007.6667   -773.3333
84          15423.6667  11007.6667
> 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/77r361258750534.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/87j1s1258750534.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/9lbcy1258750534.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/10jj1p1258750534.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/11mr3h1258750534.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/125at01258750534.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/13xtya1258750534.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/14warb1258750534.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/152nhu1258750534.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/16r8231258750534.tab") 
+ }
> 
> system("convert tmp/1jcr21258750534.ps tmp/1jcr21258750534.png")
> system("convert tmp/2uy5j1258750534.ps tmp/2uy5j1258750534.png")
> system("convert tmp/3vbl41258750534.ps tmp/3vbl41258750534.png")
> system("convert tmp/4oysb1258750534.ps tmp/4oysb1258750534.png")
> system("convert tmp/5slv71258750534.ps tmp/5slv71258750534.png")
> system("convert tmp/6959i1258750534.ps tmp/6959i1258750534.png")
> system("convert tmp/77r361258750534.ps tmp/77r361258750534.png")
> system("convert tmp/87j1s1258750534.ps tmp/87j1s1258750534.png")
> system("convert tmp/9lbcy1258750534.ps tmp/9lbcy1258750534.png")
> system("convert tmp/10jj1p1258750534.ps tmp/10jj1p1258750534.png")
> 
> 
> proc.time()
   user  system elapsed 
  2.703   1.566   3.220