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 <- c(2.05,2.11,2.09,2.05,2.08,2.06,2.06,2.08,2.07,2.06,2.07,2.06,2.09,2.07,2.09,2.28,2.33,2.35,2.52,2.63,2.58,2.70,2.81,2.97,3.04,3.28,3.33,3.50,3.56,3.57,3.69,3.82,3.79,3.96,4.06,4.05,4.03,3.94,4.02,3.88,4.02,4.03,4.09,3.99,4.01,4.01,4.19,4.30,4.27,3.82,3.15,2.49,1.81,1.26,1.06,0.84,0.78,0.70,0.36,0.35) > 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.1693670440 2.0483308 -0.167697832 Feb 1 0.1363219782 2.0543799 -0.080701919 Mar 1 0.0472769687 2.0604291 -0.017706063 Apr 1 -0.0275412354 2.0640032 0.013538074 May 1 -0.0863594630 2.0675772 0.098782234 Jun 1 -0.1728276291 2.0704530 0.162374661 Jul 1 -0.1232958371 2.0733287 0.109967130 Aug 1 -0.1036553026 2.0770174 0.106637916 Sep 1 -0.0980148491 2.0807061 0.087308783 Oct 1 -0.0009940517 2.0841023 -0.023108231 Nov 1 0.0680268124 2.0874985 -0.085525311 Dec 1 0.1916954892 2.1107185 -0.242413972 Jan 2 0.1693670440 2.1339385 -0.213305510 Feb 2 0.1363219782 2.1798157 -0.246137670 Mar 2 0.0472769687 2.2256929 -0.182969886 Apr 2 -0.0275412354 2.2901101 0.017431144 May 2 -0.0863594630 2.3545273 0.061832198 Jun 2 -0.1728276291 2.4323983 0.090429320 Jul 2 -0.1232958371 2.5102694 0.133026483 Aug 2 -0.1036553026 2.5974941 0.136161240 Sep 2 -0.0980148491 2.6847188 -0.006703921 Oct 2 -0.0009940517 2.7787600 -0.077765907 Nov 2 0.0680268124 2.8728011 -0.130827959 Dec 2 0.1916954892 2.9734413 -0.195136782 Jan 3 0.1693670440 3.0740814 -0.203448482 Feb 3 0.1363219782 3.1815739 -0.037895898 Mar 3 0.0472769687 3.2890664 -0.006343369 Apr 3 -0.0275412354 3.3977427 0.129798539 May 3 -0.0863594630 3.5064190 0.139940472 Jun 3 -0.1728276291 3.5981064 0.144721230 Jul 3 -0.1232958371 3.6897938 0.123502031 Aug 3 -0.1036553026 3.7511214 0.172533945 Sep 3 -0.0980148491 3.8124489 0.075565940 Oct 3 -0.0009940517 3.8497950 0.111199048 Nov 3 0.0680268124 3.8871411 0.104832090 Dec 3 0.1916954892 3.9156603 -0.057355828 Jan 4 0.1693670440 3.9441796 -0.083546625 Feb 4 0.1363219782 3.9664194 -0.162741411 Mar 4 0.0472769687 3.9886593 -0.015936254 Apr 4 -0.0275412354 4.0066255 -0.099084314 May 4 -0.0863594630 4.0245918 0.081767650 Jun 4 -0.1728276291 4.0391645 0.163663152 Jul 4 -0.1232958371 4.0537371 0.159558697 Aug 4 -0.1036553026 4.0320765 0.061578775 Sep 4 -0.0980148491 4.0104159 0.097598934 Oct 4 -0.0009940517 3.9010766 0.109917478 Nov 4 0.0680268124 3.7917372 0.330235956 Dec 4 0.1916954892 3.5836907 0.524613831 Jan 5 0.1693670440 3.3756441 0.724988829 Feb 5 0.1363219782 3.1012558 0.582422207 Mar 5 0.0472769687 2.8268675 0.275855528 Apr 5 -0.0275412354 2.5038005 0.013740775 May 5 -0.0863594630 2.1807334 -0.284373954 Jun 5 -0.1728276291 1.8551695 -0.422341832 Jul 5 -0.1232958371 1.5296055 -0.346309668 Aug 5 -0.1036553026 1.1989312 -0.255275940 Sep 5 -0.0980148491 0.8682570 0.009757868 Oct 5 -0.0009940517 0.5382334 0.162760618 Nov 5 0.0680268124 0.2082099 0.083763301 Dec 5 0.1916954892 -0.1188396 0.277144102 > m$win s t l 601 19 13 > m$deg s t l 0 1 1 > m$jump s t l 61 2 2 > m$inner [1] 2 > m$outer [1] 0 > postscript(file="/var/www/html/rcomp/tmp/118nf1259959051.ps",horizontal=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/rcomp/tmp/2ssjt1259959051.ps",horizontal=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/rcomp/tmp/3e0wd1259959051.ps",horizontal=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/rcomp/tmp/4vclr1259959052.ps",horizontal=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/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,'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/rcomp/tmp/5d8q41259959052.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/rcomp/tmp/6pydt1259959052.tab") > > system("convert tmp/118nf1259959051.ps tmp/118nf1259959051.png") > system("convert tmp/2ssjt1259959051.ps tmp/2ssjt1259959051.png") > system("convert tmp/3e0wd1259959051.ps tmp/3e0wd1259959051.png") > system("convert tmp/4vclr1259959052.ps tmp/4vclr1259959052.png") > > > proc.time() user system elapsed 0.958 0.603 1.117