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(2360,2214,2825,2355,2333,3016,2155,2172,2150,2533,2058,2160,2260,2498,2695,2799,2947,2930,2318,2540,2570,2669,2450,2842,3440,2678,2981,2260,2844,2546,2456,2295,2379,2479,2057,2280,2351,2276,2548,2311,2201,2725,2408,2139,1898,2537,2069,2063,2524,2437,2189,2793,2074,2622,2278,2144,2427,2139,1828,2072,1800,1758,2246,1987,1868,2514,2121) > 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 75.64979 2362.799 -78.4485919 Feb 1 -66.49498 2363.215 -82.7203355 Mar 1 207.52684 2363.632 253.8413188 Apr 1 50.49558 2365.088 -60.5836452 May 1 16.96441 2366.544 -50.5087079 Jun 1 372.42245 2369.784 273.7932357 Jul 1 -55.95226 2373.024 -162.0720685 Aug 1 -155.32252 2377.340 -50.0171052 Sep 1 -120.39149 2381.655 -111.2634217 Oct 1 74.14538 2393.398 65.4561218 Nov 1 -296.91781 2405.142 -50.2242682 Dec 1 -102.12534 2430.337 -168.2116040 Jan 2 75.64979 2455.532 -271.1815911 Feb 2 -66.49498 2482.778 81.7168429 Mar 2 207.52684 2510.024 -22.5513253 Apr 2 50.49558 2541.903 206.6010465 May 2 16.96441 2573.782 356.2533198 Jun 2 372.42245 2618.286 -60.7081627 Jul 2 -55.95226 2662.789 -288.8368930 Aug 2 -155.32252 2694.931 0.3916654 Sep 2 -120.39149 2727.073 -36.6810561 Oct 2 74.14538 2730.124 -135.2696523 Nov 2 -296.91781 2733.176 13.7418181 Dec 2 -102.12534 2720.758 223.3671341 Jan 3 75.64979 2708.340 656.0097990 Feb 3 -66.49498 2689.692 54.8033510 Mar 3 207.52684 2671.043 102.4303009 Apr 3 50.49558 2635.820 -426.3151131 May 3 16.96441 2600.596 226.4393744 Jun 3 372.42245 2552.474 -378.8964931 Jul 3 -55.95226 2504.352 7.6003915 Aug 3 -155.32252 2464.081 -13.7582745 Sep 3 -120.39149 2423.810 75.5817796 Oct 3 74.14538 2401.897 2.9572115 Nov 3 -296.91781 2379.985 -26.0672902 Dec 3 -102.12534 2366.192 15.9334604 Jan 4 75.64979 2352.399 -77.0484403 Feb 4 -66.49498 2339.386 3.1094326 Mar 4 207.52684 2326.372 14.1007034 Apr 4 50.49558 2317.077 -56.5727321 May 4 16.96441 2307.782 -123.7462662 Jun 4 372.42245 2305.952 46.6256113 Jul 4 -55.95226 2304.122 159.8302411 Aug 4 -155.32252 2307.329 -13.0062231 Sep 4 -120.39149 2310.535 -292.1439670 Oct 4 74.14538 2313.395 149.4594153 Nov 4 -296.91781 2316.255 49.6628642 Dec 4 -102.12534 2315.948 -150.8228431 Jan 5 75.64979 2315.641 132.7087985 Feb 5 -66.49498 2316.842 186.6533047 Mar 5 207.52684 2318.042 -336.5687913 Apr 5 50.49558 2313.250 429.2545365 May 5 16.96441 2308.458 -251.4222342 Jun 5 372.42245 2284.896 -35.3182225 Jul 5 -55.95226 2261.334 72.6185414 Aug 5 -155.32252 2225.198 74.1242817 Sep 5 -120.39149 2189.063 358.3287422 Oct 5 74.14538 2155.143 -90.2885858 Nov 5 -296.91781 2121.224 3.6941528 Dec 5 -102.12534 2098.547 75.5787632 Jan 6 75.64979 2075.869 -351.5192776 Feb 6 -66.49498 2052.659 -228.1638436 Mar 6 207.52684 2029.448 9.0249883 Apr 6 50.49558 2006.628 -70.1233333 May 6 16.96441 1983.807 -132.7717536 Jun 6 372.42245 1963.503 178.0748932 Jul 6 -55.95226 1943.198 233.7542922 > m$win s t l 671 19 13 > m$deg s t l 0 1 1 > m$jump s t l 68 2 2 > m$inner [1] 2 > m$outer [1] 0 > postscript(file="/var/www/html/rcomp/tmp/16gk01260385744.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/2rv5m1260385744.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/3r9081260385744.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/4y9d91260385744.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/51uyh1260385744.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/6bo6j1260385744.tab") > system("convert tmp/16gk01260385744.ps tmp/16gk01260385744.png") > system("convert tmp/2rv5m1260385744.ps tmp/2rv5m1260385744.png") > system("convert tmp/3r9081260385744.ps tmp/3r9081260385744.png") > system("convert tmp/4y9d91260385744.ps tmp/4y9d91260385744.png") > > > proc.time() user system elapsed 0.960 0.604 1.824