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(100,0,95.84395716,0,105.5073942,1,118.1540031,1,101.8612953,1,109.8419174,1,105.6348802,1,112.927078,1,133.0698623,1,125.6756757,1,146.736359,1,142.5803162,1,106.1448241,1,126.5170831,1,132.7893932,1,121.2391637,1,114.5079041,1,146.1499235,1,146.1244263,1,128.5058644,1,155.5838858,1,125.0382458,1,136.8944416,1,142.2233554,1,117.7715451,1,120.627231,1,127.7664457,1,135.1096379,1,105.7113717,1,117.9245283,1,120.754717,1,107.572667,1,130.4436512,1,107.2157063,1,105.0739419,1,130.1121877,1,109.6379398,1,116.7261601,1,97.11881693,0,140.8975013,1,108.2865885,1,97.65425803,0,112.0346762,1,123.0494646,1,112.4171341,1,116.4966854,1,104.6914839,1,122.2335543,1,99.79602244,0,96.71086181,0,112.3151453,1,102.5497195,1,104.5385008,1,122.0805711,1,80.64762876,0,91.40744518,0,99.51555329,0,106.527282,1,98.49566548,0,106.7567568,1),dim=c(2,60),dimnames=list(c('Y','X'),1:60))
> y <- array(NA,dim=c(2,60),dimnames=list(c('Y','X'),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 = 'Linear Trend'
> par2 = 'Include Monthly 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
Y X M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 t
1 100.00000 0 1 0 0 0 0 0 0 0 0 0 0 1
2 95.84396 0 0 1 0 0 0 0 0 0 0 0 0 2
3 105.50739 1 0 0 1 0 0 0 0 0 0 0 0 3
4 118.15400 1 0 0 0 1 0 0 0 0 0 0 0 4
5 101.86130 1 0 0 0 0 1 0 0 0 0 0 0 5
6 109.84192 1 0 0 0 0 0 1 0 0 0 0 0 6
7 105.63488 1 0 0 0 0 0 0 1 0 0 0 0 7
8 112.92708 1 0 0 0 0 0 0 0 1 0 0 0 8
9 133.06986 1 0 0 0 0 0 0 0 0 1 0 0 9
10 125.67568 1 0 0 0 0 0 0 0 0 0 1 0 10
11 146.73636 1 0 0 0 0 0 0 0 0 0 0 1 11
12 142.58032 1 0 0 0 0 0 0 0 0 0 0 0 12
13 106.14482 1 1 0 0 0 0 0 0 0 0 0 0 13
14 126.51708 1 0 1 0 0 0 0 0 0 0 0 0 14
15 132.78939 1 0 0 1 0 0 0 0 0 0 0 0 15
16 121.23916 1 0 0 0 1 0 0 0 0 0 0 0 16
17 114.50790 1 0 0 0 0 1 0 0 0 0 0 0 17
18 146.14992 1 0 0 0 0 0 1 0 0 0 0 0 18
19 146.12443 1 0 0 0 0 0 0 1 0 0 0 0 19
20 128.50586 1 0 0 0 0 0 0 0 1 0 0 0 20
21 155.58389 1 0 0 0 0 0 0 0 0 1 0 0 21
22 125.03825 1 0 0 0 0 0 0 0 0 0 1 0 22
23 136.89444 1 0 0 0 0 0 0 0 0 0 0 1 23
24 142.22336 1 0 0 0 0 0 0 0 0 0 0 0 24
25 117.77155 1 1 0 0 0 0 0 0 0 0 0 0 25
26 120.62723 1 0 1 0 0 0 0 0 0 0 0 0 26
27 127.76645 1 0 0 1 0 0 0 0 0 0 0 0 27
28 135.10964 1 0 0 0 1 0 0 0 0 0 0 0 28
29 105.71137 1 0 0 0 0 1 0 0 0 0 0 0 29
30 117.92453 1 0 0 0 0 0 1 0 0 0 0 0 30
31 120.75472 1 0 0 0 0 0 0 1 0 0 0 0 31
32 107.57267 1 0 0 0 0 0 0 0 1 0 0 0 32
33 130.44365 1 0 0 0 0 0 0 0 0 1 0 0 33
34 107.21571 1 0 0 0 0 0 0 0 0 0 1 0 34
35 105.07394 1 0 0 0 0 0 0 0 0 0 0 1 35
36 130.11219 1 0 0 0 0 0 0 0 0 0 0 0 36
37 109.63794 1 1 0 0 0 0 0 0 0 0 0 0 37
38 116.72616 1 0 1 0 0 0 0 0 0 0 0 0 38
39 97.11882 0 0 0 1 0 0 0 0 0 0 0 0 39
40 140.89750 1 0 0 0 1 0 0 0 0 0 0 0 40
41 108.28659 1 0 0 0 0 1 0 0 0 0 0 0 41
42 97.65426 0 0 0 0 0 0 1 0 0 0 0 0 42
43 112.03468 1 0 0 0 0 0 0 1 0 0 0 0 43
44 123.04946 1 0 0 0 0 0 0 0 1 0 0 0 44
45 112.41713 1 0 0 0 0 0 0 0 0 1 0 0 45
46 116.49669 1 0 0 0 0 0 0 0 0 0 1 0 46
47 104.69148 1 0 0 0 0 0 0 0 0 0 0 1 47
48 122.23355 1 0 0 0 0 0 0 0 0 0 0 0 48
49 99.79602 0 1 0 0 0 0 0 0 0 0 0 0 49
50 96.71086 0 0 1 0 0 0 0 0 0 0 0 0 50
51 112.31515 1 0 0 1 0 0 0 0 0 0 0 0 51
52 102.54972 1 0 0 0 1 0 0 0 0 0 0 0 52
53 104.53850 1 0 0 0 0 1 0 0 0 0 0 0 53
54 122.08057 1 0 0 0 0 0 1 0 0 0 0 0 54
55 80.64763 0 0 0 0 0 0 0 1 0 0 0 0 55
56 91.40745 0 0 0 0 0 0 0 0 1 0 0 0 56
57 99.51555 0 0 0 0 0 0 0 0 0 1 0 0 57
58 106.52728 1 0 0 0 0 0 0 0 0 0 1 0 58
59 98.49567 0 0 0 0 0 0 0 0 0 0 0 1 59
60 106.75676 1 0 0 0 0 0 0 0 0 0 0 0 60
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) X M1 M2 M3 M4
116.5972 21.2164 -16.3845 -11.5186 -11.6966 -7.1984
M5 M6 M7 M8 M9 M10
-23.5564 -7.3131 -12.7532 -12.8490 0.9154 -13.0923
M11 t
-6.4105 -0.2509
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-19.8569 -8.7934 0.4469 5.3056 25.8311
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 116.59717 8.69287 13.413 < 2e-16 ***
X 21.21640 4.77563 4.443 5.55e-05 ***
M1 -16.38449 8.01708 -2.044 0.0467 *
M2 -11.51860 7.99670 -1.440 0.1565
M3 -11.69660 7.76743 -1.506 0.1389
M4 -7.19841 7.66088 -0.940 0.3523
M5 -23.55639 7.65157 -3.079 0.0035 **
M6 -7.31310 7.72757 -0.946 0.3489
M7 -12.75318 7.71670 -1.653 0.1052
M8 -12.84904 7.70705 -1.667 0.1023
M9 0.91537 7.69862 0.119 0.9059
M10 -13.09231 7.62358 -1.717 0.0926 .
M11 -6.41047 7.68544 -0.834 0.4085
t -0.25090 0.09748 -2.574 0.0133 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 12.05 on 46 degrees of freedom
Multiple R-squared: 0.5553, Adjusted R-squared: 0.4296
F-statistic: 4.418 on 13 and 46 DF, p-value: 8.439e-05
> 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.6045529 0.7908943 0.3954471
[2,] 0.6820955 0.6358090 0.3179045
[3,] 0.7672622 0.4654757 0.2327378
[4,] 0.6776249 0.6447501 0.3223751
[5,] 0.6637239 0.6725522 0.3362761
[6,] 0.7386174 0.5227653 0.2613826
[7,] 0.8888238 0.2223523 0.1111762
[8,] 0.8948577 0.2102846 0.1051423
[9,] 0.8522701 0.2954598 0.1477299
[10,] 0.8132450 0.3735101 0.1867550
[11,] 0.7756068 0.4487863 0.2243932
[12,] 0.7070740 0.5858520 0.2929260
[13,] 0.7146004 0.5707992 0.2853996
[14,] 0.7518160 0.4963680 0.2481840
[15,] 0.7394620 0.5210760 0.2605380
[16,] 0.7899181 0.4201638 0.2100819
[17,] 0.7906074 0.4187853 0.2093926
[18,] 0.8095862 0.3808275 0.1904138
[19,] 0.8970848 0.2058305 0.1029152
[20,] 0.8518783 0.2962434 0.1481217
[21,] 0.8160656 0.3678688 0.1839344
[22,] 0.7308461 0.5383078 0.2691539
[23,] 0.6292170 0.7415659 0.3707830
[24,] 0.8877880 0.2244240 0.1122120
[25,] 0.7966864 0.4066272 0.2033136
[26,] 0.8067028 0.3865944 0.1932972
[27,] 0.7692247 0.4615507 0.2307753
> postscript(file="/var/www/html/rcomp/tmp/12ao31258719199.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/2osvi1258719199.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/3qrfo1258719199.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/4ztkm1258719199.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/5gei61258719199.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 = 60
Frequency = 1
1 2 3 4 5
0.038219088 -8.732816098 -19.856881970 -11.457558484 -11.141393264
6 7 8 9 10
-19.153159370 -17.669222596 -10.030262940 -3.400992142 3.463400176
11 12 13 14 15
18.093143520 7.777525636 -12.022581670 3.734684984 10.435895272
16 17 18 19 20
-5.361619642 4.515993778 20.165624972 25.831101746 8.559301702
21 22 23 24 25
22.123809600 5.836748518 11.262004362 10.431343078 2.614917572
26 27 28 29 30
0.855611126 8.423726014 11.519632800 -1.269760380 -5.048991986
31 32 33 34 35
3.472170688 -9.363117456 -0.005646758 -8.975012740 -17.547717096
36 37 38 39 40
1.330953620 -2.507909486 -0.034681532 2.003278586 20.318274442
41 42 43 44 45
4.316234662 -1.092080914 -2.237091870 9.124458386 -15.021385616
46 47 48 49 50
3.316744602 -14.919396854 -3.536901538 11.877354496 4.177201520
51 52 53 54 55
-1.006017902 -15.018729116 3.578925204 5.128607298 -9.396957968
56 57 58 59 60
1.709620308 -3.695785084 -3.641880556 3.111966068 -16.002920796
> postscript(file="/var/www/html/rcomp/tmp/68k0z1258719199.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 = 60
Frequency = 1
lag(myerror, k = 1) myerror
0 0.038219088 NA
1 -8.732816098 0.038219088
2 -19.856881970 -8.732816098
3 -11.457558484 -19.856881970
4 -11.141393264 -11.457558484
5 -19.153159370 -11.141393264
6 -17.669222596 -19.153159370
7 -10.030262940 -17.669222596
8 -3.400992142 -10.030262940
9 3.463400176 -3.400992142
10 18.093143520 3.463400176
11 7.777525636 18.093143520
12 -12.022581670 7.777525636
13 3.734684984 -12.022581670
14 10.435895272 3.734684984
15 -5.361619642 10.435895272
16 4.515993778 -5.361619642
17 20.165624972 4.515993778
18 25.831101746 20.165624972
19 8.559301702 25.831101746
20 22.123809600 8.559301702
21 5.836748518 22.123809600
22 11.262004362 5.836748518
23 10.431343078 11.262004362
24 2.614917572 10.431343078
25 0.855611126 2.614917572
26 8.423726014 0.855611126
27 11.519632800 8.423726014
28 -1.269760380 11.519632800
29 -5.048991986 -1.269760380
30 3.472170688 -5.048991986
31 -9.363117456 3.472170688
32 -0.005646758 -9.363117456
33 -8.975012740 -0.005646758
34 -17.547717096 -8.975012740
35 1.330953620 -17.547717096
36 -2.507909486 1.330953620
37 -0.034681532 -2.507909486
38 2.003278586 -0.034681532
39 20.318274442 2.003278586
40 4.316234662 20.318274442
41 -1.092080914 4.316234662
42 -2.237091870 -1.092080914
43 9.124458386 -2.237091870
44 -15.021385616 9.124458386
45 3.316744602 -15.021385616
46 -14.919396854 3.316744602
47 -3.536901538 -14.919396854
48 11.877354496 -3.536901538
49 4.177201520 11.877354496
50 -1.006017902 4.177201520
51 -15.018729116 -1.006017902
52 3.578925204 -15.018729116
53 5.128607298 3.578925204
54 -9.396957968 5.128607298
55 1.709620308 -9.396957968
56 -3.695785084 1.709620308
57 -3.641880556 -3.695785084
58 3.111966068 -3.641880556
59 -16.002920796 3.111966068
60 NA -16.002920796
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] -8.732816098 0.038219088
[2,] -19.856881970 -8.732816098
[3,] -11.457558484 -19.856881970
[4,] -11.141393264 -11.457558484
[5,] -19.153159370 -11.141393264
[6,] -17.669222596 -19.153159370
[7,] -10.030262940 -17.669222596
[8,] -3.400992142 -10.030262940
[9,] 3.463400176 -3.400992142
[10,] 18.093143520 3.463400176
[11,] 7.777525636 18.093143520
[12,] -12.022581670 7.777525636
[13,] 3.734684984 -12.022581670
[14,] 10.435895272 3.734684984
[15,] -5.361619642 10.435895272
[16,] 4.515993778 -5.361619642
[17,] 20.165624972 4.515993778
[18,] 25.831101746 20.165624972
[19,] 8.559301702 25.831101746
[20,] 22.123809600 8.559301702
[21,] 5.836748518 22.123809600
[22,] 11.262004362 5.836748518
[23,] 10.431343078 11.262004362
[24,] 2.614917572 10.431343078
[25,] 0.855611126 2.614917572
[26,] 8.423726014 0.855611126
[27,] 11.519632800 8.423726014
[28,] -1.269760380 11.519632800
[29,] -5.048991986 -1.269760380
[30,] 3.472170688 -5.048991986
[31,] -9.363117456 3.472170688
[32,] -0.005646758 -9.363117456
[33,] -8.975012740 -0.005646758
[34,] -17.547717096 -8.975012740
[35,] 1.330953620 -17.547717096
[36,] -2.507909486 1.330953620
[37,] -0.034681532 -2.507909486
[38,] 2.003278586 -0.034681532
[39,] 20.318274442 2.003278586
[40,] 4.316234662 20.318274442
[41,] -1.092080914 4.316234662
[42,] -2.237091870 -1.092080914
[43,] 9.124458386 -2.237091870
[44,] -15.021385616 9.124458386
[45,] 3.316744602 -15.021385616
[46,] -14.919396854 3.316744602
[47,] -3.536901538 -14.919396854
[48,] 11.877354496 -3.536901538
[49,] 4.177201520 11.877354496
[50,] -1.006017902 4.177201520
[51,] -15.018729116 -1.006017902
[52,] 3.578925204 -15.018729116
[53,] 5.128607298 3.578925204
[54,] -9.396957968 5.128607298
[55,] 1.709620308 -9.396957968
[56,] -3.695785084 1.709620308
[57,] -3.641880556 -3.695785084
[58,] 3.111966068 -3.641880556
[59,] -16.002920796 3.111966068
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 -8.732816098 0.038219088
2 -19.856881970 -8.732816098
3 -11.457558484 -19.856881970
4 -11.141393264 -11.457558484
5 -19.153159370 -11.141393264
6 -17.669222596 -19.153159370
7 -10.030262940 -17.669222596
8 -3.400992142 -10.030262940
9 3.463400176 -3.400992142
10 18.093143520 3.463400176
11 7.777525636 18.093143520
12 -12.022581670 7.777525636
13 3.734684984 -12.022581670
14 10.435895272 3.734684984
15 -5.361619642 10.435895272
16 4.515993778 -5.361619642
17 20.165624972 4.515993778
18 25.831101746 20.165624972
19 8.559301702 25.831101746
20 22.123809600 8.559301702
21 5.836748518 22.123809600
22 11.262004362 5.836748518
23 10.431343078 11.262004362
24 2.614917572 10.431343078
25 0.855611126 2.614917572
26 8.423726014 0.855611126
27 11.519632800 8.423726014
28 -1.269760380 11.519632800
29 -5.048991986 -1.269760380
30 3.472170688 -5.048991986
31 -9.363117456 3.472170688
32 -0.005646758 -9.363117456
33 -8.975012740 -0.005646758
34 -17.547717096 -8.975012740
35 1.330953620 -17.547717096
36 -2.507909486 1.330953620
37 -0.034681532 -2.507909486
38 2.003278586 -0.034681532
39 20.318274442 2.003278586
40 4.316234662 20.318274442
41 -1.092080914 4.316234662
42 -2.237091870 -1.092080914
43 9.124458386 -2.237091870
44 -15.021385616 9.124458386
45 3.316744602 -15.021385616
46 -14.919396854 3.316744602
47 -3.536901538 -14.919396854
48 11.877354496 -3.536901538
49 4.177201520 11.877354496
50 -1.006017902 4.177201520
51 -15.018729116 -1.006017902
52 3.578925204 -15.018729116
53 5.128607298 3.578925204
54 -9.396957968 5.128607298
55 1.709620308 -9.396957968
56 -3.695785084 1.709620308
57 -3.641880556 -3.695785084
58 3.111966068 -3.641880556
59 -16.002920796 3.111966068
> 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/7eum61258719199.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/8mshr1258719199.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/9pfyx1258719199.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/10t3aa1258719199.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/11gz471258719199.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/128qjx1258719199.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/133msf1258719199.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/14l4em1258719199.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/15d6481258719199.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/16ibvn1258719199.tab")
+ }
>
> system("convert tmp/12ao31258719199.ps tmp/12ao31258719199.png")
> system("convert tmp/2osvi1258719199.ps tmp/2osvi1258719199.png")
> system("convert tmp/3qrfo1258719199.ps tmp/3qrfo1258719199.png")
> system("convert tmp/4ztkm1258719199.ps tmp/4ztkm1258719199.png")
> system("convert tmp/5gei61258719199.ps tmp/5gei61258719199.png")
> system("convert tmp/68k0z1258719199.ps tmp/68k0z1258719199.png")
> system("convert tmp/7eum61258719199.ps tmp/7eum61258719199.png")
> system("convert tmp/8mshr1258719199.ps tmp/8mshr1258719199.png")
> system("convert tmp/9pfyx1258719199.ps tmp/9pfyx1258719199.png")
> system("convert tmp/10t3aa1258719199.ps tmp/10t3aa1258719199.png")
>
>
> proc.time()
user system elapsed
2.337 1.516 2.736