R version 2.13.0 (2011-04-13) Copyright (C) 2011 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i486-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(12008,9169,8788,8417,8247,8197,8236,8253,7733,8366,8626,8863,10102,8463,9114,8563,8872,8301,8301,8278,7736,7973,8268,9476,11100,8962,9173,8738,8459,8078,8411,8291,7810,8616,8312,9692,9911,8915,9452,9112,8472,8230,8384,8625,8221,8649,8625,10443,10357,8586,8892,8329,8101,7922,8120,7838,7735,8406,8209,9451,10041,9411,10405,8467,8464,8102,7627,7513,7510,8291,8064,9383,9706,8579,9474,8318,8213,8059,9111,7708,7680,8014,8007,8718,9486,9113,9025,8476,7952,7759,7835,7600,7651,8319,8812,8630) > 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 179835.6920 210.5086 12427.3516 0.0000 > m$fitted level slope sea Jan 1 12008.000 0.000000 0.000000 Feb 1 9404.590 -140.607836 -235.589921 Mar 1 8928.832 -147.265469 -140.832077 Apr 1 8578.450 -148.949751 -161.450490 May 1 8388.569 -149.310519 -141.569276 Jun 1 8326.478 -148.403971 -129.478110 Jul 1 8356.972 -146.335974 -120.971768 Aug 1 8380.111 -144.205402 -127.110873 Sep 1 7938.235 -148.239371 -205.235083 Oct 1 8392.261 -139.500744 -26.260648 Nov 1 8724.440 -132.218051 -98.440459 Dec 1 8967.931 -126.078791 -104.930566 Jan 2 8669.757 -122.142131 1432.243113 Feb 2 8631.859 -119.587314 -168.859394 Mar 2 9072.837 -104.335789 41.163023 Apr 2 8775.162 -108.208061 -212.162081 May 2 8925.051 -103.402565 -53.050835 Jun 2 8522.557 -109.218612 -221.556873 Jul 2 8406.016 -109.368082 -105.016379 Aug 2 8314.873 -108.981829 -36.872520 Sep 2 8032.749 -112.774523 -296.748699 Oct 2 8006.715 -110.806837 -33.714794 Nov 2 8309.017 -101.202860 -41.016555 Dec 2 9311.454 -77.984631 164.545683 Jan 3 9656.794 -72.898493 1443.206214 Feb 3 9338.590 -80.031055 -376.590458 Mar 3 9131.737 -83.970021 41.262813 Apr 3 9015.485 -84.866042 -277.484929 May 3 8560.304 -94.551105 -101.303992 Jun 3 8324.319 -98.258076 -246.318924 Jul 3 8453.349 -92.181123 -42.348991 Aug 3 8312.941 -93.495583 -21.941400 Sep 3 8129.017 -96.003919 -319.016763 Oct 3 8574.781 -80.748196 41.218785 Nov 3 8587.433 -78.130883 -275.433156 Dec 3 9388.655 -55.047410 303.345486 Jan 4 8660.130 -71.328960 1250.869625 Feb 4 9114.797 -55.186720 -199.797489 Mar 4 9349.946 -45.698071 102.053673 Apr 4 9316.261 -45.320821 -204.261267 May 4 8698.818 -62.639282 -226.817938 Jun 4 8515.290 -66.268037 -285.289606 Jul 4 8394.958 -67.901297 -10.958359 Aug 4 8544.546 -61.266413 80.453531 Sep 4 8609.150 -57.391696 -388.150474 Oct 4 8620.071 -55.279895 28.928756 Nov 4 8967.038 -42.977827 -342.037699 Dec 4 9792.785 -17.257658 650.215124 Jan 5 9395.689 -28.380223 961.311294 Feb 5 8962.447 -41.338044 -376.447035 Mar 5 8801.828 -45.308240 90.171889 Apr 5 8475.595 -54.538434 -146.595437 May 5 8310.132 -58.108959 -209.131892 Jun 5 8198.438 -59.818806 -276.438200 Jul 5 8159.820 -59.141747 -39.819815 Aug 5 7843.454 -67.392731 -5.454465 Sep 5 8049.317 -58.592380 -314.316768 Oct 5 8389.759 -45.743281 16.241038 Nov 5 8676.847 -35.119612 -467.847343 Dec 5 8719.140 -32.685432 731.860152 Jan 6 8928.673 -25.052656 1112.326604 Feb 6 9586.138 -2.706419 -175.137854 Mar 6 10140.504 15.938702 264.496120 Apr 6 8958.064 -24.034918 -491.063902 May 6 8688.098 -32.149780 -224.097821 Jun 6 8419.162 -39.912346 -317.161697 Jul 6 7756.699 -60.299814 -129.699383 Aug 6 7572.690 -64.357427 -59.689582 Sep 6 7808.686 -54.492761 -298.685828 Oct 6 8216.441 -39.332544 74.558978 Nov 6 8482.013 -29.390268 -418.012831 Dec 6 8677.932 -22.092359 705.067791 Jan 7 8756.613 -18.816348 949.386843 Feb 7 8841.567 -15.381508 -262.566929 Mar 7 8930.259 -11.894961 543.741110 Apr 7 8798.988 -15.892883 -480.988434 May 7 8458.043 -26.719011 -245.043380 Jun 7 8260.893 -32.371292 -201.892652 Jul 7 8987.621 -7.232075 123.378711 Aug 7 8106.878 -36.166470 -398.878258 Sep 7 8053.242 -36.745218 -373.241863 Oct 7 8007.172 -37.053683 6.828433 Nov 7 8361.948 -24.137357 -354.948268 Dec 7 8114.182 -31.487774 603.818194 Jan 8 8452.183 -19.307729 1033.816593 Feb 8 9201.233 6.251572 -88.232917 Mar 8 8639.985 -12.750128 385.015228 Apr 8 8818.007 -6.359365 -342.006966 May 8 8330.319 -22.439491 -378.318753 Jun 8 8139.663 -28.044607 -380.662947 Jul 8 7637.863 -43.810832 197.137036 Aug 8 7884.652 -34.142910 -284.652445 Sep 8 7991.588 -29.451000 -340.588491 Oct 8 8297.085 -18.324851 21.915056 Nov 8 8951.117 3.964216 -139.117093 Dec 8 8334.215 -16.593253 295.784817 > m$resid Jan Feb Mar Apr May Jun 1 0.00000000 -3.75533775 -0.78093774 -0.47931765 -0.09646214 0.20538009 2 -0.46822251 0.17628826 1.27989710 -0.45196665 0.60274943 -0.69791337 3 1.03671563 -0.54672463 -0.28803453 -0.07490760 -0.85958474 -0.32802075 4 -1.59569228 1.19501533 0.66010928 0.02774464 -1.32343748 -0.27939022 5 -0.88791618 -0.92617376 -0.27187385 -0.64715023 -0.25611498 -0.12362123 6 0.56267564 1.56576898 1.27263450 -2.75707684 -0.56723834 -0.54581870 7 0.23333763 0.23836908 0.23817523 -0.27447616 -0.74924883 -0.39271677 8 0.85395589 1.76612860 -1.30039384 0.43852086 -1.10898150 -0.38754149 Jul Aug Sep Oct Nov Dec 1 0.42096796 0.39852661 -0.69951401 1.41436816 1.10698799 0.88119482 2 -0.01707513 0.04247783 -0.40337662 0.20202394 0.96181132 2.56517640 3 0.52693277 -0.11177493 -0.20955585 1.25538741 0.21622627 2.03685484 4 -0.12490921 0.50246748 0.29081204 0.15780137 0.92812044 2.00879351 5 0.04889291 -0.59323650 0.63025228 0.92000497 0.76668170 0.17881621 6 -1.43442256 -0.28503904 0.69202719 1.06450044 0.70183267 0.52007525 7 1.74827376 -2.01151435 -0.04022344 -0.02145926 0.90165362 -0.51592112 8 -1.09087317 0.66895659 0.32466920 0.77046682 1.54702066 -1.43175920 > 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/1smr61322575166.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/275vh1322575166.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/3d9201322575166.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/4rx1n1322575166.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/56az41322575166.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/6qcbz1322575166.tab") > > try(system("convert tmp/1smr61322575166.ps tmp/1smr61322575166.png",intern=TRUE)) character(0) > try(system("convert tmp/275vh1322575166.ps tmp/275vh1322575166.png",intern=TRUE)) character(0) > try(system("convert tmp/3d9201322575166.ps tmp/3d9201322575166.png",intern=TRUE)) character(0) > try(system("convert tmp/4rx1n1322575166.ps tmp/4rx1n1322575166.png",intern=TRUE)) character(0) > try(system("convert tmp/56az41322575166.ps tmp/56az41322575166.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 2.028 0.273 2.316