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(103.2
+ ,123297
+ ,116476
+ ,109375
+ ,106370
+ ,103.7
+ ,114813
+ ,123297
+ ,116476
+ ,109375
+ ,106.2
+ ,117925
+ ,114813
+ ,123297
+ ,116476
+ ,107.7
+ ,126466
+ ,117925
+ ,114813
+ ,123297
+ ,109.9
+ ,131235
+ ,126466
+ ,117925
+ ,114813
+ ,111.7
+ ,120546
+ ,131235
+ ,126466
+ ,117925
+ ,114.9
+ ,123791
+ ,120546
+ ,131235
+ ,126466
+ ,116
+ ,129813
+ ,123791
+ ,120546
+ ,131235
+ ,118.3
+ ,133463
+ ,129813
+ ,123791
+ ,120546
+ ,120.4
+ ,122987
+ ,133463
+ ,129813
+ ,123791
+ ,126
+ ,125418
+ ,122987
+ ,133463
+ ,129813
+ ,128.1
+ ,130199
+ ,125418
+ ,122987
+ ,133463
+ ,130.1
+ ,133016
+ ,130199
+ ,125418
+ ,122987
+ ,130.8
+ ,121454
+ ,133016
+ ,130199
+ ,125418
+ ,133.6
+ ,122044
+ ,121454
+ ,133016
+ ,130199
+ ,134.2
+ ,128313
+ ,122044
+ ,121454
+ ,133016
+ ,135.5
+ ,131556
+ ,128313
+ ,122044
+ ,121454
+ ,136.2
+ ,120027
+ ,131556
+ ,128313
+ ,122044
+ ,139.1
+ ,123001
+ ,120027
+ ,131556
+ ,128313
+ ,139
+ ,130111
+ ,123001
+ ,120027
+ ,131556
+ ,139.6
+ ,132524
+ ,130111
+ ,123001
+ ,120027
+ ,138.7
+ ,123742
+ ,132524
+ ,130111
+ ,123001
+ ,140.9
+ ,124931
+ ,123742
+ ,132524
+ ,130111
+ ,141.3
+ ,133646
+ ,124931
+ ,123742
+ ,132524
+ ,141.8
+ ,136557
+ ,133646
+ ,124931
+ ,123742
+ ,142
+ ,127509
+ ,136557
+ ,133646
+ ,124931
+ ,144.5
+ ,128945
+ ,127509
+ ,136557
+ ,133646
+ ,144.6
+ ,137191
+ ,128945
+ ,127509
+ ,136557
+ ,145.5
+ ,139716
+ ,137191
+ ,128945
+ ,127509
+ ,146.8
+ ,129083
+ ,139716
+ ,137191
+ ,128945
+ ,149.5
+ ,131604
+ ,129083
+ ,139716
+ ,137191
+ ,149.9
+ ,139413
+ ,131604
+ ,129083
+ ,139716
+ ,150.1
+ ,143125
+ ,139413
+ ,131604
+ ,129083
+ ,150.9
+ ,133948
+ ,143125
+ ,139413
+ ,131604
+ ,152.8
+ ,137116
+ ,133948
+ ,143125
+ ,139413
+ ,153.1
+ ,144864
+ ,137116
+ ,133948
+ ,143125
+ ,154
+ ,149277
+ ,144864
+ ,137116
+ ,133948
+ ,154.9
+ ,138796
+ ,149277
+ ,144864
+ ,137116
+ ,156.9
+ ,143258
+ ,138796
+ ,149277
+ ,144864
+ ,158.4
+ ,150034
+ ,143258
+ ,138796
+ ,149277
+ ,159.7
+ ,154708
+ ,150034
+ ,143258
+ ,138796
+ ,160.2
+ ,144888
+ ,154708
+ ,150034
+ ,143258
+ ,163.2
+ ,148762
+ ,144888
+ ,154708
+ ,150034
+ ,163.7
+ ,156500
+ ,148762
+ ,144888
+ ,154708
+ ,164.4
+ ,161088
+ ,156500
+ ,148762
+ ,144888
+ ,163.7
+ ,152772
+ ,161088
+ ,156500
+ ,148762
+ ,165.5
+ ,158011
+ ,152772
+ ,161088
+ ,156500
+ ,165.6
+ ,163318
+ ,158011
+ ,152772
+ ,161088
+ ,166.8
+ ,169969
+ ,163318
+ ,158011
+ ,152772
+ ,167.5
+ ,162269
+ ,169969
+ ,163318
+ ,158011
+ ,170.6
+ ,165765
+ ,162269
+ ,169969
+ ,163318
+ ,170.9
+ ,170600
+ ,165765
+ ,162269
+ ,169969
+ ,172
+ ,174681
+ ,170600
+ ,165765
+ ,162269
+ ,171.8
+ ,166364
+ ,174681
+ ,170600
+ ,165765
+ ,173.9
+ ,170240
+ ,166364
+ ,174681
+ ,170600
+ ,174
+ ,176150
+ ,170240
+ ,166364
+ ,174681
+ ,173.8
+ ,182056
+ ,176150
+ ,170240
+ ,166364
+ ,173.9
+ ,172218
+ ,182056
+ ,176150
+ ,170240
+ ,176
+ ,177856
+ ,172218
+ ,182056
+ ,176150
+ ,176.6
+ ,182253
+ ,177856
+ ,172218
+ ,182056
+ ,178.2
+ ,188090
+ ,182253
+ ,177856
+ ,172218
+ ,179.2
+ ,176863
+ ,188090
+ ,182253
+ ,177856
+ ,181.3
+ ,183273
+ ,176863
+ ,188090
+ ,182253
+ ,181.8
+ ,187969
+ ,183273
+ ,176863
+ ,188090
+ ,182.9
+ ,194650
+ ,187969
+ ,183273
+ ,176863
+ ,183.8
+ ,183036
+ ,194650
+ ,187969
+ ,183273
+ ,186.3
+ ,189516
+ ,183036
+ ,194650
+ ,187969
+ ,187.4
+ ,193805
+ ,189516
+ ,183036
+ ,194650
+ ,189.2
+ ,200499
+ ,193805
+ ,189516
+ ,183036
+ ,189.7
+ ,188142
+ ,200499
+ ,193805
+ ,189516
+ ,191.9
+ ,193732
+ ,188142
+ ,200499
+ ,193805
+ ,192.6
+ ,197126
+ ,193732
+ ,188142
+ ,200499
+ ,193.7
+ ,205140
+ ,197126
+ ,193732
+ ,188142
+ ,194.2
+ ,191751
+ ,205140
+ ,197126
+ ,193732
+ ,197.6
+ ,196700
+ ,191751
+ ,205140
+ ,197126
+ ,199.3
+ ,199784
+ ,196700
+ ,191751
+ ,205140
+ ,201.4
+ ,207360
+ ,199784
+ ,196700
+ ,191751
+ ,203
+ ,196101
+ ,207360
+ ,199784
+ ,196700
+ ,206.3
+ ,200824
+ ,196101
+ ,207360
+ ,199784
+ ,207.1
+ ,205743
+ ,200824
+ ,196101
+ ,207360
+ ,209.8
+ ,212489
+ ,205743
+ ,200824
+ ,196101
+ ,211.1
+ ,200810
+ ,212489
+ ,205743
+ ,200824
+ ,215.3
+ ,203683
+ ,200810
+ ,212489
+ ,205743
+ ,217.4
+ ,207286
+ ,203683
+ ,200810
+ ,212489
+ ,215.5
+ ,210910
+ ,207286
+ ,203683
+ ,200810
+ ,210.9
+ ,194915
+ ,210910
+ ,207286
+ ,203683
+ ,212.6
+ ,217920
+ ,194915
+ ,210910
+ ,207286)
+ ,dim=c(5
+ ,87)
+ ,dimnames=list(c('RPI'
+ ,'HFCE'
+ ,'HFCE-1'
+ ,'HFCE-2'
+ ,'HFCE-3')
+ ,1:87))
> y <- array(NA,dim=c(5,87),dimnames=list(c('RPI','HFCE','HFCE-1','HFCE-2','HFCE-3'),1:87))
> 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 Quarterly Dummies'
> par1 = '2'
> #'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
HFCE RPI HFCE-1 HFCE-2 HFCE-3 Q1 Q2 Q3 t
1 123297 103.2 116476 109375 106370 1 0 0 1
2 114813 103.7 123297 116476 109375 0 1 0 2
3 117925 106.2 114813 123297 116476 0 0 1 3
4 126466 107.7 117925 114813 123297 0 0 0 4
5 131235 109.9 126466 117925 114813 1 0 0 5
6 120546 111.7 131235 126466 117925 0 1 0 6
7 123791 114.9 120546 131235 126466 0 0 1 7
8 129813 116.0 123791 120546 131235 0 0 0 8
9 133463 118.3 129813 123791 120546 1 0 0 9
10 122987 120.4 133463 129813 123791 0 1 0 10
11 125418 126.0 122987 133463 129813 0 0 1 11
12 130199 128.1 125418 122987 133463 0 0 0 12
13 133016 130.1 130199 125418 122987 1 0 0 13
14 121454 130.8 133016 130199 125418 0 1 0 14
15 122044 133.6 121454 133016 130199 0 0 1 15
16 128313 134.2 122044 121454 133016 0 0 0 16
17 131556 135.5 128313 122044 121454 1 0 0 17
18 120027 136.2 131556 128313 122044 0 1 0 18
19 123001 139.1 120027 131556 128313 0 0 1 19
20 130111 139.0 123001 120027 131556 0 0 0 20
21 132524 139.6 130111 123001 120027 1 0 0 21
22 123742 138.7 132524 130111 123001 0 1 0 22
23 124931 140.9 123742 132524 130111 0 0 1 23
24 133646 141.3 124931 123742 132524 0 0 0 24
25 136557 141.8 133646 124931 123742 1 0 0 25
26 127509 142.0 136557 133646 124931 0 1 0 26
27 128945 144.5 127509 136557 133646 0 0 1 27
28 137191 144.6 128945 127509 136557 0 0 0 28
29 139716 145.5 137191 128945 127509 1 0 0 29
30 129083 146.8 139716 137191 128945 0 1 0 30
31 131604 149.5 129083 139716 137191 0 0 1 31
32 139413 149.9 131604 129083 139716 0 0 0 32
33 143125 150.1 139413 131604 129083 1 0 0 33
34 133948 150.9 143125 139413 131604 0 1 0 34
35 137116 152.8 133948 143125 139413 0 0 1 35
36 144864 153.1 137116 133948 143125 0 0 0 36
37 149277 154.0 144864 137116 133948 1 0 0 37
38 138796 154.9 149277 144864 137116 0 1 0 38
39 143258 156.9 138796 149277 144864 0 0 1 39
40 150034 158.4 143258 138796 149277 0 0 0 40
41 154708 159.7 150034 143258 138796 1 0 0 41
42 144888 160.2 154708 150034 143258 0 1 0 42
43 148762 163.2 144888 154708 150034 0 0 1 43
44 156500 163.7 148762 144888 154708 0 0 0 44
45 161088 164.4 156500 148762 144888 1 0 0 45
46 152772 163.7 161088 156500 148762 0 1 0 46
47 158011 165.5 152772 161088 156500 0 0 1 47
48 163318 165.6 158011 152772 161088 0 0 0 48
49 169969 166.8 163318 158011 152772 1 0 0 49
50 162269 167.5 169969 163318 158011 0 1 0 50
51 165765 170.6 162269 169969 163318 0 0 1 51
52 170600 170.9 165765 162269 169969 0 0 0 52
53 174681 172.0 170600 165765 162269 1 0 0 53
54 166364 171.8 174681 170600 165765 0 1 0 54
55 170240 173.9 166364 174681 170600 0 0 1 55
56 176150 174.0 170240 166364 174681 0 0 0 56
57 182056 173.8 176150 170240 166364 1 0 0 57
58 172218 173.9 182056 176150 170240 0 1 0 58
59 177856 176.0 172218 182056 176150 0 0 1 59
60 182253 176.6 177856 172218 182056 0 0 0 60
61 188090 178.2 182253 177856 172218 1 0 0 61
62 176863 179.2 188090 182253 177856 0 1 0 62
63 183273 181.3 176863 188090 182253 0 0 1 63
64 187969 181.8 183273 176863 188090 0 0 0 64
65 194650 182.9 187969 183273 176863 1 0 0 65
66 183036 183.8 194650 187969 183273 0 1 0 66
67 189516 186.3 183036 194650 187969 0 0 1 67
68 193805 187.4 189516 183036 194650 0 0 0 68
69 200499 189.2 193805 189516 183036 1 0 0 69
70 188142 189.7 200499 193805 189516 0 1 0 70
71 193732 191.9 188142 200499 193805 0 0 1 71
72 197126 192.6 193732 188142 200499 0 0 0 72
73 205140 193.7 197126 193732 188142 1 0 0 73
74 191751 194.2 205140 197126 193732 0 1 0 74
75 196700 197.6 191751 205140 197126 0 0 1 75
76 199784 199.3 196700 191751 205140 0 0 0 76
77 207360 201.4 199784 196700 191751 1 0 0 77
78 196101 203.0 207360 199784 196700 0 1 0 78
79 200824 206.3 196101 207360 199784 0 0 1 79
80 205743 207.1 200824 196101 207360 0 0 0 80
81 212489 209.8 205743 200824 196101 1 0 0 81
82 200810 211.1 212489 205743 200824 0 1 0 82
83 203683 215.3 200810 212489 205743 0 0 1 83
84 207286 217.4 203683 200810 212489 0 0 0 84
85 210910 215.5 207286 203683 200810 1 0 0 85
86 194915 210.9 210910 207286 203683 0 1 0 86
87 217920 212.6 194915 210910 207286 0 0 1 87
> k <- length(x[1,])
> df <- as.data.frame(x)
> (mylm <- lm(df))
Call:
lm(formula = df)
Coefficients:
(Intercept) RPI `HFCE-1` `HFCE-2` `HFCE-3` Q1
7.148e+04 -4.217e+02 -1.049e-01 8.700e-01 7.424e-02 2.800e+03
Q2 Q3 t
-1.322e+04 -1.400e+04 6.958e+02
> (mysum <- summary(mylm))
Call:
lm(formula = df)
Residuals:
Min 1Q Median 3Q Max
-7599.77 -693.23 -35.29 628.55 11109.84
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.148e+04 1.108e+04 6.449 8.61e-09 ***
RPI -4.217e+02 7.935e+01 -5.314 9.89e-07 ***
`HFCE-1` -1.049e-01 1.584e-01 -0.662 0.510
`HFCE-2` 8.700e-01 1.556e-01 5.592 3.20e-07 ***
`HFCE-3` 7.424e-02 1.618e-01 0.459 0.648
Q1 2.800e+03 2.549e+03 1.099 0.275
Q2 -1.322e+04 2.825e+03 -4.681 1.18e-05 ***
Q3 -1.400e+04 2.381e+03 -5.882 9.59e-08 ***
t 6.958e+02 1.243e+02 5.596 3.14e-07 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 2004 on 78 degrees of freedom
Multiple R-squared: 0.9961, Adjusted R-squared: 0.9957
F-statistic: 2466 on 8 and 78 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,] 7.251928e-02 1.450386e-01 0.9274807
[2,] 2.269881e-02 4.539762e-02 0.9773012
[3,] 6.533617e-03 1.306723e-02 0.9934664
[4,] 1.825253e-03 3.650506e-03 0.9981747
[5,] 1.994634e-03 3.989269e-03 0.9980054
[6,] 8.689201e-04 1.737840e-03 0.9991311
[7,] 2.744541e-04 5.489082e-04 0.9997255
[8,] 2.722542e-04 5.445084e-04 0.9997277
[9,] 1.873613e-04 3.747226e-04 0.9998126
[10,] 9.980750e-05 1.996150e-04 0.9999002
[11,] 7.557508e-05 1.511502e-04 0.9999244
[12,] 7.067704e-05 1.413541e-04 0.9999293
[13,] 4.747544e-05 9.495088e-05 0.9999525
[14,] 2.005423e-05 4.010845e-05 0.9999799
[15,] 7.397436e-06 1.479487e-05 0.9999926
[16,] 3.309133e-06 6.618267e-06 0.9999967
[17,] 1.176294e-06 2.352587e-06 0.9999988
[18,] 5.218922e-07 1.043784e-06 0.9999995
[19,] 4.831923e-07 9.663845e-07 0.9999995
[20,] 1.647896e-07 3.295793e-07 0.9999998
[21,] 7.311984e-08 1.462397e-07 0.9999999
[22,] 2.493291e-08 4.986583e-08 1.0000000
[23,] 1.421472e-08 2.842944e-08 1.0000000
[24,] 1.218743e-08 2.437486e-08 1.0000000
[25,] 4.606489e-09 9.212979e-09 1.0000000
[26,] 2.358403e-09 4.716805e-09 1.0000000
[27,] 8.178078e-10 1.635616e-09 1.0000000
[28,] 1.650296e-09 3.300592e-09 1.0000000
[29,] 6.139971e-10 1.227994e-09 1.0000000
[30,] 3.058713e-10 6.117426e-10 1.0000000
[31,] 1.497578e-10 2.995157e-10 1.0000000
[32,] 3.191874e-10 6.383748e-10 1.0000000
[33,] 1.934289e-10 3.868577e-10 1.0000000
[34,] 9.243104e-11 1.848621e-10 1.0000000
[35,] 8.473274e-11 1.694655e-10 1.0000000
[36,] 3.793501e-10 7.587001e-10 1.0000000
[37,] 1.008831e-09 2.017661e-09 1.0000000
[38,] 3.864115e-10 7.728230e-10 1.0000000
[39,] 2.159279e-09 4.318558e-09 1.0000000
[40,] 9.832501e-10 1.966500e-09 1.0000000
[41,] 7.591562e-09 1.518312e-08 1.0000000
[42,] 9.886893e-09 1.977379e-08 1.0000000
[43,] 4.104507e-09 8.209013e-09 1.0000000
[44,] 1.920297e-09 3.840595e-09 1.0000000
[45,] 1.556140e-09 3.112280e-09 1.0000000
[46,] 6.765026e-10 1.353005e-09 1.0000000
[47,] 3.486307e-10 6.972615e-10 1.0000000
[48,] 2.719471e-10 5.438942e-10 1.0000000
[49,] 7.316613e-10 1.463323e-09 1.0000000
[50,] 3.641899e-10 7.283798e-10 1.0000000
[51,] 2.279662e-10 4.559324e-10 1.0000000
[52,] 5.654772e-10 1.130954e-09 1.0000000
[53,] 4.104136e-10 8.208273e-10 1.0000000
[54,] 2.520596e-10 5.041191e-10 1.0000000
[55,] 1.256520e-10 2.513040e-10 1.0000000
[56,] 7.550070e-10 1.510014e-09 1.0000000
[57,] 1.781351e-09 3.562702e-09 1.0000000
[58,] 1.625140e-09 3.250280e-09 1.0000000
[59,] 5.976613e-10 1.195323e-09 1.0000000
[60,] 4.689215e-10 9.378429e-10 1.0000000
[61,] 3.007718e-10 6.015436e-10 1.0000000
[62,] 2.166302e-10 4.332605e-10 1.0000000
[63,] 1.977530e-10 3.955060e-10 1.0000000
[64,] 1.128662e-10 2.257325e-10 1.0000000
> postscript(file="/var/www/html/rcomp/tmp/1yij11258726369.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/2kl4a1258726369.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/3by4t1258726370.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/4lduo1258726370.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/5y4ti1258726370.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 = 87
Frequency = 1
1 2 3 4 5
998.3980900 2365.7261118 -733.8824458 941.7344520 1960.3243804
6 7 8 9 10
194.8689975 -1029.3870912 43.2543883 -231.1536998 407.2833818
11 12 13 14 15
563.9427114 629.6556740 -41.9950287 -27.1145710 -2189.1567038
16 17 18 19 20
-454.3109632 843.2436694 -222.1385956 -435.6495231 2034.6764521
21 22 23 24 25
218.6487707 229.8807283 -1116.0589182 654.5748575 811.7542309
26 27 28 29 30
-190.8367413 -1743.5915878 -348.1014217 -652.5824975 -2427.1020190
31 32 33 34 35
-2606.0778295 0.3390062 -284.4534907 -389.6844576 -1106.6049047
36 37 38 39 40
109.6082749 143.5931888 -1144.9226850 -1267.7756196 700.6125012
41 42 43 44 45
33.3296006 14.1188464 -360.6676720 1491.9244823 1049.0199465
46 47 48 49 50
1225.4667718 1870.7893146 -35.2924398 241.4825051 3854.3059737
51 52 53 54 55
1754.9974128 -1410.6095060 -2324.7654203 562.1091225 627.4425141
56 57 58 59 60
-780.0559571 -589.4689809 130.7867599 131.0202745 -1206.0455616
61 62 63 64 65
-1904.1776211 -1015.1027216 -216.1682204 -1.8566571 -604.0195841
66 67 68 69 70
-373.0591749 -132.5375517 209.1635576 -159.6787692 -490.1192454
71 72 73 74 75
-1325.1069035 -1494.9344333 -103.3144649 -482.5926728 -2642.7847972
76 77 78 79 80
-1968.2772948 8.8098956 2494.7379814 693.8172588 979.4745721
81 82 83 84 85
2610.4869077 2883.1533583 153.5967710 -95.5339837 -2023.4816286
86 87
-7599.7651493 11109.8435115
> postscript(file="/var/www/html/rcomp/tmp/612671258726370.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 = 87
Frequency = 1
lag(myerror, k = 1) myerror
0 998.3980900 NA
1 2365.7261118 998.3980900
2 -733.8824458 2365.7261118
3 941.7344520 -733.8824458
4 1960.3243804 941.7344520
5 194.8689975 1960.3243804
6 -1029.3870912 194.8689975
7 43.2543883 -1029.3870912
8 -231.1536998 43.2543883
9 407.2833818 -231.1536998
10 563.9427114 407.2833818
11 629.6556740 563.9427114
12 -41.9950287 629.6556740
13 -27.1145710 -41.9950287
14 -2189.1567038 -27.1145710
15 -454.3109632 -2189.1567038
16 843.2436694 -454.3109632
17 -222.1385956 843.2436694
18 -435.6495231 -222.1385956
19 2034.6764521 -435.6495231
20 218.6487707 2034.6764521
21 229.8807283 218.6487707
22 -1116.0589182 229.8807283
23 654.5748575 -1116.0589182
24 811.7542309 654.5748575
25 -190.8367413 811.7542309
26 -1743.5915878 -190.8367413
27 -348.1014217 -1743.5915878
28 -652.5824975 -348.1014217
29 -2427.1020190 -652.5824975
30 -2606.0778295 -2427.1020190
31 0.3390062 -2606.0778295
32 -284.4534907 0.3390062
33 -389.6844576 -284.4534907
34 -1106.6049047 -389.6844576
35 109.6082749 -1106.6049047
36 143.5931888 109.6082749
37 -1144.9226850 143.5931888
38 -1267.7756196 -1144.9226850
39 700.6125012 -1267.7756196
40 33.3296006 700.6125012
41 14.1188464 33.3296006
42 -360.6676720 14.1188464
43 1491.9244823 -360.6676720
44 1049.0199465 1491.9244823
45 1225.4667718 1049.0199465
46 1870.7893146 1225.4667718
47 -35.2924398 1870.7893146
48 241.4825051 -35.2924398
49 3854.3059737 241.4825051
50 1754.9974128 3854.3059737
51 -1410.6095060 1754.9974128
52 -2324.7654203 -1410.6095060
53 562.1091225 -2324.7654203
54 627.4425141 562.1091225
55 -780.0559571 627.4425141
56 -589.4689809 -780.0559571
57 130.7867599 -589.4689809
58 131.0202745 130.7867599
59 -1206.0455616 131.0202745
60 -1904.1776211 -1206.0455616
61 -1015.1027216 -1904.1776211
62 -216.1682204 -1015.1027216
63 -1.8566571 -216.1682204
64 -604.0195841 -1.8566571
65 -373.0591749 -604.0195841
66 -132.5375517 -373.0591749
67 209.1635576 -132.5375517
68 -159.6787692 209.1635576
69 -490.1192454 -159.6787692
70 -1325.1069035 -490.1192454
71 -1494.9344333 -1325.1069035
72 -103.3144649 -1494.9344333
73 -482.5926728 -103.3144649
74 -2642.7847972 -482.5926728
75 -1968.2772948 -2642.7847972
76 8.8098956 -1968.2772948
77 2494.7379814 8.8098956
78 693.8172588 2494.7379814
79 979.4745721 693.8172588
80 2610.4869077 979.4745721
81 2883.1533583 2610.4869077
82 153.5967710 2883.1533583
83 -95.5339837 153.5967710
84 -2023.4816286 -95.5339837
85 -7599.7651493 -2023.4816286
86 11109.8435115 -7599.7651493
87 NA 11109.8435115
> dum1 <- dum[2:length(myerror),]
> dum1
lag(myerror, k = 1) myerror
[1,] 2365.7261118 998.3980900
[2,] -733.8824458 2365.7261118
[3,] 941.7344520 -733.8824458
[4,] 1960.3243804 941.7344520
[5,] 194.8689975 1960.3243804
[6,] -1029.3870912 194.8689975
[7,] 43.2543883 -1029.3870912
[8,] -231.1536998 43.2543883
[9,] 407.2833818 -231.1536998
[10,] 563.9427114 407.2833818
[11,] 629.6556740 563.9427114
[12,] -41.9950287 629.6556740
[13,] -27.1145710 -41.9950287
[14,] -2189.1567038 -27.1145710
[15,] -454.3109632 -2189.1567038
[16,] 843.2436694 -454.3109632
[17,] -222.1385956 843.2436694
[18,] -435.6495231 -222.1385956
[19,] 2034.6764521 -435.6495231
[20,] 218.6487707 2034.6764521
[21,] 229.8807283 218.6487707
[22,] -1116.0589182 229.8807283
[23,] 654.5748575 -1116.0589182
[24,] 811.7542309 654.5748575
[25,] -190.8367413 811.7542309
[26,] -1743.5915878 -190.8367413
[27,] -348.1014217 -1743.5915878
[28,] -652.5824975 -348.1014217
[29,] -2427.1020190 -652.5824975
[30,] -2606.0778295 -2427.1020190
[31,] 0.3390062 -2606.0778295
[32,] -284.4534907 0.3390062
[33,] -389.6844576 -284.4534907
[34,] -1106.6049047 -389.6844576
[35,] 109.6082749 -1106.6049047
[36,] 143.5931888 109.6082749
[37,] -1144.9226850 143.5931888
[38,] -1267.7756196 -1144.9226850
[39,] 700.6125012 -1267.7756196
[40,] 33.3296006 700.6125012
[41,] 14.1188464 33.3296006
[42,] -360.6676720 14.1188464
[43,] 1491.9244823 -360.6676720
[44,] 1049.0199465 1491.9244823
[45,] 1225.4667718 1049.0199465
[46,] 1870.7893146 1225.4667718
[47,] -35.2924398 1870.7893146
[48,] 241.4825051 -35.2924398
[49,] 3854.3059737 241.4825051
[50,] 1754.9974128 3854.3059737
[51,] -1410.6095060 1754.9974128
[52,] -2324.7654203 -1410.6095060
[53,] 562.1091225 -2324.7654203
[54,] 627.4425141 562.1091225
[55,] -780.0559571 627.4425141
[56,] -589.4689809 -780.0559571
[57,] 130.7867599 -589.4689809
[58,] 131.0202745 130.7867599
[59,] -1206.0455616 131.0202745
[60,] -1904.1776211 -1206.0455616
[61,] -1015.1027216 -1904.1776211
[62,] -216.1682204 -1015.1027216
[63,] -1.8566571 -216.1682204
[64,] -604.0195841 -1.8566571
[65,] -373.0591749 -604.0195841
[66,] -132.5375517 -373.0591749
[67,] 209.1635576 -132.5375517
[68,] -159.6787692 209.1635576
[69,] -490.1192454 -159.6787692
[70,] -1325.1069035 -490.1192454
[71,] -1494.9344333 -1325.1069035
[72,] -103.3144649 -1494.9344333
[73,] -482.5926728 -103.3144649
[74,] -2642.7847972 -482.5926728
[75,] -1968.2772948 -2642.7847972
[76,] 8.8098956 -1968.2772948
[77,] 2494.7379814 8.8098956
[78,] 693.8172588 2494.7379814
[79,] 979.4745721 693.8172588
[80,] 2610.4869077 979.4745721
[81,] 2883.1533583 2610.4869077
[82,] 153.5967710 2883.1533583
[83,] -95.5339837 153.5967710
[84,] -2023.4816286 -95.5339837
[85,] -7599.7651493 -2023.4816286
[86,] 11109.8435115 -7599.7651493
> z <- as.data.frame(dum1)
> z
lag(myerror, k = 1) myerror
1 2365.7261118 998.3980900
2 -733.8824458 2365.7261118
3 941.7344520 -733.8824458
4 1960.3243804 941.7344520
5 194.8689975 1960.3243804
6 -1029.3870912 194.8689975
7 43.2543883 -1029.3870912
8 -231.1536998 43.2543883
9 407.2833818 -231.1536998
10 563.9427114 407.2833818
11 629.6556740 563.9427114
12 -41.9950287 629.6556740
13 -27.1145710 -41.9950287
14 -2189.1567038 -27.1145710
15 -454.3109632 -2189.1567038
16 843.2436694 -454.3109632
17 -222.1385956 843.2436694
18 -435.6495231 -222.1385956
19 2034.6764521 -435.6495231
20 218.6487707 2034.6764521
21 229.8807283 218.6487707
22 -1116.0589182 229.8807283
23 654.5748575 -1116.0589182
24 811.7542309 654.5748575
25 -190.8367413 811.7542309
26 -1743.5915878 -190.8367413
27 -348.1014217 -1743.5915878
28 -652.5824975 -348.1014217
29 -2427.1020190 -652.5824975
30 -2606.0778295 -2427.1020190
31 0.3390062 -2606.0778295
32 -284.4534907 0.3390062
33 -389.6844576 -284.4534907
34 -1106.6049047 -389.6844576
35 109.6082749 -1106.6049047
36 143.5931888 109.6082749
37 -1144.9226850 143.5931888
38 -1267.7756196 -1144.9226850
39 700.6125012 -1267.7756196
40 33.3296006 700.6125012
41 14.1188464 33.3296006
42 -360.6676720 14.1188464
43 1491.9244823 -360.6676720
44 1049.0199465 1491.9244823
45 1225.4667718 1049.0199465
46 1870.7893146 1225.4667718
47 -35.2924398 1870.7893146
48 241.4825051 -35.2924398
49 3854.3059737 241.4825051
50 1754.9974128 3854.3059737
51 -1410.6095060 1754.9974128
52 -2324.7654203 -1410.6095060
53 562.1091225 -2324.7654203
54 627.4425141 562.1091225
55 -780.0559571 627.4425141
56 -589.4689809 -780.0559571
57 130.7867599 -589.4689809
58 131.0202745 130.7867599
59 -1206.0455616 131.0202745
60 -1904.1776211 -1206.0455616
61 -1015.1027216 -1904.1776211
62 -216.1682204 -1015.1027216
63 -1.8566571 -216.1682204
64 -604.0195841 -1.8566571
65 -373.0591749 -604.0195841
66 -132.5375517 -373.0591749
67 209.1635576 -132.5375517
68 -159.6787692 209.1635576
69 -490.1192454 -159.6787692
70 -1325.1069035 -490.1192454
71 -1494.9344333 -1325.1069035
72 -103.3144649 -1494.9344333
73 -482.5926728 -103.3144649
74 -2642.7847972 -482.5926728
75 -1968.2772948 -2642.7847972
76 8.8098956 -1968.2772948
77 2494.7379814 8.8098956
78 693.8172588 2494.7379814
79 979.4745721 693.8172588
80 2610.4869077 979.4745721
81 2883.1533583 2610.4869077
82 153.5967710 2883.1533583
83 -95.5339837 153.5967710
84 -2023.4816286 -95.5339837
85 -7599.7651493 -2023.4816286
86 11109.8435115 -7599.7651493
> 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/7u6re1258726370.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/8dpfb1258726370.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/962yq1258726370.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/10as8m1258726370.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/11jpuy1258726370.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/12bcf31258726370.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/13i63l1258726370.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/14yw0q1258726370.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/15ng831258726370.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/16tnfj1258726370.tab")
+ }
>
> system("convert tmp/1yij11258726369.ps tmp/1yij11258726369.png")
> system("convert tmp/2kl4a1258726369.ps tmp/2kl4a1258726369.png")
> system("convert tmp/3by4t1258726370.ps tmp/3by4t1258726370.png")
> system("convert tmp/4lduo1258726370.ps tmp/4lduo1258726370.png")
> system("convert tmp/5y4ti1258726370.ps tmp/5y4ti1258726370.png")
> system("convert tmp/612671258726370.ps tmp/612671258726370.png")
> system("convert tmp/7u6re1258726370.ps tmp/7u6re1258726370.png")
> system("convert tmp/8dpfb1258726370.ps tmp/8dpfb1258726370.png")
> system("convert tmp/962yq1258726370.ps tmp/962yq1258726370.png")
> system("convert tmp/10as8m1258726370.ps tmp/10as8m1258726370.png")
>
>
> proc.time()
user system elapsed
2.885 1.639 3.550