R version 2.8.1 (2008-12-22) 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. 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. > y <- c(1169,2154,2249,2687,4359,5382,4459,6398,4596,3024,1887,2070,1351,2218,2461,3028,4784,4975,4607,6249,4809,3157,1910,2228,1594,2467,2222,3607,4685,4962,5770,5480,5000,3228,1993,2288,1580,2111,2192,3601,4665,4876,5813,5589,5331,3075,2002,2306,1507,1992,2487,3490,4647,5594,5611,5788,6204,3013,1931,2549) > x <- c(1.00,1.00,0.99,0.97,0.96,0.96,1.09,1.10,1.10,1.09,1.07,1.08,1.09,1.09,1.07,1.07,1.05,1.06,1.20,1.22,1.21,1.19,1.15,1.17,1.18,1.18,1.17,1.14,1.13,1.14,1.25,1.28,1.29,1.28,1.25,1.25,1.24,1.24,1.23,1.20,1.19,1.20,1.30,1.32,1.32,1.29,1.25,1.25,1.25,1.24,1.22,1.21,1.20,1.20,1.30,1.31,1.30,1.23,1.19,1.17) > par2 = '36' > par1 = '0' > par1 <- as.numeric(par1) > par2 <- as.numeric(par2) > x <- as.ts(x) > y <- as.ts(y) > mylm <- lm(y~x) > cbind(mylm$resid) [,1] 1 -1772.81652 2 -787.81652 3 -657.18296 4 -147.91582 5 1559.71774 6 2582.71774 7 1196.48138 8 3099.84781 9 1297.84781 10 -238.51862 11 -1304.25149 12 -1156.88506 13 -1911.51862 14 -1044.51862 15 -730.25149 16 -163.25149 17 1664.01564 18 1819.38208 19 952.51214 20 2523.24501 21 1118.87858 22 -461.85429 23 -1566.32002 24 -1319.58716 25 -1989.22072 26 -1116.22072 27 -1325.58716 28 166.31354 29 1279.94711 30 1521.31354 31 1937.34431 32 1540.44361 33 1024.81004 34 -711.55639 35 -1839.65569 36 -1544.65569 37 -2217.02212 38 -1686.02212 39 -1569.38856 40 -53.48786 41 1046.14571 42 1221.51214 43 1802.17648 44 1506.90934 45 1248.90934 46 -900.18996 47 -1830.65569 48 -1526.65569 49 -2325.65569 50 -1805.02212 51 -1238.75499 52 -200.12142 53 992.51214 54 1939.51214 55 1600.17648 56 1741.54291 57 2193.17648 58 -748.38856 59 -1687.85429 60 -998.58716 > library(lattice) > postscript(file="/var/www/rcomp/tmp/1n10z1257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(y,type='l',main='Run Sequence Plot of Y[t]',xlab='time or index',ylab='value') > grid() > dev.off() null device 1 > postscript(file="/var/www/rcomp/tmp/2uvf11257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(x,type='l',main='Run Sequence Plot of X[t]',xlab='time or index',ylab='value') > grid() > dev.off() null device 1 > postscript(file="/var/www/rcomp/tmp/3tl281257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(x,y,main='Scatter Plot',xlab='X[t]',ylab='Y[t]') > grid() > dev.off() null device 1 > postscript(file="/var/www/rcomp/tmp/4mz621257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(mylm$resid,type='l',main='Run Sequence Plot of e[t]',xlab='time or index',ylab='value') > grid() > dev.off() null device 1 > postscript(file="/var/www/rcomp/tmp/5s8og1257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > hist(mylm$resid,main='Histogram of e[t]') > dev.off() null device 1 > postscript(file="/var/www/rcomp/tmp/677ai1257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > if (par1 > 0) + { + densityplot(~mylm$resid,col='black',main=paste('Density Plot of e[t] bw = ',par1),bw=par1) + } else { + densityplot(~mylm$resid,col='black',main='Density Plot of e[t]') + } > dev.off() null device 1 > postscript(file="/var/www/rcomp/tmp/7uqcd1257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > qqnorm(mylm$resid,main='QQ plot of e[t]') > qqline(mylm$resid) > grid() > dev.off() null device 1 > if (par2 > 0) + { + postscript(file="/var/www/rcomp/tmp/8k69m1257360534.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) + acf(mylm$resid,lag.max=par2,main='Residual Autocorrelation Function') + grid() + dev.off() + } null device 1 > summary(x) Min. 1st Qu. Median Mean 3rd Qu. Max. 0.960 1.090 1.195 1.173 1.250 1.320 > > #Note: the /var/www/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/rcomp/createtable") > > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Model: Y[t] = c + b X[t] + e[t]',2,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'c',1,TRUE) > a<-table.element(a,mylm$coeff[[1]]) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'b',1,TRUE) > a<-table.element(a,mylm$coeff[[2]]) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/www/rcomp/tmp/9cvrr1257360534.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Descriptive Statistics about e[t]',2,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'# observations',header=TRUE) > a<-table.element(a,length(mylm$resid)) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'minimum',header=TRUE) > a<-table.element(a,min(mylm$resid)) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Q1',header=TRUE) > a<-table.element(a,quantile(mylm$resid,0.25)) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'median',header=TRUE) > a<-table.element(a,median(mylm$resid)) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'mean',header=TRUE) > a<-table.element(a,mean(mylm$resid)) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Q3',header=TRUE) > a<-table.element(a,quantile(mylm$resid,0.75)) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'maximum',header=TRUE) > a<-table.element(a,max(mylm$resid)) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/www/rcomp/tmp/10iukt1257360534.tab") > > system("convert tmp/1n10z1257360534.ps tmp/1n10z1257360534.png") > system("convert tmp/2uvf11257360534.ps tmp/2uvf11257360534.png") > system("convert tmp/3tl281257360534.ps tmp/3tl281257360534.png") > system("convert tmp/4mz621257360534.ps tmp/4mz621257360534.png") > system("convert tmp/5s8og1257360534.ps tmp/5s8og1257360534.png") > system("convert tmp/677ai1257360534.ps tmp/677ai1257360534.png") > system("convert tmp/7uqcd1257360534.ps tmp/7uqcd1257360534.png") > system("convert tmp/8k69m1257360534.ps tmp/8k69m1257360534.png") > > > proc.time() user system elapsed 1.990 1.700 3.517