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(867.887509505211 + ,-2250.28069676838 + ,33618.3570412959 + ,9954.34468238836 + ,354.191730842355 + ,18882.406400463 + ,20229.4310915672 + ,268402.416151187 + ,-113346.926055862 + ,-45016.394227939 + ,35069.861367254 + ,58531.0957290091 + ,-77256.3771198791 + ,-31473.594568955 + ,-52391.0075132882 + ,32854.9847569661 + ,101107.732845397 + ,-176275.960398033 + ,79531.884415102 + ,-176414.251561376 + ,151290.579462589 + ,167731.594163443 + ,143237.122434691 + ,80251.9665265577 + ,118735.726273623 + ,75035.8259037494 + ,19198.3085437346 + ,-36364.5639276314 + ,-36170.5787440905 + ,-109567.395064155 + ,-100783.336097857 + ,-149267.403931369 + ,38947.3510583149 + ,58613.0600994635 + ,16074.4602044407 + ,-41563.0049150659 + ,-15970.5964777959 + ,-47563.9548420802 + ,59595.3577179922 + ,65897.8405390448 + ,-166489.283203891 + ,46312.3269884632 + ,-15952.8722863516 + ,-87780.6523566012 + ,134744.172737777 + ,75232.8122408289 + ,24408.7558471514 + ,-15406.1403381955 + ,-3766.75348767364 + ,27197.2239951006 + ,-46777.2890031503 + ,-82472.8212609495 + ,-35154.7184801844 + ,-46946.870011609 + ,-43641.5364941684 + ,-54920.7084991732 + ,54905.4038222157 + ,-10509.5840707596 + ,-13706.8046976985 + ,-42347.6087628011 + ,-28990.4687708701) > par1 = '12' > par1 <- '12' > #'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) > nx <- length(x) > x <- ts(x,frequency=par1) > m <- StructTS(x,type='BSM') > m$coef level slope seas epsilon 0 0 0 7166614985 > m$fitted level slope sea Jan 1 867.8875 0.000000 0.000000 Feb 1 482.1351 -128.584124 -128.584124 Mar 1 5827.8728 1239.996323 1239.996323 Apr 1 7380.4916 1302.520811 1302.520811 May 1 6728.8281 976.823431 976.823431 Jun 1 9813.3857 1277.928308 1277.928308 Jul 1 12699.0794 1478.898987 1478.898987 Aug 1 64503.4064 7070.613215 7070.613215 Sep 1 34717.0078 3384.912026 3384.912026 Oct 1 22176.6419 1937.159490 1937.159490 Nov 1 25701.7792 2069.490972 2069.490972 Dec 1 32597.1855 2440.715231 2440.715231 Jan 2 31748.4884 2205.757203 -24263.329242 Feb 2 21903.7169 1402.388623 1402.388624 Mar 2 10487.6072 601.232482 601.232482 Apr 2 14393.1575 795.604118 795.604118 May 2 27731.7131 1492.434754 1492.434754 Jun 2 316.7008 -29.009829 -29.009829 Jul 2 10817.3381 497.472530 497.472530 Aug 2 -12536.1649 -638.288208 -638.288208 Sep 2 6835.6735 271.262998 271.262998 Oct 2 25737.7662 1081.299073 1081.299073 Nov 2 39691.4711 1617.649315 1617.649315 Dec 2 45317.9402 1778.002107 1778.002107 Jan 3 52762.2071 1995.935369 -21955.289055 Feb 3 56715.5142 2068.430616 2068.430609 Mar 3 54503.6231 1915.561985 1915.561985 Apr 3 47058.9943 1592.796784 1592.796784 May 3 40420.0893 1318.406726 1318.406726 Jun 3 27702.7769 865.641594 865.641594 Jul 3 16989.6404 503.804777 503.804777 Aug 3 3092.6670 67.417562 67.417562 Sep 3 6141.5407 155.107448 155.107448 Oct 3 10522.4080 275.843444 275.843444 Nov 3 11191.8507 286.776756 286.776756 Dec 3 7396.0977 176.438112 176.438112 Jan 4 6411.8145 145.892815 -1604.820967 Feb 4 2466.9904 41.002639 41.002636 Mar 4 6689.0814 145.529848 145.529848 Apr 4 11037.4397 248.037859 248.037859 May 4 -1085.4937 -46.509077 -46.509077 Jun 4 2085.6397 28.319818 28.319818 Jul 4 917.6857 1.131777 1.131777 Aug 4 -4807.1928 -126.112895 -126.112895 Sep 4 3885.2316 65.594436 65.594436 Oct 4 8343.0778 159.046603 159.046603 Nov 4 9452.3122 178.842181 178.842181 Dec 4 8141.7731 148.446645 148.446645 Jan 5 7817.5117 138.992485 -1528.917330 Feb 5 9064.8683 160.725114 160.725113 Mar 5 6035.0319 99.368161 99.368161 Apr 5 1198.0155 6.228830 6.228830 May 5 -779.1295 -30.500314 -30.500314 Jun 5 -3275.4449 -75.333316 -75.333316 Jul 5 -5460.6999 -113.010490 -113.010490 Aug 5 -8108.8727 -157.487022 -157.487022 Sep 5 -5064.5766 -102.283863 -102.283863 Oct 5 -5427.0997 -106.694699 -106.694699 Nov 5 -5927.3991 -113.254777 -113.254777 Dec 5 -7776.6662 -141.713996 -141.713996 Jan 6 -8796.9335 -155.884210 1714.726302 > m$resid Jan Feb Mar Apr May Jun 1 0.00000000 -0.03365884 0.35173005 0.01709034 -0.09939697 0.10530301 2 -1.00516658 -0.71900851 -0.82639311 0.22841418 0.92402182 -2.25824452 3 1.05779316 0.20361213 -0.46511995 -1.05990664 -0.96924144 -1.71517306 4 -0.25025956 -0.61575838 0.64805581 0.67005637 -2.02671385 0.54119154 5 -0.12093722 0.21886334 -0.64394129 -1.01771020 -0.41746323 -0.52960119 6 -0.26302573 Jul Aug Sep Oct Nov Dec 1 0.08146842 2.63468434 -2.01426689 -0.91346731 0.09583834 0.30667743 2 0.86866577 -2.07060096 1.82254754 1.77564788 1.28077720 0.41555187 3 -1.46597720 -1.88615204 0.40340756 0.58992717 0.05662731 -0.60467240 4 -0.20639699 -1.01261778 1.59733112 0.81432121 0.18019881 -0.28875071 5 -0.46218765 -0.56613699 0.72858821 -0.06031951 -0.09289006 -0.41698756 6 > mylevel <- as.numeric(m$fitted[,'level']) > myslope <- as.numeric(m$fitted[,'slope']) > myseas <- as.numeric(m$fitted[,'sea']) > myresid <- as.numeric(m$resid) > myfit <- mylevel+myseas > mylagmax <- nx/2 > postscript(file="/var/fisher/rcomp/tmp/1gygf1353966046.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(mylevel,na.action=na.pass,lag.max = mylagmax,main='Level') > acf(myseas,na.action=na.pass,lag.max = mylagmax,main='Seasonal') > acf(myresid,na.action=na.pass,lag.max = mylagmax,main='Standardized Residals') > par(op) > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/2qfov1353966046.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(mylevel,main='Level') > spectrum(myseas,main='Seasonal') > spectrum(myresid,main='Standardized Residals') > par(op) > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/3u1w21353966046.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(mylevel,main='Level') > cpgram(myseas,main='Seasonal') > cpgram(myresid,main='Standardized Residals') > par(op) > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/4r0c51353966046.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > plot(as.numeric(m$resid),main='Standardized Residuals',ylab='Residuals',xlab='time',type='b') > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/5i7vj1353966046.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > op <- par(mfrow = c(2,2)) > hist(m$resid,main='Residual Histogram') > plot(density(m$resid),main='Residual Kernel Density') > qqnorm(m$resid,main='Residual Normal QQ Plot') > qqline(m$resid) > plot(m$resid^2, myfit^2,main='Sq.Resid vs. Sq.Fit',xlab='Squared residuals',ylab='Squared Fit') > 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,'Structural Time Series Model',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,'Level',header=TRUE) > a<-table.element(a,'Slope',header=TRUE) > a<-table.element(a,'Seasonal',header=TRUE) > a<-table.element(a,'Stand. Residuals',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,mylevel[i]) + a<-table.element(a,myslope[i]) + a<-table.element(a,myseas[i]) + a<-table.element(a,myresid[i]) + a<-table.row.end(a) + } > a<-table.end(a) > table.save(a,file="/var/fisher/rcomp/tmp/6xlaq1353966046.tab") > > try(system("convert tmp/1gygf1353966046.ps tmp/1gygf1353966046.png",intern=TRUE)) character(0) > try(system("convert tmp/2qfov1353966046.ps tmp/2qfov1353966046.png",intern=TRUE)) character(0) > try(system("convert tmp/3u1w21353966046.ps tmp/3u1w21353966046.png",intern=TRUE)) character(0) > try(system("convert tmp/4r0c51353966046.ps tmp/4r0c51353966046.png",intern=TRUE)) character(0) > try(system("convert tmp/5i7vj1353966046.ps tmp/5i7vj1353966046.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 3.130 0.696 3.810