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(4716.99,4926.65,4920.10,5170.09,5246.24,5283.61,4979.05,4825.20,4695.12,4711.54,4727.22,4384.96,4378.75,4472.93,4564.07,4310.54,4171.38,4049.38,3591.37,3720.46,4107.23,4101.71,4162.34,4136.22,4125.88,4031.48,3761.36,3408.56,3228.47,3090.45,2741.14,2980.44,3104.33,3181.57,2863.86,2898.01,3112.33,3254.33,3513.47,3587.61,3727.45,3793.34,3817.58,3845.13,3931.86,4197.52,4307.13,4229.43,4362.28,4217.34,4361.28,4327.74,4417.65,4557.68,4650.35,4967.18,5123.42,5290.85,5535.66,5514.06,5493.88,5694.83,5850.41,6116.64,6175.00,6513.58,6383.78,6673.66,6936.61,7300.68,7392.93,7497.31,7584.71,7160.79,7196.19,7245.63,7347.51,7425.75,7778.51,7822.33,8181.22,8371.47,8347.71,8672.11,8802.79,9138.46,9123.29,9023.21,8850.41,8864.58,9163.74,8516.66,8553.44,7555.20,7851.22,7442.00,7992.53,8264.04,7517.39,7200.40,7193.69,6193.58,5104.21,4800.46,4461.61,4398.59,4243.63,4293.82) > 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 55.428865 5040.131 -378.569422 Feb 1 125.204255 5014.603 -213.157009 Mar 1 90.091139 4989.075 -159.066091 Apr 1 53.826094 4957.067 159.196531 May 1 60.285475 4925.060 260.894728 Jun 1 1.376052 4890.268 391.965624 Jul 1 -166.018898 4855.477 289.592047 Aug 1 -165.566276 4819.976 170.789842 Sep 1 -53.832603 4784.476 -35.523412 Oct 1 -33.845172 4718.555 26.830076 Nov 1 20.395447 4652.634 54.190376 Dec 1 12.655447 4558.902 -186.597406 Jan 2 55.428865 4465.170 -141.848606 Feb 2 125.204255 4381.473 -33.746768 Mar 2 90.091139 4297.775 176.203576 Apr 2 53.826094 4242.506 14.207945 May 2 60.285475 4187.237 -76.142111 Jun 2 1.376052 4152.639 -104.635550 Jul 2 -166.018898 4118.042 -360.653460 Aug 2 -165.566276 4081.197 -195.170532 Sep 2 -53.832603 4044.351 116.711346 Oct 2 -33.845172 3988.595 146.960057 Nov 2 20.395447 3932.839 209.105580 Dec 2 12.655447 3859.262 264.302257 Jan 3 55.428865 3785.686 284.765515 Feb 3 125.204255 3698.779 207.496596 Mar 3 90.091139 3611.873 59.396183 Apr 3 53.826094 3513.843 -159.109463 May 3 60.285475 3415.814 -247.629536 Jun 3 1.376052 3324.432 -235.358165 Jul 3 -166.018898 3233.050 -325.891265 Aug 3 -165.566276 3180.911 -34.904239 Sep 3 -53.832603 3128.771 29.391737 Oct 3 -33.845172 3136.767 78.647802 Nov 3 20.395447 3144.764 -301.299322 Dec 3 12.655447 3200.206 -314.851796 Jan 4 55.428865 3255.649 -198.747688 Feb 4 125.204255 3336.594 -207.468045 Mar 4 90.091139 3417.539 5.840104 Apr 4 53.826094 3517.852 15.932119 May 4 60.285475 3618.165 48.999708 Jun 4 1.376052 3725.192 66.771973 Jul 4 -166.018898 3832.219 151.379766 Aug 4 -165.566276 3918.659 92.037215 Sep 4 -53.832603 4005.099 -19.406386 Oct 4 -33.845172 4068.177 163.188535 Nov 4 20.395447 4131.254 155.480268 Dec 4 12.655447 4191.728 25.046241 Jan 5 55.428865 4252.202 54.648796 Feb 5 125.204255 4328.930 -236.794455 Mar 5 90.091139 4405.658 -134.469201 Apr 5 53.826094 4503.174 -229.260324 May 5 60.285475 4600.690 -243.325873 Jun 5 1.376052 4713.236 -156.932051 Jul 5 -166.018898 4825.782 -9.412702 Aug 5 -165.566276 4951.600 181.146119 Sep 5 -53.832603 5077.419 99.833888 Oct 5 -33.845172 5214.111 110.584550 Nov 5 20.395447 5350.803 164.462022 Dec 5 12.655447 5492.890 8.514786 Jan 6 55.428865 5634.977 -196.525868 Feb 6 125.204255 5781.460 -211.834363 Mar 6 90.091139 5927.943 -167.624352 Apr 6 53.826094 6088.872 -26.058324 May 6 60.285475 6249.801 -135.086723 Jun 6 1.376052 6418.016 94.188033 Jul 6 -166.018898 6586.231 -36.431683 Aug 6 -165.566276 6730.225 109.001551 Sep 6 -53.832603 6874.219 116.223735 Oct 6 -33.845172 6978.276 356.248673 Nov 6 20.395447 7082.334 290.200422 Dec 6 12.655447 7165.837 318.817178 Jan 7 55.428865 7249.341 279.940515 Feb 7 125.204255 7332.116 -296.530028 Mar 7 90.091139 7414.891 -308.792066 Apr 7 53.826094 7504.119 -312.314636 May 7 60.285475 7593.346 -306.121632 Jun 7 1.376052 7706.167 -281.792706 Jul 7 -166.018898 7818.987 125.541749 Aug 7 -165.566276 7965.427 22.468900 Sep 7 -53.832603 8111.868 123.185001 Oct 7 -33.845172 8259.723 145.592504 Nov 7 20.395447 8407.578 -80.263180 Dec 7 12.655447 8524.308 135.146368 Jan 8 55.428865 8641.039 106.322497 Feb 8 125.204255 8702.473 310.782606 Mar 8 90.091139 8763.908 269.291219 Apr 8 53.826094 8747.773 221.610728 May 8 60.285475 8731.639 58.485811 Jun 8 1.376052 8652.791 210.412607 Jul 8 -166.018898 8573.944 755.814931 Aug 8 -165.566276 8460.995 221.231465 Sep 8 -53.832603 8348.046 259.226948 Oct 8 -33.845172 8197.904 -608.858540 Nov 8 20.395447 8047.762 -216.937217 Dec 8 12.655447 7829.067 -399.722458 Jan 9 55.428865 7610.372 326.728882 Feb 9 125.204255 7326.887 811.948552 Mar 9 90.091139 7043.402 383.896728 Apr 9 53.826094 6700.000 446.573584 May 9 60.285475 6356.599 776.806014 Jun 9 1.376052 6020.747 171.457084 Jul 9 -166.018898 5684.895 -414.666318 Aug 9 -165.566276 5337.497 -371.470774 Sep 9 -53.832603 4990.099 -474.656280 Oct 9 -33.845172 4631.680 -199.244613 Nov 9 20.395447 4273.261 -50.026133 Dec 9 12.655447 3910.210 370.955052 > m$win s t l 1081 19 13 > m$deg s t l 0 1 1 > m$jump s t l 109 2 2 > m$inner [1] 2 > m$outer [1] 0 > postscript(file="/var/www/html/rcomp/tmp/1pe9b1259862884.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/2y7sl1259862884.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/3iyvn1259862884.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/4ngy61259862884.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/5tvrs1259862884.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/6yf0p1259862884.tab") > system("convert tmp/1pe9b1259862884.ps tmp/1pe9b1259862884.png") > system("convert tmp/2y7sl1259862884.ps tmp/2y7sl1259862884.png") > system("convert tmp/3iyvn1259862884.ps tmp/3iyvn1259862884.png") > system("convert tmp/4ngy61259862884.ps tmp/4ngy61259862884.png") > > > proc.time() user system elapsed 1.192 0.631 1.400