R version 2.13.0 (2011-04-13) Copyright (C) 2011 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i486-pc-linux-gnu (32-bit) 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 <- c(127,128.7,127.3,136.7,133.8,137.2,147.4,137.6,123.6,117.4,113.7,106.8,103.3,96.3,96.2,94.7,94.6,95.7,106.7,100.2,94.2,97.6,94.3,98,93.6,86.3,90.7,81.8,87.6,73.8,63.3,59.1,52.9,54.6,52.4,67.5,90.4,126,144.3,167.8,166.2,156,137,129.3,118,114.7,112.8,115.7,103.9,96.9,88.8,93,86.3,82.3,82.4,76.6,72.7,67.5,77.3,73.7,73,78.2,90.7,91.5,86.3,86.8,86.1,77.1,75.7,78.7,71.5,69.6,73.6,78.1,78.3,71.5,68.7,61.2,64.7,64.6,56.3,54.5,49.5,54,59.2,52.4,52.8,47.8,45.2,47.1,42.6,42.1,39.4,39.6,37.8,36.8,36.5,34,37.5,36.1,35.1,32.8,32.2,38.1,41.4,38.5,34.6,31.2,34.7,35.8,33.8,32.5,31.2,32.5,32.3,30.1,25.3,24.5,23.7,23.8,26.9,32.4,32.6,31.7,34.1,34.9,33.1,32.3,34.6,32.7,32.6,41,40) > 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 -2.2844634 135.22673 -5.94226288 Feb 1 0.8541024 133.93251 -6.08661725 Mar 1 4.3203676 132.63830 -9.65867106 Apr 1 6.1895201 131.03296 -0.52247702 May 1 5.5041255 129.42761 -1.13173593 Jun 1 3.4693076 127.56939 6.16130016 Jul 1 2.9163069 125.71117 18.77251910 Aug 1 -0.2971697 123.78708 14.11008497 Sep 1 -4.6288280 121.86300 6.36583254 Oct 1 -5.2997081 118.93382 3.76589111 Nov 1 -6.5433159 116.00464 4.23867734 Dec 1 -4.2002438 112.16202 -1.16177770 Jan 2 -2.2844634 108.31940 -2.73494104 Feb 2 0.8541024 105.20403 -9.75812786 Mar 2 4.3203676 102.08865 -10.20901411 Apr 2 6.1895201 100.32584 -11.81536279 May 2 5.5041255 98.56304 -9.46716443 Jun 2 3.4693076 97.86457 -5.63387747 Jul 2 2.9163069 97.16610 6.61759235 Aug 2 -0.2971697 96.82133 3.67584097 Sep 2 -4.6288280 96.47656 2.35227129 Oct 2 -5.2997081 95.60214 7.29756642 Nov 2 -6.5433159 94.72773 6.11558921 Dec 2 -4.2002438 92.55745 9.64278944 Jan 3 -2.2844634 90.38718 5.49728137 Feb 3 0.8541024 86.97545 -1.52954925 Mar 3 4.3203676 83.56371 2.81592070 Apr 3 6.1895201 79.90345 -4.29297031 May 3 5.5041255 76.24319 5.85268572 Jun 3 3.4693076 74.18377 -3.85307718 Jul 3 2.9163069 72.12435 -11.74065723 Aug 3 -0.2971697 73.88230 -14.48512977 Sep 3 -4.6288280 75.64025 -18.11142061 Oct 3 -5.2997081 81.43154 -21.53183524 Nov 3 -6.5433159 87.22284 -28.27952221 Dec 3 -4.2002438 94.88427 -23.18402692 Jan 4 -2.2844634 102.54570 -9.86123994 Feb 4 0.8541024 109.81281 15.33308791 Mar 4 4.3203676 117.07992 22.89971632 Apr 4 6.1895201 122.47604 39.13444231 May 4 5.5041255 127.87216 32.82371535 Jun 4 3.4693076 129.93904 22.59165318 Jul 4 2.9163069 132.00592 2.07777386 Aug 4 -0.2971697 129.51952 0.07765328 Sep 4 -4.6288280 127.03311 -4.40428560 Oct 4 -5.2997081 121.54812 -1.54841129 Nov 4 -6.5433159 116.06313 3.28019067 Dec 4 -4.2002438 110.40016 9.50008156 Jan 5 -2.2844634 104.73720 1.44726415 Feb 5 0.8541024 100.05632 -4.01042540 Mar 5 4.3203676 95.37545 -10.89581438 Apr 5 6.1895201 91.52494 -4.71446458 May 5 5.5041255 87.67444 -6.87856774 Jun 5 3.4693076 84.87099 -6.04029475 Jul 5 2.9163069 82.06753 -2.58383891 Aug 5 -0.2971697 80.76599 -3.86882407 Sep 5 -4.6288280 79.46446 -2.13562753 Oct 5 -5.2997081 79.31757 -6.51786527 Nov 5 -6.5433159 79.17069 4.67262465 Dec 5 -4.2002438 79.48251 -1.58226855 Jan 6 -2.2844634 79.79433 -4.50987005 Feb 6 0.8541024 80.15364 -2.80774449 Mar 6 4.3203676 80.51295 5.86668162 Apr 6 6.1895201 80.74613 4.56435298 May 6 5.5041255 80.97930 -0.18342862 Jun 6 3.4693076 80.86299 2.46770454 Jul 6 2.9163069 80.74667 2.43702056 Aug 6 -0.2971697 80.19389 -2.79671725 Sep 6 -4.6288280 79.64110 0.68772664 Oct 6 -5.2997081 78.39871 5.60099848 Nov 6 -6.5433159 77.15632 0.88699798 Dec 6 -4.2002438 75.48566 -1.68541696 Jan 7 -2.2844634 73.81500 2.06945981 Feb 7 0.8541024 72.05172 5.19418172 Mar 7 4.3203676 70.28843 3.69120420 Apr 7 6.1895201 68.53000 -3.21951601 May 7 5.5041255 66.77156 -3.57568918 Jun 7 3.4693076 65.13165 -7.40096002 Jul 7 2.9163069 63.49174 -1.70804801 Aug 7 -0.2971697 61.85178 3.04538866 Sep 7 -4.6288280 60.21182 0.71700703 Oct 7 -5.2997081 58.48796 1.31174584 Nov 7 -6.5433159 56.76410 -0.72078769 Dec 7 -4.2002438 54.96013 3.24011248 Jan 8 -2.2844634 53.15616 8.32830434 Feb 8 0.8541024 51.37599 0.16990856 Mar 8 4.3203676 49.59582 -1.11618666 Apr 8 6.1895201 48.08554 -6.47505824 May 8 5.5041255 46.57526 -6.87938278 Jun 8 3.4693076 45.25423 -1.62353409 Jul 8 2.9163069 43.93320 -4.24950255 Aug 8 -0.2971697 42.74806 -0.35089032 Sep 8 -4.6288280 41.56292 2.46590360 Oct 8 -5.2997081 40.52869 4.37102033 Nov 8 -6.5433159 39.49445 4.84886472 Dec 8 -4.2002438 38.45039 2.54985741 Jan 9 -2.2844634 37.40632 1.37814179 Feb 9 0.8541024 36.68270 -3.53679887 Mar 9 4.3203676 35.95907 -2.77943897 Apr 9 6.1895201 35.71219 -5.80170556 May 9 5.5041255 35.46530 -5.86942511 Jun 9 3.4693076 35.53145 -6.20075534 Jul 9 2.9163069 35.59760 -6.31390273 Aug 9 -0.2971697 35.73537 2.66180255 Sep 9 -4.6288280 35.87314 10.15568953 Oct 9 -5.2997081 35.78323 8.01647846 Nov 9 -6.5433159 35.69332 5.44999505 Dec 9 -4.2002438 35.20869 0.19155097 Jan 10 -2.2844634 34.72406 2.26039859 Feb 10 0.8541024 33.80762 1.13827884 Mar 10 4.3203676 32.89117 -3.41154035 Apr 10 6.1895201 31.88024 -5.56976280 May 10 5.5041255 30.86931 -5.17343820 Jun 10 3.4693076 30.24274 -1.21204490 Jul 10 2.9163069 29.61616 -0.23246874 Aug 10 -0.2971697 29.39720 0.99996870 Sep 10 -4.6288280 29.17824 0.75058784 Oct 10 -5.2997081 29.19150 0.60820884 Nov 10 -6.5433159 29.20476 1.03855750 Dec 10 -4.2002438 29.27862 -1.27837929 Jan 11 -2.2844634 29.35249 -0.16802439 Feb 11 0.8541024 29.63028 1.91561863 Mar 11 4.3203676 29.90807 -1.62843780 Apr 11 6.1895201 30.78269 -5.27221202 May 11 5.5041255 31.65731 -3.06143921 Jun 11 3.4693076 32.88105 -1.45035490 Jul 11 2.9163069 34.10478 -3.92108775 Aug 11 -0.2971697 35.36209 -2.76491950 Sep 11 -4.6288280 36.61940 2.60943044 Oct 11 -5.2997081 37.95524 0.04447070 Nov 11 -6.5433159 39.29108 -0.14776137 Dec 11 -4.2002438 40.71172 4.48852087 Jan 12 -2.2844634 42.13237 0.15209480 > m$win s t l 1331 19 13 > m$deg s t l 0 1 1 > m$jump s t l 134 2 2 > m$inner [1] 2 > m$outer [1] 0 > postscript(file="/var/wessaorg/rcomp/tmp/14xux1322607689.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/wessaorg/rcomp/tmp/2hp301322607689.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/wessaorg/rcomp/tmp/3s8g51322607690.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/wessaorg/rcomp/tmp/4y5s91322607690.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/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/wessaorg/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/wessaorg/rcomp/tmp/5nfkl1322607690.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/wessaorg/rcomp/tmp/6iti21322607690.tab") > > try(system("convert tmp/14xux1322607689.ps tmp/14xux1322607689.png",intern=TRUE)) character(0) > try(system("convert tmp/2hp301322607689.ps tmp/2hp301322607689.png",intern=TRUE)) character(0) > try(system("convert tmp/3s8g51322607690.ps tmp/3s8g51322607690.png",intern=TRUE)) character(0) > try(system("convert tmp/4y5s91322607690.ps tmp/4y5s91322607690.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 1.635 0.213 1.857