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(-820.8,993.3,741.7,603.6,-145.8,-35.1,395.1,523.1,462.3,183.4,791.5,344.8,-217.0,406.7,228.6,-580.1,-1550.4,-1447.5,-40.1,-1033.5,-925.6,-347.8,-447.7,-102.6,-2062.2,-929.7,-720.7,-1541.8,-1432.3,-1216.2,-212.8,-378.2,76.9,-101.3,220.4,495.6,-1035.2,61.8,-734.8,-6.9,-1061.1,-854.6,-186.5,244.0,-992.6,-335.2,316.8,477.6,-572.1,1115.2) > 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 -778.04136 515.47858 -558.237218 Feb 1 501.21604 489.67782 2.406139 Mar 1 127.34117 463.87706 150.481774 Apr 1 -115.42844 435.40277 283.625668 May 1 -764.29737 406.92848 211.568888 Jun 1 -603.45075 377.47169 190.879054 Jul 1 275.62115 348.01491 -228.536059 Aug 1 126.49857 314.13604 82.465388 Sep 1 -56.14825 280.25718 238.191071 Oct 1 146.05976 212.40524 -175.064999 Nov 1 524.21755 144.55329 122.729158 Dec 1 616.41203 50.71503 -322.327062 Jan 2 -778.04136 -43.12323 604.164590 Feb 2 501.21604 -141.78087 47.264825 Mar 2 127.34117 -240.43851 341.697338 Apr 2 -115.42844 -334.02792 -130.643646 May 2 -764.29737 -427.61733 -358.485304 Jun 2 -603.45075 -521.67481 -322.374446 Jul 2 275.62115 -615.73228 300.011133 Aug 2 126.49857 -706.16408 -453.834486 Sep 2 -56.14825 -796.59588 -72.855868 Oct 2 146.05976 -850.94268 357.082919 Nov 2 524.21755 -905.28948 -66.628069 Dec 2 616.41203 -912.25886 193.246828 Jan 3 -778.04136 -919.22824 -364.930401 Feb 3 501.21604 -888.31578 -542.600262 Mar 3 127.34117 -857.40333 9.362155 Apr 3 -115.42844 -798.80074 -627.570826 May 3 -764.29737 -740.19815 72.195518 Jun 3 -603.45075 -664.28215 51.532893 Jul 3 275.62115 -588.36614 99.944989 Aug 3 126.49857 -516.69306 11.994491 Sep 3 -56.14825 -445.01998 578.068229 Oct 3 146.05976 -393.67668 146.316920 Nov 3 524.21755 -342.33339 38.515837 Dec 3 616.41203 -318.47215 197.660114 Jan 4 -778.04136 -294.61091 37.452265 Feb 4 501.21604 -298.99693 -140.419109 Mar 4 127.34117 -303.38296 -558.758206 Apr 4 -115.42844 -318.87826 427.406695 May 4 -764.29737 -334.37355 37.570920 Jun 4 -603.45075 -302.98286 51.833604 Jul 4 275.62115 -271.59216 -190.528993 Aug 4 126.49857 -225.51934 343.020774 Sep 4 -56.14825 -179.44653 -757.005222 Oct 4 146.05976 -132.68276 -348.577000 Nov 4 524.21755 -85.91900 -121.498552 Dec 4 616.41203 -36.28091 -102.531123 Jan 5 -778.04136 13.35718 192.584178 Feb 5 501.21604 67.53882 546.445133 > 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/www/html/rcomp/tmp/1lzm91291927104.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/rcomp/tmp/2lzm91291927104.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/rcomp/tmp/3w8mu1291927104.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/rcomp/tmp/4w8mu1291927104.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/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/52amg1291927105.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/66a341291927105.tab") > > try(system("convert tmp/1lzm91291927104.ps tmp/1lzm91291927104.png",intern=TRUE)) character(0) > try(system("convert tmp/2lzm91291927104.ps tmp/2lzm91291927104.png",intern=TRUE)) character(0) > try(system("convert tmp/3w8mu1291927104.ps tmp/3w8mu1291927104.png",intern=TRUE)) character(0) > try(system("convert tmp/4w8mu1291927104.ps tmp/4w8mu1291927104.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 0.916 0.661 2.341