R version 2.10.1 (2009-12-14)
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(31245,30951,30872,30752,30967,30781,30681,31356,31434,31594,31949,32396,32441,32447,32288,32418,32346,32091,31855,31683,31615,31840,31536,31383,31638,31626,31720,31472,31372,31419,31341,31171,31036,30532,30666,30571,30173,30032,29874,30018,29911,29963,30050,29901,29544,29451,29293,29334,29389,29563)
> 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  158.923126 30574.51  511.562048
Feb 1  119.080612 30690.29  141.625109
Mar 1   42.212320 30806.07   23.713947
Apr 1   35.852488 30927.18 -211.029548
May 1   36.992627 31048.28 -118.273012
Jun 1  -24.242738 31170.46 -365.219210
Jul 1  -81.728152 31292.64 -529.915360
Aug 1   -4.365376 31413.85  -53.481415
Sep 1  -93.502399 31535.05   -7.547671
Oct 1 -117.952224 31670.51   41.439126
Nov 1  -82.652054 31805.98  225.675929
Dec 1   11.381570 31913.02  471.602402
Jan 2  158.923126 32020.06  262.020945
Feb 2  119.080612 32067.42  260.503881
Mar 2   42.212320 32114.78  131.012595
Apr 2   35.852488 32100.51  281.640713
May 2   36.992627 32086.24  222.768860
Jun 2  -24.242738 32027.76   87.486514
Jul 2  -81.728152 31969.27  -32.545784
Aug 2   -4.365376 31896.84 -209.474480
Sep 2  -93.502399 31824.41 -115.903376
Oct 2 -117.952224 31755.65  202.301759
Nov 2  -82.652054 31686.90  -68.243100
Dec 2   11.381570 31631.09 -259.472254
Jan 3  158.923126 31575.29  -96.209339
Feb 3  119.080612 31524.20  -17.282179
Mar 3   42.212320 31473.12  204.670757
Apr 3   35.852488 31407.30   28.849289
May 3   36.992627 31341.48   -6.472149
Jun 3  -24.242738 31248.96  194.283865
Jul 3  -81.728152 31156.44  266.289927
Aug 3   -4.365376 31030.01  145.357991
Sep 3  -93.502399 30903.58  225.925854
Oct 3 -117.952224 30762.63 -112.676249
Nov 3  -82.652054 30621.68  126.971652
Dec 3   11.381570 30493.33   66.290007
Jan 4  158.923126 30364.98 -350.899569
Feb 4  119.080612 30254.03 -341.114962
Mar 4   42.212320 30143.09 -311.304579
Apr 4   35.852488 30047.70  -65.554797
May 4   36.992627 29952.31  -78.304985
Jun 4  -24.242738 29874.28  112.963171
Jul 4  -81.728152 29796.25  335.481376
Aug 4   -4.365376 29724.98  180.387392
Sep 4  -93.502399 29653.71  -16.206793
Oct 4 -117.952224 29582.87  -13.918271
Nov 4  -82.652054 29512.03 -136.379745
Dec 4   11.381570 29439.24 -116.619141
Jan 5  158.923126 29366.44 -136.366468
Feb 5  119.080612 29292.20  151.718341
> m$win
  s   t   l 
501  19  13 
> m$deg
s t l 
0 1 1 
> m$jump
 s  t  l 
51  2  2 
> m$inner
[1] 2
> m$outer
[1] 0
> postscript(file="/var/yougetitorg/rcomp/tmp/1p7ft1293638423.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/yougetitorg/rcomp/tmp/2p7ft1293638423.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/yougetitorg/rcomp/tmp/3igwe1293638423.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/yougetitorg/rcomp/tmp/4speh1293638423.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/yougetitorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/yougetitorg/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/yougetitorg/rcomp/tmp/5ohbq1293638423.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/yougetitorg/rcomp/tmp/6azsw1293638423.tab") 
> 
> try(system("convert tmp/1p7ft1293638423.ps tmp/1p7ft1293638423.png",intern=TRUE))
character(0)
> try(system("convert tmp/2p7ft1293638423.ps tmp/2p7ft1293638423.png",intern=TRUE))
character(0)
> try(system("convert tmp/3igwe1293638423.ps tmp/3igwe1293638423.png",intern=TRUE))
character(0)
> try(system("convert tmp/4speh1293638423.ps tmp/4speh1293638423.png",intern=TRUE))
character(0)
> 
> 
> proc.time()
   user  system elapsed 
  1.220   0.800   1.628