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(563278,571477,583775,598714,581744,530455,527676,536352,545732,558685,561649,540735,538712,547612,563280,581302,572273,518654,520579,530577,540324,547970,555654,551174,548604,563668,586111,604378,600991,544686,537034,551531,563250,574761,580112,575093,557560,564478,580523,596594,586570,536214,523597,536535,536322,532638,528222,516141,501866,506174,517945,533590,528379,477580,469357,490243,492622,507561,516922,514258,509846,527070,541657,564591,555362,498662,511038,525919,531673,548854,560576,557274,565742,587625,619916,625809,619567,572942,572775,574205,579799,590072,593408,597141,595404,612117,628232,628884,620735,569028,567456,573100,584428,589379,590865,595454,594167,611324,612613,610763,593530,542722,536662,543599,555332,560854,562325,554788,547344,565464,577992,579714,569323,506971,500857,509127,509933,517009,519164,512238,509239) > 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 10092094 47490304 15037130 0 > m$fitted level slope sea Jan 1 563278.0 0.00000 0.0000 Feb 1 569468.7 4018.05152 2008.3327 Mar 1 582157.0 9365.42657 1618.0101 Apr 1 597541.0 13261.64103 1173.0323 May 1 588798.4 -771.19503 -7054.3694 Jun 1 542180.6 -29959.55967 -11725.6384 Jul 1 520239.7 -24835.89443 7436.3498 Aug 1 526852.3 -4722.90698 9499.7368 Sep 1 542122.2 8061.14086 3609.7573 Oct 1 557813.7 12939.04403 871.2533 Nov 1 564043.7 8650.46058 -2394.7425 Dec 1 547498.1 -7455.74953 -6763.1412 Jan 2 539216.1 -7977.71306 -504.1031 Feb 2 546456.7 1482.00476 1155.3294 Mar 2 563167.2 10634.86093 112.7765 Apr 2 574364.0 10980.46001 6937.9780 May 2 568652.5 666.56031 3620.4944 Jun 2 538969.3 -17976.63922 -20315.3434 Jul 2 519394.9 -18959.11114 1184.1482 Aug 2 520422.0 -6643.84368 10155.0218 Sep 2 534660.4 6232.11136 5663.6146 Oct 2 545128.4 8843.54990 2841.6083 Nov 2 550622.2 6780.93060 5031.7773 Dec 2 556919.7 6483.81749 -5745.6860 Jan 3 557530.0 2866.13010 -8926.0230 Feb 3 567641.8 7301.10818 -3973.7883 Mar 3 585422.0 13623.33477 688.9943 Apr 3 593677.3 10364.47263 10700.6710 May 3 589805.0 1668.28591 11185.9697 Jun 3 569387.7 -11780.70275 -24701.7392 Jul 3 544912.8 -19501.65462 -7878.8006 Aug 3 542019.6 -9384.11130 9511.4034 Sep 3 553065.1 3074.07290 10184.9181 Oct 3 567726.4 10138.42474 7034.6074 Nov 3 575364.8 8616.06656 4747.2318 Dec 3 580029.9 6210.98743 -4936.9478 Jan 4 574344.1 -1042.04740 -16784.0924 Feb 4 573468.8 -940.76964 -8990.8185 Mar 4 578022.8 2372.99819 2500.1814 Apr 4 580897.8 2676.36838 15696.1959 May 4 570221.5 -5428.42332 16348.4723 Jun 4 556922.9 -10200.89781 -20708.9499 Jul 4 538797.2 -14998.91418 -15200.2277 Aug 4 531356.0 -10420.80462 5179.0472 Sep 4 528428.4 -5878.00402 7893.6401 Oct 4 524714.9 -4565.91379 7923.1278 Nov 4 521338.4 -3845.33363 6883.6324 Dec 4 517189.6 -4029.26114 -1048.5613 Jan 5 517953.1 -1122.62699 -16087.0730 Feb 5 518152.5 -323.08312 -11978.4952 Mar 5 516563.0 -1085.96725 1381.9664 Apr 5 513856.9 -2062.84171 19733.0522 May 5 509335.6 -3549.61686 19043.3849 Jun 5 497778.2 -8392.35447 -20198.2025 Jul 5 487304.8 -9649.13018 -17947.8395 Aug 5 485033.5 -5194.57135 5209.4897 Sep 5 484227.0 -2544.15038 8395.0274 Oct 5 495378.8 5728.66218 12182.1926 Nov 5 508331.7 10092.18752 8590.3428 Dec 5 517662.9 9632.28010 -3404.8821 Jan 6 526823.5 9347.17625 -16977.5012 Feb 6 538481.3 10740.92379 -11411.2973 Mar 6 542532.9 6716.01841 -875.8722 Apr 6 544777.4 4024.19764 19813.6401 May 6 536630.7 -3317.90024 18731.2944 Jun 6 520868.2 -10827.90585 -22206.2357 Jul 6 525158.6 -1714.14623 -14120.6476 Aug 6 523095.1 -1924.61875 2823.8763 Sep 6 526138.9 1068.45142 5534.0569 Oct 6 536331.1 6565.48093 12522.9264 Nov 6 550066.3 10886.68449 10509.7280 Dec 6 561454.8 11189.33911 -4180.8134 Jan 7 581373.8 16453.40246 -15631.8274 Feb 7 598035.9 16579.04922 -10410.9484 Mar 7 618043.2 18639.22322 1872.7652 Apr 7 609704.9 2421.93667 16104.0660 May 7 599341.6 -5275.88057 20225.4132 Jun 7 598022.1 -2892.61707 -25080.1107 Jul 7 588616.7 -6812.50655 -15841.7312 Aug 7 574511.0 -11198.28387 -306.0162 Sep 7 573372.0 -5150.82038 6427.0074 Oct 7 577194.4 244.74048 12877.6357 Nov 7 582281.4 3157.96414 11126.5791 Dec 7 600500.9 12224.41869 -3359.8573 Jan 8 613456.6 12664.54537 -18052.5966 Feb 8 624227.4 11526.42968 -12110.3818 Mar 8 621525.8 2986.82464 6706.1859 Apr 8 612235.6 -4383.61938 16648.3802 May 8 602785.3 -7429.16999 17949.7361 Jun 8 591840.2 -9543.62675 -22812.1861 Jul 8 580578.3 -10576.28855 -13122.3176 Aug 8 574625.7 -7800.07117 -1525.7086 Sep 8 577330.7 -1495.23494 7097.3409 Oct 8 578424.2 58.77881 10954.8226 Nov 8 583920.2 3324.76642 6944.7859 Dec 8 597280.9 9356.18647 -1826.8617 Jan 9 611144.4 12064.25437 -16977.3934 Feb 9 619879.6 10066.51187 -8555.5956 Mar 9 608345.5 -2882.48933 4267.4764 Apr 9 594454.8 -9482.96071 16308.1989 May 9 576211.8 -14740.74371 17318.1834 Jun 9 563445.8 -13554.94422 -20723.7584 Jul 9 550974.7 -12904.40728 -14312.7211 Aug 9 545668.3 -8348.46796 -2069.3441 Sep 9 546261.5 -2989.63095 9070.5329 Oct 9 550399.9 1283.14504 10454.1116 Nov 9 558231.1 5210.79977 4093.8920 Dec 9 560726.6 3581.30512 -5938.5708 Jan 10 563460.0 3072.65797 -16116.0232 Feb 10 566612.7 3120.63281 -1148.7327 Mar 10 569554.7 3013.62953 8437.3331 Apr 10 562994.1 -2720.32222 16719.9445 May 10 553512.8 -6772.83804 15810.1598 Jun 10 531942.7 -15646.53304 -24971.6943 Jul 10 516665.4 -15425.18490 -15808.4343 Aug 10 511292.3 -9405.10037 -2165.3209 Sep 10 503629.6 -8362.15226 6303.3874 Oct 10 505491.8 -2241.08946 11517.2249 Nov 10 512010.5 3006.50069 7153.5193 Dec 10 517842.9 4700.12254 -5604.9023 Jan 11 525399.0 6411.24285 -16160.0029 > m$resid Jan Feb Mar Apr May 1 0.000000000 0.826887946 0.816217836 0.561469521 -2.009620636 2 -0.078694640 1.354027127 1.333570784 0.050751465 -1.487723619 3 -0.529265287 0.640833367 0.917071389 -0.476368203 -1.261133732 4 -1.055206019 0.014660767 0.480616923 0.044202223 -1.177244457 5 0.422132710 0.115828512 -0.110663080 -0.142108592 -0.216023242 6 -0.041376953 0.201992504 -0.583914529 -0.391274842 -1.066740279 7 0.763733435 0.018214162 0.298901385 -2.356269949 -1.118298234 8 0.063846595 -0.165011567 -1.239026149 -1.070599125 -0.442389322 9 0.392821576 -0.289679447 -1.878844047 -0.958590583 -0.763649428 10 -0.073781452 0.006957138 -0.015525988 -0.832641940 -0.588536380 11 0.248208208 Jun Jul Aug Sep Oct 1 -4.233195726 0.743748332 2.918454906 1.855089981 0.707834783 2 -2.695178380 -0.142593690 1.787293494 1.868327431 0.378970774 3 -1.945010573 -1.119588215 1.468358151 1.807802347 1.025107338 4 -0.691102453 -0.695435625 0.664262405 0.659231332 0.190415841 5 -0.701970189 -0.182158151 0.646176824 0.384619680 1.200807052 6 -1.089245155 1.321151190 -0.030526200 0.434337163 0.798025283 7 0.345780937 -0.568330183 -0.636049284 0.877554807 0.783369078 8 -0.306834821 -0.149742961 0.402614449 0.914885278 0.225632387 9 0.172090231 0.094342869 0.660724511 0.777601892 0.620377783 10 -1.287856587 0.032103438 0.873091628 0.151337135 0.888713676 11 Nov Dec 1 -0.622316089 -2.337172958 2 -0.299312322 -0.043224770 3 -0.220998469 -0.349871956 4 0.104636157 -0.026736766 5 0.633708968 -0.066816192 6 0.627582317 0.043953586 7 0.423091810 1.316381532 8 0.474307629 0.875588479 9 0.570371500 -0.236533724 10 0.762008899 0.245827315 11 > 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/wessaorg/rcomp/tmp/1iitk1355677960.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/wessaorg/rcomp/tmp/2z3xy1355677960.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/wessaorg/rcomp/tmp/3hei11355677960.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/wessaorg/rcomp/tmp/4li1z1355677960.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/wessaorg/rcomp/tmp/53fni1355677960.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/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/wessaorg/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/wessaorg/rcomp/tmp/62r9i1355677960.tab") > > try(system("convert tmp/1iitk1355677960.ps tmp/1iitk1355677960.png",intern=TRUE)) character(0) > try(system("convert tmp/2z3xy1355677960.ps tmp/2z3xy1355677960.png",intern=TRUE)) character(0) > try(system("convert tmp/3hei11355677960.ps tmp/3hei11355677960.png",intern=TRUE)) character(0) > try(system("convert tmp/4li1z1355677960.ps tmp/4li1z1355677960.png",intern=TRUE)) character(0) > try(system("convert tmp/53fni1355677960.ps tmp/53fni1355677960.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 3.375 0.583 3.944