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(54.64 + ,4606.07 + ,14.36 + ,52.39 + ,4176.60 + ,14.62 + ,52.51 + ,4331.53 + ,13.51 + ,52.92 + ,3987.30 + ,14.95 + ,55.22 + ,4205.86 + ,16.72 + ,55.41 + ,4331.37 + ,16.33 + ,57.02 + ,4106.30 + ,15.21 + ,58.55 + ,4009.33 + ,16.69 + ,57.49 + ,3857.64 + ,15.81 + ,55.52 + ,3929.00 + ,16.02 + ,57.84 + ,4210.47 + ,16.7 + ,58.69 + ,4445.82 + ,15.99 + ,59.74 + ,4497.07 + ,17.68 + ,60.7 + ,4443.71 + ,18.89 + ,60.74 + ,4529.25 + ,18.72 + ,64.32 + ,4634.22 + ,21.14 + ,66.9 + ,4772.67 + ,20.97 + ,70.93 + ,4881.52 + ,23.75 + ,75.89 + ,5153.13 + ,23.05 + ,80.6 + ,5324.19 + ,23.45 + ,81.39 + ,5209.36 + ,21.74 + ,81.33 + ,5108.92 + ,19.37 + ,77.04 + ,5130.88 + ,21.1 + ,79.54 + ,5195.68 + ,21.2 + ,81.93 + ,5050.42 + ,22.67 + ,80.79 + ,5101.03 + ,22.24 + ,81.98 + ,5139.84 + ,23.78 + ,85.94 + ,5234.31 + ,23.27 + ,86.6 + ,5435.73 + ,25.74 + ,87.42 + ,5633.57 + ,26.1 + ,93.14 + ,5498.28 + ,27.49 + ,95.76 + ,5668.13 + ,31.41 + ,99.75 + ,5537.64 + ,28.79 + ,97.71 + ,5442.78 + ,26.76 + ,94.99 + ,5491.50 + ,26.41 + ,96.41 + ,5501.20 + ,27.05 + ,96.28 + ,5658.08 + ,29.43 + ,100.14 + ,5686.16 + ,32.1 + ,99.9 + ,5801.24 + ,36.84 + ,102.87 + ,5678.40 + ,34.22 + ,107.37 + ,5793.68 + ,36.53 + ,115.68 + ,5866.10 + ,40.99 + ,124.33 + ,6087.27 + ,45.97 + ,128.44 + ,6058.70 + ,43.6 + ,130.19 + ,6171.63 + ,47.84 + ,148.4 + ,6385.55 + ,51.47 + ,169.14 + ,6180.05 + ,51.31 + ,153.98 + ,6159.65 + ,48.47 + ,163.13 + ,6271.59 + ,48.28 + ,165.4 + ,6365.59 + ,46.56 + ,166.35 + ,6420.76 + ,43.83 + ,173.73 + ,6628.97 + ,51.17 + ,174.23 + ,6731.85 + ,49.59 + ,177.04 + ,6884.40 + ,49.11 + ,170.78 + ,6927.25 + ,49.97 + ,174.01 + ,6796.18 + ,50.07 + ,183.76 + ,6903.69 + ,53.3 + ,201.95 + ,7189.46 + ,57.08 + ,205.38 + ,7435.08 + ,68.54 + ,197.36 + ,7379.99 + ,71.62 + ,196.53 + ,7174.80 + ,67.64 + ,179.94 + ,7233.46 + ,64.79 + ,174.84 + ,7597.58 + ,80.97 + ,179.86 + ,7790.04 + ,88.42 + ,172.77 + ,7466.77 + ,110.22 + ,162.56 + ,7350.45 + ,99 + ,178.4 + ,6850.64 + ,95.95 + ,190.83 + ,6717.18 + ,107.94 + ,201.07 + ,6600.54 + ,97.82 + ,198.95 + ,6979.58 + ,111.64 + ,190.46 + ,6971.45 + ,114.73 + ,186.27 + ,6411.48 + ,117.58 + ,187.96 + ,6276.01 + ,99.19 + ,174.99 + ,6190.56 + ,90.19 + ,164.1 + ,5618.64 + ,59.74 + ,131.48 + ,4595.00 + ,44.51 + ,116.14 + ,4287.38 + ,23.94 + ,103.43 + ,4385.35 + ,21.29 + ,96.87 + ,3927.01 + ,20.77 + ,93.68 + ,3495.28 + ,25.07 + ,96.49 + ,3757.45 + ,32.95 + ,105.22 + ,4076.75 + ,40.05 + ,110.11 + ,4475.36 + ,44.59 + ,118.47 + ,4396.44 + ,40.28 + ,122.15 + ,4766.91 + ,41.19 + ,137.35 + ,4907.38 + ,38.14 + ,134.83 + ,5069.67 + ,41.85 + ,138.34 + ,4983.56 + ,43.76 + ,141.98 + ,5245.81 + ,50.16 + ,149.45 + ,5241.07 + ,52.94 + ,154.68 + ,5003.74 + ,47.69 + ,145.98 + ,5075.38 + ,51.52 + ,156.33 + ,5358.33 + ,58.69 + ,176.28 + ,5298.80 + ,50.44 + ,159.08 + ,4784.64 + ,45.72 + ,151.18 + ,4579.82 + ,43.24 + ,162.63 + ,4982.42 + ,51.49 + ,174.2 + ,4776.39 + ,50.43 + ,180.51 + ,5157.48 + ,58.73 + ,185.31 + ,5305.21 + ,65.12 + ,186.33 + ,5187.20 + ,64.13) + ,dim=c(3 + ,101) + ,dimnames=list(c('Commodity' + ,'WorldLeaders' + ,'RioTinto') + ,1:101)) > y <- array(NA,dim=c(3,101),dimnames=list(c('Commodity','WorldLeaders','RioTinto'),1:101)) > 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 = '3' > #'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 RioTinto Commodity WorldLeaders 1 14.36 54.64 4606.07 2 14.62 52.39 4176.60 3 13.51 52.51 4331.53 4 14.95 52.92 3987.30 5 16.72 55.22 4205.86 6 16.33 55.41 4331.37 7 15.21 57.02 4106.30 8 16.69 58.55 4009.33 9 15.81 57.49 3857.64 10 16.02 55.52 3929.00 11 16.70 57.84 4210.47 12 15.99 58.69 4445.82 13 17.68 59.74 4497.07 14 18.89 60.70 4443.71 15 18.72 60.74 4529.25 16 21.14 64.32 4634.22 17 20.97 66.90 4772.67 18 23.75 70.93 4881.52 19 23.05 75.89 5153.13 20 23.45 80.60 5324.19 21 21.74 81.39 5209.36 22 19.37 81.33 5108.92 23 21.10 77.04 5130.88 24 21.20 79.54 5195.68 25 22.67 81.93 5050.42 26 22.24 80.79 5101.03 27 23.78 81.98 5139.84 28 23.27 85.94 5234.31 29 25.74 86.60 5435.73 30 26.10 87.42 5633.57 31 27.49 93.14 5498.28 32 31.41 95.76 5668.13 33 28.79 99.75 5537.64 34 26.76 97.71 5442.78 35 26.41 94.99 5491.50 36 27.05 96.41 5501.20 37 29.43 96.28 5658.08 38 32.10 100.14 5686.16 39 36.84 99.90 5801.24 40 34.22 102.87 5678.40 41 36.53 107.37 5793.68 42 40.99 115.68 5866.10 43 45.97 124.33 6087.27 44 43.60 128.44 6058.70 45 47.84 130.19 6171.63 46 51.47 148.40 6385.55 47 51.31 169.14 6180.05 48 48.47 153.98 6159.65 49 48.28 163.13 6271.59 50 46.56 165.40 6365.59 51 43.83 166.35 6420.76 52 51.17 173.73 6628.97 53 49.59 174.23 6731.85 54 49.11 177.04 6884.40 55 49.97 170.78 6927.25 56 50.07 174.01 6796.18 57 53.30 183.76 6903.69 58 57.08 201.95 7189.46 59 68.54 205.38 7435.08 60 71.62 197.36 7379.99 61 67.64 196.53 7174.80 62 64.79 179.94 7233.46 63 80.97 174.84 7597.58 64 88.42 179.86 7790.04 65 110.22 172.77 7466.77 66 99.00 162.56 7350.45 67 95.95 178.40 6850.64 68 107.94 190.83 6717.18 69 97.82 201.07 6600.54 70 111.64 198.95 6979.58 71 114.73 190.46 6971.45 72 117.58 186.27 6411.48 73 99.19 187.96 6276.01 74 90.19 174.99 6190.56 75 59.74 164.10 5618.64 76 44.51 131.48 4595.00 77 23.94 116.14 4287.38 78 21.29 103.43 4385.35 79 20.77 96.87 3927.01 80 25.07 93.68 3495.28 81 32.95 96.49 3757.45 82 40.05 105.22 4076.75 83 44.59 110.11 4475.36 84 40.28 118.47 4396.44 85 41.19 122.15 4766.91 86 38.14 137.35 4907.38 87 41.85 134.83 5069.67 88 43.76 138.34 4983.56 89 50.16 141.98 5245.81 90 52.94 149.45 5241.07 91 47.69 154.68 5003.74 92 51.52 145.98 5075.38 93 58.69 156.33 5358.33 94 50.44 176.28 5298.80 95 45.72 159.08 4784.64 96 43.24 151.18 4579.82 97 51.49 162.63 4982.42 98 50.43 174.20 4776.39 99 58.73 180.51 5157.48 100 65.12 185.31 5305.21 101 64.13 186.33 5187.20 > k <- length(x[1,]) > df <- as.data.frame(x) > (mylm <- lm(df)) Call: lm(formula = df) Coefficients: (Intercept) Commodity WorldLeaders -32.524900 0.377868 0.005352 > (mysum <- summary(mylm)) Call: lm(formula = df) Residuals: Min 1Q Median 3Q Max -25.185 -6.257 -2.518 4.508 45.404 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -32.524900 7.551553 -4.307 3.93e-05 *** Commodity 0.377868 0.041498 9.106 1.06e-14 *** WorldLeaders 0.005352 0.001911 2.800 0.00615 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 13.19 on 98 degrees of freedom Multiple R-squared: 0.7498, Adjusted R-squared: 0.7447 F-statistic: 146.9 on 2 and 98 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,] 3.313819e-05 6.627639e-05 9.999669e-01 [2,] 7.029529e-05 1.405906e-04 9.999297e-01 [3,] 4.767847e-06 9.535693e-06 9.999952e-01 [4,] 3.887130e-07 7.774261e-07 9.999996e-01 [5,] 2.472582e-08 4.945164e-08 1.000000e+00 [6,] 1.592963e-09 3.185926e-09 1.000000e+00 [7,] 9.422899e-11 1.884580e-10 1.000000e+00 [8,] 9.776076e-12 1.955215e-11 1.000000e+00 [9,] 1.860884e-12 3.721767e-12 1.000000e+00 [10,] 1.615346e-13 3.230693e-13 1.000000e+00 [11,] 2.160524e-14 4.321047e-14 1.000000e+00 [12,] 1.541657e-15 3.083314e-15 1.000000e+00 [13,] 9.199890e-17 1.839978e-16 1.000000e+00 [14,] 1.505909e-16 3.011819e-16 1.000000e+00 [15,] 2.290212e-16 4.580425e-16 1.000000e+00 [16,] 6.745642e-16 1.349128e-15 1.000000e+00 [17,] 4.282791e-15 8.565582e-15 1.000000e+00 [18,] 5.596170e-16 1.119234e-15 1.000000e+00 [19,] 8.435459e-17 1.687092e-16 1.000000e+00 [20,] 8.822909e-18 1.764582e-17 1.000000e+00 [21,] 9.139401e-19 1.827880e-18 1.000000e+00 [22,] 1.068056e-19 2.136111e-19 1.000000e+00 [23,] 1.120277e-20 2.240554e-20 1.000000e+00 [24,] 1.647569e-21 3.295138e-21 1.000000e+00 [25,] 2.043050e-22 4.086100e-22 1.000000e+00 [26,] 2.813608e-23 5.627216e-23 1.000000e+00 [27,] 6.502548e-23 1.300510e-22 1.000000e+00 [28,] 6.583907e-24 1.316781e-23 1.000000e+00 [29,] 8.531816e-25 1.706363e-24 1.000000e+00 [30,] 9.400363e-26 1.880073e-25 1.000000e+00 [31,] 9.372514e-27 1.874503e-26 1.000000e+00 [32,] 1.596922e-27 3.193845e-27 1.000000e+00 [33,] 1.087314e-27 2.174629e-27 1.000000e+00 [34,] 3.354834e-25 6.709669e-25 1.000000e+00 [35,] 1.784203e-25 3.568407e-25 1.000000e+00 [36,] 1.226821e-25 2.453643e-25 1.000000e+00 [37,] 1.692800e-25 3.385600e-25 1.000000e+00 [38,] 3.358679e-25 6.717359e-25 1.000000e+00 [39,] 4.781036e-26 9.562073e-26 1.000000e+00 [40,] 2.160926e-26 4.321851e-26 1.000000e+00 [41,] 3.265628e-27 6.531255e-27 1.000000e+00 [42,] 5.163398e-26 1.032680e-25 1.000000e+00 [43,] 1.333113e-26 2.666225e-26 1.000000e+00 [44,] 1.441593e-26 2.883186e-26 1.000000e+00 [45,] 4.680573e-26 9.361147e-26 1.000000e+00 [46,] 7.445330e-25 1.489066e-24 1.000000e+00 [47,] 3.685284e-25 7.370568e-25 1.000000e+00 [48,] 4.374518e-25 8.749036e-25 1.000000e+00 [49,] 1.402874e-24 2.805748e-24 1.000000e+00 [50,] 1.684258e-24 3.368516e-24 1.000000e+00 [51,] 2.619817e-24 5.239634e-24 1.000000e+00 [52,] 5.497247e-24 1.099449e-23 1.000000e+00 [53,] 6.665158e-23 1.333032e-22 1.000000e+00 [54,] 2.638359e-21 5.276718e-21 1.000000e+00 [55,] 1.295298e-18 2.590595e-18 1.000000e+00 [56,] 1.006219e-16 2.012437e-16 1.000000e+00 [57,] 3.254948e-14 6.509895e-14 1.000000e+00 [58,] 7.380214e-09 1.476043e-08 1.000000e+00 [59,] 8.430267e-05 1.686053e-04 9.999157e-01 [60,] 9.810429e-02 1.962086e-01 9.018957e-01 [61,] 3.447498e-01 6.894996e-01 6.552502e-01 [62,] 6.211660e-01 7.576679e-01 3.788340e-01 [63,] 8.686567e-01 2.626865e-01 1.313433e-01 [64,] 9.023770e-01 1.952461e-01 9.762305e-02 [65,] 9.406714e-01 1.186573e-01 5.932865e-02 [66,] 9.693417e-01 6.131655e-02 3.065828e-02 [67,] 9.993565e-01 1.287079e-03 6.435393e-04 [68,] 9.999262e-01 1.475364e-04 7.376822e-05 [69,] 9.999980e-01 4.029549e-06 2.014774e-06 [70,] 9.999953e-01 9.371532e-06 4.685766e-06 [71,] 9.999898e-01 2.048481e-05 1.024240e-05 [72,] 9.999951e-01 9.745862e-06 4.872931e-06 [73,] 9.999988e-01 2.344644e-06 1.172322e-06 [74,] 9.999994e-01 1.104959e-06 5.524795e-07 [75,] 9.999987e-01 2.660552e-06 1.330276e-06 [76,] 9.999960e-01 8.028711e-06 4.014356e-06 [77,] 9.999947e-01 1.059607e-05 5.298035e-06 [78,] 9.999981e-01 3.798737e-06 1.899368e-06 [79,] 9.999988e-01 2.375789e-06 1.187894e-06 [80,] 9.999979e-01 4.181655e-06 2.090827e-06 [81,] 9.999960e-01 7.980048e-06 3.990024e-06 [82,] 9.999896e-01 2.080733e-05 1.040367e-05 [83,] 9.999620e-01 7.591210e-05 3.795605e-05 [84,] 9.998542e-01 2.916547e-04 1.458274e-04 [85,] 9.994569e-01 1.086121e-03 5.430603e-04 [86,] 9.983855e-01 3.228904e-03 1.614452e-03 [87,] 9.951310e-01 9.737951e-03 4.868976e-03 [88,] 9.934229e-01 1.315422e-02 6.577111e-03 [89,] 9.993507e-01 1.298599e-03 6.492994e-04 [90,] 9.966024e-01 6.795201e-03 3.397600e-03 > postscript(file="/var/www/html/freestat/rcomp/tmp/1ik2i1292693115.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/2ik2i1292693115.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/3ik2i1292693115.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/4tt1l1292693115.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/5tt1l1292693115.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 = 101 Frequency = 1 1 2 3 4 5 6 1.5857390 4.9945360 3.0099813 6.1374302 5.8685644 4.7350195 7 8 9 10 11 12 4.2112627 5.6321234 5.9645334 6.5370037 4.8338759 2.5430564 13 14 15 16 17 18 3.5619963 4.6948341 4.0518954 4.5573103 2.6714031 3.3460109 19 20 21 22 23 24 -0.6819169 -2.9772176 -4.3711448 -6.1809015 -2.9473799 -4.1388708 25 26 27 28 29 30 -2.7945209 -3.0646240 -2.1820047 -4.6939821 -3.5514078 -4.5601317 31 32 33 34 35 36 -4.6074444 -2.5865242 -6.0158151 -6.7672576 -6.3502131 -6.2987021 37 38 39 40 41 42 -4.7092265 -3.6480869 0.5666746 -2.5181347 -2.5255394 -1.5932288 43 44 45 46 47 48 -1.0655276 -4.8356551 -1.8613444 -6.2572611 -13.1543807 -10.1567128 49 50 51 52 53 54 -14.4033290 -17.4841934 -20.8684471 -17.4314891 -19.7510538 -22.1093362 55 56 57 58 59 60 -19.1132206 -19.5322274 -20.5618545 -25.1847668 -16.3354533 -9.9300987 61 62 63 64 65 66 -12.4982578 -9.3933802 6.7649188 11.2879426 37.4972226 30.7578217 67 68 69 70 71 72 24.3974519 32.4048484 19.0397530 31.6321501 37.9737652 45.4040840 73 74 75 76 77 78 27.1005440 23.4588384 0.1848333 2.7595857 -10.3674820 -8.7391270 79 80 81 82 83 84 -4.3272006 3.4888885 8.9039021 10.9961661 11.5549645 4.5083780 85 86 87 88 89 90 2.0450069 -7.5004099 -3.7067842 -2.6622273 0.9587274 0.9414203 91 92 93 94 95 96 -5.0146023 1.7194232 3.4640916 -12.0057671 -7.4745640 -5.8731744 97 98 99 100 101 -4.1045473 -8.4337777 -4.5577825 -0.7722254 -1.5160424 > postscript(file="/var/www/html/freestat/rcomp/tmp/6tt1l1292693115.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 = 101 Frequency = 1 lag(myerror, k = 1) myerror 0 1.5857390 NA 1 4.9945360 1.5857390 2 3.0099813 4.9945360 3 6.1374302 3.0099813 4 5.8685644 6.1374302 5 4.7350195 5.8685644 6 4.2112627 4.7350195 7 5.6321234 4.2112627 8 5.9645334 5.6321234 9 6.5370037 5.9645334 10 4.8338759 6.5370037 11 2.5430564 4.8338759 12 3.5619963 2.5430564 13 4.6948341 3.5619963 14 4.0518954 4.6948341 15 4.5573103 4.0518954 16 2.6714031 4.5573103 17 3.3460109 2.6714031 18 -0.6819169 3.3460109 19 -2.9772176 -0.6819169 20 -4.3711448 -2.9772176 21 -6.1809015 -4.3711448 22 -2.9473799 -6.1809015 23 -4.1388708 -2.9473799 24 -2.7945209 -4.1388708 25 -3.0646240 -2.7945209 26 -2.1820047 -3.0646240 27 -4.6939821 -2.1820047 28 -3.5514078 -4.6939821 29 -4.5601317 -3.5514078 30 -4.6074444 -4.5601317 31 -2.5865242 -4.6074444 32 -6.0158151 -2.5865242 33 -6.7672576 -6.0158151 34 -6.3502131 -6.7672576 35 -6.2987021 -6.3502131 36 -4.7092265 -6.2987021 37 -3.6480869 -4.7092265 38 0.5666746 -3.6480869 39 -2.5181347 0.5666746 40 -2.5255394 -2.5181347 41 -1.5932288 -2.5255394 42 -1.0655276 -1.5932288 43 -4.8356551 -1.0655276 44 -1.8613444 -4.8356551 45 -6.2572611 -1.8613444 46 -13.1543807 -6.2572611 47 -10.1567128 -13.1543807 48 -14.4033290 -10.1567128 49 -17.4841934 -14.4033290 50 -20.8684471 -17.4841934 51 -17.4314891 -20.8684471 52 -19.7510538 -17.4314891 53 -22.1093362 -19.7510538 54 -19.1132206 -22.1093362 55 -19.5322274 -19.1132206 56 -20.5618545 -19.5322274 57 -25.1847668 -20.5618545 58 -16.3354533 -25.1847668 59 -9.9300987 -16.3354533 60 -12.4982578 -9.9300987 61 -9.3933802 -12.4982578 62 6.7649188 -9.3933802 63 11.2879426 6.7649188 64 37.4972226 11.2879426 65 30.7578217 37.4972226 66 24.3974519 30.7578217 67 32.4048484 24.3974519 68 19.0397530 32.4048484 69 31.6321501 19.0397530 70 37.9737652 31.6321501 71 45.4040840 37.9737652 72 27.1005440 45.4040840 73 23.4588384 27.1005440 74 0.1848333 23.4588384 75 2.7595857 0.1848333 76 -10.3674820 2.7595857 77 -8.7391270 -10.3674820 78 -4.3272006 -8.7391270 79 3.4888885 -4.3272006 80 8.9039021 3.4888885 81 10.9961661 8.9039021 82 11.5549645 10.9961661 83 4.5083780 11.5549645 84 2.0450069 4.5083780 85 -7.5004099 2.0450069 86 -3.7067842 -7.5004099 87 -2.6622273 -3.7067842 88 0.9587274 -2.6622273 89 0.9414203 0.9587274 90 -5.0146023 0.9414203 91 1.7194232 -5.0146023 92 3.4640916 1.7194232 93 -12.0057671 3.4640916 94 -7.4745640 -12.0057671 95 -5.8731744 -7.4745640 96 -4.1045473 -5.8731744 97 -8.4337777 -4.1045473 98 -4.5577825 -8.4337777 99 -0.7722254 -4.5577825 100 -1.5160424 -0.7722254 101 NA -1.5160424 > dum1 <- dum[2:length(myerror),] > dum1 lag(myerror, k = 1) myerror [1,] 4.9945360 1.5857390 [2,] 3.0099813 4.9945360 [3,] 6.1374302 3.0099813 [4,] 5.8685644 6.1374302 [5,] 4.7350195 5.8685644 [6,] 4.2112627 4.7350195 [7,] 5.6321234 4.2112627 [8,] 5.9645334 5.6321234 [9,] 6.5370037 5.9645334 [10,] 4.8338759 6.5370037 [11,] 2.5430564 4.8338759 [12,] 3.5619963 2.5430564 [13,] 4.6948341 3.5619963 [14,] 4.0518954 4.6948341 [15,] 4.5573103 4.0518954 [16,] 2.6714031 4.5573103 [17,] 3.3460109 2.6714031 [18,] -0.6819169 3.3460109 [19,] -2.9772176 -0.6819169 [20,] -4.3711448 -2.9772176 [21,] -6.1809015 -4.3711448 [22,] -2.9473799 -6.1809015 [23,] -4.1388708 -2.9473799 [24,] -2.7945209 -4.1388708 [25,] -3.0646240 -2.7945209 [26,] -2.1820047 -3.0646240 [27,] -4.6939821 -2.1820047 [28,] -3.5514078 -4.6939821 [29,] -4.5601317 -3.5514078 [30,] -4.6074444 -4.5601317 [31,] -2.5865242 -4.6074444 [32,] -6.0158151 -2.5865242 [33,] -6.7672576 -6.0158151 [34,] -6.3502131 -6.7672576 [35,] -6.2987021 -6.3502131 [36,] -4.7092265 -6.2987021 [37,] -3.6480869 -4.7092265 [38,] 0.5666746 -3.6480869 [39,] -2.5181347 0.5666746 [40,] -2.5255394 -2.5181347 [41,] -1.5932288 -2.5255394 [42,] -1.0655276 -1.5932288 [43,] -4.8356551 -1.0655276 [44,] -1.8613444 -4.8356551 [45,] -6.2572611 -1.8613444 [46,] -13.1543807 -6.2572611 [47,] -10.1567128 -13.1543807 [48,] -14.4033290 -10.1567128 [49,] -17.4841934 -14.4033290 [50,] -20.8684471 -17.4841934 [51,] -17.4314891 -20.8684471 [52,] -19.7510538 -17.4314891 [53,] -22.1093362 -19.7510538 [54,] -19.1132206 -22.1093362 [55,] -19.5322274 -19.1132206 [56,] -20.5618545 -19.5322274 [57,] -25.1847668 -20.5618545 [58,] -16.3354533 -25.1847668 [59,] -9.9300987 -16.3354533 [60,] -12.4982578 -9.9300987 [61,] -9.3933802 -12.4982578 [62,] 6.7649188 -9.3933802 [63,] 11.2879426 6.7649188 [64,] 37.4972226 11.2879426 [65,] 30.7578217 37.4972226 [66,] 24.3974519 30.7578217 [67,] 32.4048484 24.3974519 [68,] 19.0397530 32.4048484 [69,] 31.6321501 19.0397530 [70,] 37.9737652 31.6321501 [71,] 45.4040840 37.9737652 [72,] 27.1005440 45.4040840 [73,] 23.4588384 27.1005440 [74,] 0.1848333 23.4588384 [75,] 2.7595857 0.1848333 [76,] -10.3674820 2.7595857 [77,] -8.7391270 -10.3674820 [78,] -4.3272006 -8.7391270 [79,] 3.4888885 -4.3272006 [80,] 8.9039021 3.4888885 [81,] 10.9961661 8.9039021 [82,] 11.5549645 10.9961661 [83,] 4.5083780 11.5549645 [84,] 2.0450069 4.5083780 [85,] -7.5004099 2.0450069 [86,] -3.7067842 -7.5004099 [87,] -2.6622273 -3.7067842 [88,] 0.9587274 -2.6622273 [89,] 0.9414203 0.9587274 [90,] -5.0146023 0.9414203 [91,] 1.7194232 -5.0146023 [92,] 3.4640916 1.7194232 [93,] -12.0057671 3.4640916 [94,] -7.4745640 -12.0057671 [95,] -5.8731744 -7.4745640 [96,] -4.1045473 -5.8731744 [97,] -8.4337777 -4.1045473 [98,] -4.5577825 -8.4337777 [99,] -0.7722254 -4.5577825 [100,] -1.5160424 -0.7722254 > z <- as.data.frame(dum1) > z lag(myerror, k = 1) myerror 1 4.9945360 1.5857390 2 3.0099813 4.9945360 3 6.1374302 3.0099813 4 5.8685644 6.1374302 5 4.7350195 5.8685644 6 4.2112627 4.7350195 7 5.6321234 4.2112627 8 5.9645334 5.6321234 9 6.5370037 5.9645334 10 4.8338759 6.5370037 11 2.5430564 4.8338759 12 3.5619963 2.5430564 13 4.6948341 3.5619963 14 4.0518954 4.6948341 15 4.5573103 4.0518954 16 2.6714031 4.5573103 17 3.3460109 2.6714031 18 -0.6819169 3.3460109 19 -2.9772176 -0.6819169 20 -4.3711448 -2.9772176 21 -6.1809015 -4.3711448 22 -2.9473799 -6.1809015 23 -4.1388708 -2.9473799 24 -2.7945209 -4.1388708 25 -3.0646240 -2.7945209 26 -2.1820047 -3.0646240 27 -4.6939821 -2.1820047 28 -3.5514078 -4.6939821 29 -4.5601317 -3.5514078 30 -4.6074444 -4.5601317 31 -2.5865242 -4.6074444 32 -6.0158151 -2.5865242 33 -6.7672576 -6.0158151 34 -6.3502131 -6.7672576 35 -6.2987021 -6.3502131 36 -4.7092265 -6.2987021 37 -3.6480869 -4.7092265 38 0.5666746 -3.6480869 39 -2.5181347 0.5666746 40 -2.5255394 -2.5181347 41 -1.5932288 -2.5255394 42 -1.0655276 -1.5932288 43 -4.8356551 -1.0655276 44 -1.8613444 -4.8356551 45 -6.2572611 -1.8613444 46 -13.1543807 -6.2572611 47 -10.1567128 -13.1543807 48 -14.4033290 -10.1567128 49 -17.4841934 -14.4033290 50 -20.8684471 -17.4841934 51 -17.4314891 -20.8684471 52 -19.7510538 -17.4314891 53 -22.1093362 -19.7510538 54 -19.1132206 -22.1093362 55 -19.5322274 -19.1132206 56 -20.5618545 -19.5322274 57 -25.1847668 -20.5618545 58 -16.3354533 -25.1847668 59 -9.9300987 -16.3354533 60 -12.4982578 -9.9300987 61 -9.3933802 -12.4982578 62 6.7649188 -9.3933802 63 11.2879426 6.7649188 64 37.4972226 11.2879426 65 30.7578217 37.4972226 66 24.3974519 30.7578217 67 32.4048484 24.3974519 68 19.0397530 32.4048484 69 31.6321501 19.0397530 70 37.9737652 31.6321501 71 45.4040840 37.9737652 72 27.1005440 45.4040840 73 23.4588384 27.1005440 74 0.1848333 23.4588384 75 2.7595857 0.1848333 76 -10.3674820 2.7595857 77 -8.7391270 -10.3674820 78 -4.3272006 -8.7391270 79 3.4888885 -4.3272006 80 8.9039021 3.4888885 81 10.9961661 8.9039021 82 11.5549645 10.9961661 83 4.5083780 11.5549645 84 2.0450069 4.5083780 85 -7.5004099 2.0450069 86 -3.7067842 -7.5004099 87 -2.6622273 -3.7067842 88 0.9587274 -2.6622273 89 0.9414203 0.9587274 90 -5.0146023 0.9414203 91 1.7194232 -5.0146023 92 3.4640916 1.7194232 93 -12.0057671 3.4640916 94 -7.4745640 -12.0057671 95 -5.8731744 -7.4745640 96 -4.1045473 -5.8731744 97 -8.4337777 -4.1045473 98 -4.5577825 -8.4337777 99 -0.7722254 -4.5577825 100 -1.5160424 -0.7722254 > 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/7ml161292693115.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/8ml161292693115.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/9fc091292693115.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/10fc091292693115.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/110vhx1292693115.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/12mvxl1292693115.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/13awuw1292693115.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/14l5bz1292693115.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/15o6sn1292693115.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/16ky7w1292693115.tab") + } > > try(system("convert tmp/1ik2i1292693115.ps tmp/1ik2i1292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/2ik2i1292693115.ps tmp/2ik2i1292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/3ik2i1292693115.ps tmp/3ik2i1292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/4tt1l1292693115.ps tmp/4tt1l1292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/5tt1l1292693115.ps tmp/5tt1l1292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/6tt1l1292693115.ps tmp/6tt1l1292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/7ml161292693115.ps tmp/7ml161292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/8ml161292693115.ps tmp/8ml161292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/9fc091292693115.ps tmp/9fc091292693115.png",intern=TRUE)) character(0) > try(system("convert tmp/10fc091292693115.ps tmp/10fc091292693115.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 4.422 2.538 4.828