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 <- c(8.7000,8.9000,8.9000,8.1000,8.0000,8.3000,8.5000,8.7000,8.6000,8.3000,7.9000,7.9000,8.1000,8.3000,8.1000,7.4000,7.3000,7.7000,8.0000,8.0000,7.7000,6.9000,6.6000,6.9000,7.5000,7.9000,7.7000,6.5000,6.1000,6.4000,6.8000,7.1000,7.3000,7.2000,7.0000,7.0000,7.0000,7.3000,7.5000,7.2000,7.7000,8.0000,7.9000,8.0000,8.0000,7.9000,7.9000,8.0000,8.1000,8.1000,8.2000,8.0000,8.3000,8.5000,8.6000,8.7000,8.7000,8.5000,8.4000,8.5000,8.7000,8.7000,8.6000,7.9000,8.1000,8.2000,8.5000,8.6000,8.5000,8.3000,8.2000,8.7000,9.3000,9.3000,8.8000,7.4000,7.2000,7.5000,8.3000,8.8000,8.9000,8.6000,8.4000,8.4000,8.4000,8.4000,8.3000,7.6000,7.6000,7.9000,8.0000,8.2000,8.3000,8.2000,8.1000,8.0000,7.8000,7.6000,7.5000,6.8000,6.9000,7.1000,7.3000,7.4000,7.6000,7.6000,7.5000,7.5000,6.8000,6.4000,6.2000,6.0000,6.3000,6.3000,6.1000,6.1000,6.3000) > par8 = 'FALSE' > par7 = '1' > par6 = '' > par5 = '1' > par4 = '' > par3 = '0' > par2 = 'periodic' > par1 = '12' > main = 'Seasonal Decomposition by Loess' > #'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!) > par1 <- as.numeric(par1) #seasonal period > if (par2 != 'periodic') par2 <- as.numeric(par2) #s.window > par3 <- as.numeric(par3) #s.degree > if (par4 == '') par4 <- NULL else par4 <- as.numeric(par4)#t.window > par5 <- as.numeric(par5)#t.degree > if (par6 != '') par6 <- as.numeric(par6)#l.window > par7 <- as.numeric(par7)#l.degree > if (par8 == 'FALSE') par8 <- FALSE else par9 <- TRUE #robust > nx <- length(x) > x <- ts(x,frequency=par1) > if (par6 != '') { + m <- stl(x,s.window=par2, s.degree=par3, t.window=par4, t.degre=par5, l.window=par6, l.degree=par7, robust=par8) + } else { + m <- stl(x,s.window=par2, s.degree=par3, t.window=par4, t.degre=par5, l.degree=par7, robust=par8) + } > m$time.series seasonal trend remainder Jan 1 0.15117940 8.746560 -0.197739403 Feb 1 0.22500546 8.686979 -0.011984504 Mar 1 0.13883144 8.627398 0.133770470 Apr 1 -0.53032407 8.566986 0.063338404 May 1 -0.44947946 8.506573 -0.057093786 Jun 1 -0.18627643 8.446141 0.040135719 Jul 1 0.04692670 8.385708 0.067365121 Aug 1 0.23572918 8.325330 0.138940644 Sep 1 0.29453158 8.264952 0.040516242 Oct 1 0.07947905 8.203032 0.017488607 Nov 1 -0.06431358 8.141113 -0.176798928 Dec 1 0.05871075 8.087073 -0.245783861 Jan 2 0.15117940 8.033034 -0.084213110 Feb 2 0.22500546 7.974739 0.100255814 Mar 2 0.13883144 7.916444 0.044724812 Apr 2 -0.53032407 7.836415 0.093908802 May 2 -0.44947946 7.756387 -0.006907331 Jun 2 -0.18627643 7.670954 0.215322453 Jul 2 0.04692670 7.585521 0.367552134 Aug 2 0.23572918 7.517885 0.246385338 Sep 2 0.29453158 7.450250 -0.044781382 Oct 2 0.07947905 7.379371 -0.558849874 Nov 2 -0.06431358 7.308492 -0.644178265 Dec 2 0.05871075 7.225803 -0.384514191 Jan 3 0.15117940 7.143115 0.205705566 Feb 3 0.22500546 7.087588 0.587406871 Mar 3 0.13883144 7.032060 0.529108252 Apr 3 -0.53032407 7.019736 0.010587901 May 3 -0.44947946 7.007412 -0.457932572 Jun 3 -0.18627643 6.993741 -0.407464693 Jul 3 0.04692670 6.980070 -0.226996917 Aug 3 0.23572918 6.969128 -0.104857014 Sep 3 0.29453158 6.958185 0.047282964 Oct 3 0.07947905 7.010727 0.109793927 Nov 3 -0.06431358 7.063269 0.001044991 Dec 3 0.05871075 7.165127 -0.223838210 Jan 4 0.15117940 7.266986 -0.418165727 Feb 4 0.22500546 7.359313 -0.284318683 Mar 4 0.13883144 7.451640 -0.090471564 Apr 4 -0.53032407 7.531999 0.198324580 May 4 -0.44947946 7.612359 0.537120600 Jun 4 -0.18627643 7.689669 0.496607029 Jul 4 0.04692670 7.766980 0.086093354 Aug 4 0.23572918 7.826305 -0.062034417 Sep 4 0.29453158 7.885631 -0.180162113 Oct 4 0.07947905 7.932534 -0.112013369 Nov 4 -0.06431358 7.979438 -0.015124524 Dec 4 0.05871075 8.035718 -0.094428310 Jan 5 0.15117940 8.091997 -0.143176413 Feb 5 0.22500546 8.157856 -0.282861092 Mar 5 0.13883144 8.223714 -0.162545696 Apr 5 -0.53032407 8.281227 0.249096740 May 5 -0.44947946 8.338740 0.410739053 Jun 5 -0.18627643 8.384952 0.301324272 Jul 5 0.04692670 8.431164 0.121909388 Aug 5 0.23572918 8.459250 0.005020712 Sep 5 0.29453158 8.487336 -0.081867889 Oct 5 0.07947905 8.488116 -0.067595057 Nov 5 -0.06431358 8.488896 -0.024582124 Dec 5 0.05871075 8.477617 -0.036328056 Jan 6 0.15117940 8.466339 0.082481695 Feb 6 0.22500546 8.456078 0.018916203 Mar 6 0.13883144 8.445818 0.015350787 Apr 6 -0.53032407 8.434056 -0.003731912 May 6 -0.44947946 8.422294 0.127185265 Jun 6 -0.18627643 8.429299 -0.043022481 Jul 6 0.04692670 8.436304 0.016769669 Aug 6 0.23572918 8.458254 -0.093983567 Sep 6 0.29453158 8.480205 -0.274736729 Oct 6 0.07947905 8.472911 -0.252389642 Nov 6 -0.06431358 8.465616 -0.201302455 Dec 6 0.05871075 8.435614 0.205674861 Jan 7 0.15117940 8.405613 0.743207860 Feb 7 0.22500546 8.395051 0.679943534 Mar 7 0.13883144 8.384489 0.276679284 Apr 7 -0.53032407 8.382244 -0.451920251 May 7 -0.44947946 8.379999 -0.730519908 Jun 7 -0.18627643 8.354722 -0.668445516 Jul 7 0.04692670 8.329445 -0.076371227 Aug 7 0.23572918 8.303191 0.261080257 Sep 7 0.29453158 8.276937 0.328531816 Oct 7 0.07947905 8.279281 0.241240206 Nov 7 -0.06431358 8.281625 0.182688696 Dec 7 0.05871075 8.276546 0.064742905 Jan 8 0.15117940 8.271468 -0.022647202 Feb 8 0.22500546 8.234331 -0.059336003 Mar 8 0.13883144 8.197193 -0.036024729 Apr 8 -0.53032407 8.158611 -0.028286759 May 8 -0.44947946 8.120028 -0.070548912 Jun 8 -0.18627643 8.084565 0.001711847 Jul 8 0.04692670 8.049101 -0.096027497 Aug 8 0.23572918 7.998973 -0.034702647 Sep 8 0.29453158 7.948846 0.056622278 Oct 8 0.07947905 7.885010 0.235511045 Nov 8 -0.06431358 7.821174 0.343139912 Dec 8 0.05871075 7.752047 0.189242199 Jan 9 0.15117940 7.682920 -0.034099830 Feb 9 0.22500546 7.612681 -0.237686879 Mar 9 0.13883144 7.542442 -0.181273854 Apr 9 -0.53032407 7.484696 -0.154371998 May 9 -0.44947946 7.426950 -0.077470266 Jun 9 -0.18627643 7.376169 -0.089892561 Jul 9 0.04692670 7.325388 -0.072314959 Aug 9 0.23572918 7.255850 -0.091579216 Sep 9 0.29453158 7.186312 0.119156603 Oct 9 0.07947905 7.105393 0.415127504 Nov 9 -0.06431358 7.024475 0.539838507 Dec 9 0.05871075 6.929764 0.511525110 Jan 10 0.15117940 6.835053 -0.186232603 Feb 10 0.22500546 6.722320 -0.547325326 Mar 10 0.13883144 6.609587 -0.548417974 Apr 10 -0.53032407 6.491871 0.038453106 May 10 -0.44947946 6.374155 0.375324063 Jun 10 -0.18627643 6.256518 0.229758330 Jul 10 0.04692670 6.138881 -0.085807505 Aug 10 0.23572918 6.022146 -0.157875183 Sep 10 0.29453158 5.905411 0.100057215 > m$win s t l 1171 19 13 > m$deg s t l 0 1 1 > m$jump s t l 118 2 2 > m$inner [1] 2 > m$outer [1] 0 > postscript(file="/var/www/html/freestat/rcomp/tmp/1puj21292786869.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(m,main=main) > dev.off() null device 1 > mylagmax <- nx/2 > postscript(file="/var/www/html/freestat/rcomp/tmp/2puj21292786869.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow = c(2,2)) > acf(as.numeric(x),lag.max = mylagmax,main='Observed') > acf(as.numeric(m$time.series[,'trend']),na.action=na.pass,lag.max = mylagmax,main='Trend') > acf(as.numeric(m$time.series[,'seasonal']),na.action=na.pass,lag.max = mylagmax,main='Seasonal') > acf(as.numeric(m$time.series[,'remainder']),na.action=na.pass,lag.max = mylagmax,main='Remainder') > par(op) > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/30lin1292786869.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow = c(2,2)) > spectrum(as.numeric(x),main='Observed') > spectrum(as.numeric(m$time.series[!is.na(m$time.series[,'trend']),'trend']),main='Trend') > spectrum(as.numeric(m$time.series[!is.na(m$time.series[,'seasonal']),'seasonal']),main='Seasonal') > spectrum(as.numeric(m$time.series[!is.na(m$time.series[,'remainder']),'remainder']),main='Remainder') > par(op) > dev.off() null device 1 > postscript(file="/var/www/html/freestat/rcomp/tmp/40lin1292786869.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow = c(2,2)) > cpgram(as.numeric(x),main='Observed') > cpgram(as.numeric(m$time.series[!is.na(m$time.series[,'trend']),'trend']),main='Trend') > cpgram(as.numeric(m$time.series[!is.na(m$time.series[,'seasonal']),'seasonal']),main='Seasonal') > cpgram(as.numeric(m$time.series[!is.na(m$time.series[,'remainder']),'remainder']),main='Remainder') > par(op) > 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,'Seasonal Decomposition by Loess - Parameters',4,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Component',header=TRUE) > a<-table.element(a,'Window',header=TRUE) > a<-table.element(a,'Degree',header=TRUE) > a<-table.element(a,'Jump',header=TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Seasonal',header=TRUE) > a<-table.element(a,m$win['s']) > a<-table.element(a,m$deg['s']) > a<-table.element(a,m$jump['s']) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Trend',header=TRUE) > a<-table.element(a,m$win['t']) > a<-table.element(a,m$deg['t']) > a<-table.element(a,m$jump['t']) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'Low-pass',header=TRUE) > a<-table.element(a,m$win['l']) > a<-table.element(a,m$deg['l']) > a<-table.element(a,m$jump['l']) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/www/html/freestat/rcomp/tmp/5vdxv1292786869.tab") > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Seasonal Decomposition by Loess - Time Series Components',6,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'t',header=TRUE) > a<-table.element(a,'Observed',header=TRUE) > a<-table.element(a,'Fitted',header=TRUE) > a<-table.element(a,'Seasonal',header=TRUE) > a<-table.element(a,'Trend',header=TRUE) > a<-table.element(a,'Remainder',header=TRUE) > a<-table.row.end(a) > for (i in 1:nx) { + a<-table.row.start(a) + a<-table.element(a,i,header=TRUE) + a<-table.element(a,x[i]) + a<-table.element(a,x[i]+m$time.series[i,'remainder']) + a<-table.element(a,m$time.series[i,'seasonal']) + a<-table.element(a,m$time.series[i,'trend']) + a<-table.element(a,m$time.series[i,'remainder']) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/www/html/freestat/rcomp/tmp/6zvw11292786869.tab") > > try(system("convert tmp/1puj21292786869.ps tmp/1puj21292786869.png",intern=TRUE)) character(0) > try(system("convert tmp/2puj21292786869.ps tmp/2puj21292786869.png",intern=TRUE)) character(0) > try(system("convert tmp/30lin1292786869.ps tmp/30lin1292786869.png",intern=TRUE)) character(0) > try(system("convert tmp/40lin1292786869.ps tmp/40lin1292786869.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 1.704 0.899 1.831