R version 2.15.2 (2012-10-26) -- "Trick or Treat" Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i686-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(88.1,101.7,114.8,103.4,96.4,110,71.1,79.4,119.2,99.1,113.2,103.6,97.5,102.4,120.8,89.5,101.7,112.5,72.4,84.7,117.2,112.8,111.3,102.3,95.2,103,116.4,95.1,100.7,112.4,75.3,93.3,118.6,118.7,110.7,113.3,89.5,106.3,115.1,105.7,95.8,114.7,79.6,80.6,125,127.5,99.5,104.3,90,96,108.9,95.8,87.2,108.4,74.9,80.8,119.1,107.9,106.9,96.8,93.7,95.2,112.7,98.5,91.5,112,76.7,84.7,114.9,108.4,104.6,111.3,90.8,109.1,121,95.2,110.5,102.4,86.7,99.1,126,110.3,104.6,103.1,102) > par8 = 'FALSE' > par7 = '1' > par6 = '' > par5 = '1' > par4 = '' > par3 = '0' > par2 = 'periodic' > par1 = '12' > main = 'Seasonal Decomposition by Loess' > 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 -8.6857782 100.20079 -3.41501073 Feb 1 0.2705805 100.18611 1.24331219 Mar 1 13.8797806 100.17143 0.74879377 Apr 1 -4.2382342 100.21184 7.42639698 May 1 -4.1991128 100.25225 0.34686399 Jun 1 8.4180773 100.35620 1.22571883 Jul 1 -25.2933151 100.46016 -4.06684379 Aug 1 -15.9341316 100.54414 -5.21000859 Sep 1 17.9250551 100.62812 0.64682353 Oct 1 9.9923863 100.64332 -11.53570330 Nov 1 5.1168438 100.65851 7.42464359 Dec 1 2.7478481 100.89242 -0.04027161 Jan 2 -8.6857782 101.12633 5.05944387 Feb 2 0.2705805 101.40025 0.72916958 Mar 2 13.8797806 101.67417 5.24605394 Apr 2 -4.2382342 101.81213 -8.07389895 May 2 -4.1991128 101.95010 3.94901195 Jun 2 8.4180773 101.91290 2.16902636 Jul 2 -25.2933151 101.87569 -4.18237668 Aug 2 -15.9341316 101.84598 -1.21185320 Sep 2 17.9250551 101.81628 -2.54133279 Oct 2 9.9923863 101.88389 0.92372691 Nov 2 5.1168438 101.95150 4.23166032 Dec 2 2.7478481 102.12032 -2.56816760 Jan 3 -8.6857782 102.28914 1.59663514 Feb 3 0.2705805 102.56958 0.15984093 Mar 3 13.8797806 102.85001 -0.32979462 Apr 3 -4.2382342 103.15701 -3.81877210 May 3 -4.1991128 103.46400 1.43511422 Jun 3 8.4180773 103.73794 0.24398571 Jul 3 -25.2933151 104.01188 -3.41856025 Aug 3 -15.9341316 104.22285 5.01127936 Sep 3 17.9250551 104.43383 -3.75888411 Oct 3 9.9923863 104.58478 4.12283341 Nov 3 5.1168438 104.73573 0.84742465 Dec 3 2.7478481 104.74774 5.80441557 Jan 4 -8.6857782 104.75974 -6.57396284 Feb 4 0.2705805 104.67803 1.35138639 Mar 4 13.8797806 104.59633 -3.37610572 Apr 4 -4.2382342 104.49710 5.44112958 May 4 -4.1991128 104.39788 -4.39877132 Jun 4 8.4180773 104.11363 2.16828888 Jul 4 -25.2933151 103.82938 1.06393163 Aug 4 -15.9341316 103.30193 -6.76779872 Sep 4 17.9250551 102.77448 4.30046784 Oct 4 9.9923863 102.06276 15.44485218 Nov 4 5.1168438 101.35105 -6.96788977 Dec 4 2.7478481 100.69982 0.85233677 Jan 5 -8.6857782 100.04858 -1.36280603 Feb 5 0.2705805 99.48621 -3.75679243 Mar 5 13.8797806 98.92384 -3.90362017 Apr 5 -4.2382342 98.51498 1.52325688 May 5 -4.1991128 98.10612 -6.70700228 Jun 5 8.4180773 98.06111 1.92081645 Jul 5 -25.2933151 98.01610 2.17721774 Aug 5 -15.9341316 98.15934 -1.42520402 Sep 5 17.9250551 98.30257 2.87237114 Oct 5 9.9923863 98.50498 -0.59736912 Nov 5 5.1168438 98.70739 3.07576434 Dec 5 2.7478481 98.91811 -4.86595412 Jan 6 -8.6857782 99.12882 3.25695809 Feb 6 0.2705805 99.23770 -4.30828342 Mar 6 13.8797806 99.34659 -0.52636628 Apr 6 -4.2382342 99.45546 3.28277433 May 6 -4.1991128 99.56433 -3.86522127 Jun 6 8.4180773 99.88466 3.69726059 Jul 6 -25.2933151 100.20499 1.78832501 Aug 6 -15.9341316 100.72147 -0.08733548 Sep 6 17.9250551 101.23794 -4.26299906 Oct 6 9.9923863 101.74708 -3.33947095 Nov 6 5.1168438 102.25623 -2.77306913 Dec 6 2.7478481 102.79867 5.75348470 Jan 7 -8.6857782 103.34111 -3.85533080 Feb 7 0.2705805 104.05345 4.77597331 Mar 7 13.8797806 104.76578 2.35443607 Apr 7 -4.2382342 105.02361 -5.58537618 May 7 -4.1991128 105.28144 9.41767537 Jun 7 8.4180773 105.29681 -11.31488547 Jul 7 -25.2933151 105.31218 6.68113623 Aug 7 -15.9341316 105.36829 9.66584043 Sep 7 17.9250551 105.42440 2.65054156 Oct 7 9.9923863 105.44241 -5.13479504 Nov 7 5.1168438 105.46041 -5.97725793 Dec 7 2.7478481 105.40993 -5.05777730 Jan 8 -8.6857782 105.35944 5.32633400 > m$win s t l 851 19 13 > m$deg s t l 0 1 1 > m$jump s t l 86 2 2 > m$inner [1] 2 > m$outer [1] 0 > postscript(file="/var/fisher/rcomp/tmp/19foa1353704160.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/fisher/rcomp/tmp/239ax1353704160.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/fisher/rcomp/tmp/3to2e1353704160.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/fisher/rcomp/tmp/4a9kd1353704160.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/fisher/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/fisher/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/fisher/rcomp/tmp/53gs01353704160.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/fisher/rcomp/tmp/6r7x71353704160.tab") > > try(system("convert tmp/19foa1353704160.ps tmp/19foa1353704160.png",intern=TRUE)) character(0) > try(system("convert tmp/239ax1353704160.ps tmp/239ax1353704160.png",intern=TRUE)) character(0) > try(system("convert tmp/3to2e1353704160.ps tmp/3to2e1353704160.png",intern=TRUE)) character(0) > try(system("convert tmp/4a9kd1353704160.ps tmp/4a9kd1353704160.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 2.812 0.627 3.431