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(116,111,104,100,93,91,119,139,134,124,113,109,109,106,101,98,93,91,122,139,140,132,117,114,113,110,107,103,98,98,137,148,147,139,130,128,127,123,118,114,108,111,151,159,158,148,138,137,136,133,126,120,114,116,153,162,161,149,139,135,130,127,122,117,112,113,149,157,157,147,137,132,125,123,117,114,111,112,144,150,149,134,123,116,117,111,105,102,95,93,124,130,124,115,106,105,105,101,95,93,84,87,116,120,117,109,105,107,109,109,108,107,99,103,131,137,135,124,118,121,121,118,113,107,100,102,130,136,133,120,112,109,110,106,102,98,92,92,120,127,124,114,108,106,111,110,104,100,96,98,122,134,133) > par1 = '12' > par1 <- as.numeric(par1) > nx <- length(x) > x <- ts(x,frequency=par1) > m <- StructTS(x,type='BSM') > m$coef level slope seas epsilon 0.000000 25.917387 2.629528 0.000000 > m$fitted level slope sea Jan 1 116.00000 0.00000000 0.0000000 Feb 1 111.66898 -4.38630593 -0.6689818 Mar 1 104.19855 -6.77087578 -0.1985527 Apr 1 99.44268 -5.15173394 0.5573243 May 1 93.29390 -5.93831711 -0.2938962 Jun 1 90.38557 -3.54380151 0.6144333 Jul 1 114.10138 18.02455377 4.8986233 Aug 1 139.15291 23.58449015 -0.1529128 Sep 1 138.81447 4.65877045 -4.8144658 Oct 1 125.96551 -9.19102185 -1.9655085 Nov 1 112.79785 -12.33684787 0.2021470 Dec 1 107.48147 -6.78304180 1.5185315 Jan 2 108.00291 -1.01439390 0.9970930 Feb 2 106.59515 -1.32595552 -0.5951504 Mar 2 102.26341 -3.60267046 -1.2634147 Apr 2 96.30632 -5.40215311 1.6936809 May 2 91.69739 -4.79764784 1.3026121 Jun 2 94.86368 1.24820505 -3.8636820 Jul 2 116.67799 16.88365129 5.3220137 Aug 2 135.75827 18.55445581 3.2417292 Sep 2 143.61466 10.41737645 -3.6146636 Oct 2 135.91187 -3.36443324 -3.9118699 Nov 2 119.77348 -13.07897123 -2.7734771 Dec 2 112.37216 -8.76393616 1.6278407 Jan 3 110.75805 -3.33109918 2.2419501 Feb 3 109.69862 -1.60281151 0.3013767 Mar 3 107.36251 -2.15363758 -0.3625113 Apr 3 100.93628 -5.36760575 2.0637193 May 3 97.15912 -4.16896632 0.8408821 Jun 3 104.95019 4.81353492 -6.9501920 Jul 3 129.95290 19.97942060 7.0471024 Aug 3 145.15268 16.38619248 2.8473205 Sep 3 149.03533 6.98483559 -2.0353334 Oct 3 141.94164 -3.60023141 -2.9416410 Nov 3 134.23899 -6.68363957 -4.2389863 Dec 3 127.88396 -6.43680687 0.1160351 Jan 4 124.77213 -3.93796173 2.2278656 Feb 4 122.09890 -2.98784833 0.9011005 Mar 4 116.69818 -4.78988436 1.3018190 Apr 4 111.01451 -5.45760085 2.9854878 May 4 109.11639 -2.79294361 -1.1163862 Jun 4 121.10614 8.25113655 -10.1061389 Jul 4 141.81836 17.55555578 9.1816416 Aug 4 154.90582 14.21704863 4.0941778 Sep 4 158.67662 6.40946026 -0.6766181 Oct 4 152.25032 -3.18356436 -4.2503228 Nov 4 143.24358 -7.53391614 -5.2435830 Dec 4 137.42438 -6.25335275 -0.4243775 Jan 5 133.59705 -4.44050001 2.4029467 Feb 5 130.60393 -3.35980697 2.3960670 Mar 5 123.90004 -5.84824865 2.0999645 Apr 5 116.98703 -6.64064517 3.0129742 May 5 117.48114 -1.32235724 -3.4811379 Jun 5 127.98206 7.48161453 -11.9820637 Jul 5 142.54885 12.75345755 10.4511514 Aug 5 155.88738 13.18897874 6.1126181 Sep 5 160.27838 6.63766978 0.7216182 Oct 5 154.11005 -2.89785525 -5.1100480 Nov 5 145.67068 -7.02281663 -6.6706758 Dec 5 136.61344 -8.53702713 -1.6134392 Jan 6 127.89013 -8.67573968 2.1098721 Feb 6 122.37680 -6.32285724 4.6232043 Mar 6 118.16354 -4.75671277 3.8364641 Apr 6 114.76706 -3.74684272 2.2329409 May 6 117.22761 0.86728478 -5.2276064 Jun 6 126.03217 6.76440199 -13.0321701 Jul 6 138.66502 11.12134169 10.3349771 Aug 6 149.47235 10.88816709 7.5276461 Sep 6 154.46098 6.50601678 2.5390214 Oct 6 152.40248 0.14437904 -5.4024833 Nov 6 144.55064 -5.79388247 -7.5506410 Dec 6 134.21488 -9.16699186 -2.2148837 Jan 7 123.64291 -10.21074078 1.3570928 Feb 7 117.48147 -7.20481467 5.5185301 Mar 7 112.35925 -5.66118114 4.6407491 Apr 7 111.76918 -1.90233375 2.2308152 May 7 116.89909 3.31500147 -5.8990859 Jun 7 125.98271 7.59401549 -13.9827083 Jul 7 134.27541 8.11197011 9.7245884 Aug 7 141.81249 7.68578845 8.1875121 Sep 7 145.55699 4.76338218 3.4430086 Oct 7 139.58623 -3.19596488 -5.5862345 Nov 7 130.27949 -7.72662749 -7.2794927 Dec 7 118.35404 -10.84012076 -2.3540401 Jan 8 114.73889 -5.48174032 2.2611081 Feb 8 106.10047 -7.82132845 4.8995337 Mar 8 100.50883 -6.17061257 4.4911687 Apr 8 100.25004 -1.79379766 1.7499559 May 8 101.91350 0.76764019 -6.9134965 Jun 8 106.82877 3.84047398 -13.8287742 Jul 8 114.25430 6.49518653 9.7457004 Aug 8 121.35502 6.94352630 8.6449763 Sep 8 119.71545 0.58781783 4.2845463 Oct 8 118.99984 -0.37739593 -3.9998372 Nov 8 113.27620 -4.33631398 -7.2762041 Dec 8 108.98889 -4.30002435 -3.9888888 Jan 9 102.03356 -6.26678295 2.9664400 Feb 9 96.09371 -6.02477085 4.9062921 Mar 9 91.38509 -5.05130527 3.6149122 Apr 9 90.97564 -1.61797471 2.0243559 May 9 91.65280 0.08060427 -7.6528014 Jun 9 100.48066 6.55459574 -13.4806616 Jul 9 106.93392 6.47962683 9.0660773 Aug 9 110.03038 3.97723273 9.9696157 Sep 9 112.78961 3.07628689 4.2103925 Oct 9 112.48001 0.57168356 -3.4800056 Nov 9 112.42240 0.10616345 -7.4224024 Dec 9 110.54159 -1.36396794 -3.5415937 Jan 10 106.28643 -3.50317484 2.7135697 Feb 10 103.77928 -2.76659776 5.2207192 Mar 10 104.49329 -0.19434498 3.5067111 Apr 10 105.29555 0.54218348 1.7044524 May 10 108.55240 2.54940115 -9.5524009 Jun 10 115.77387 6.00424742 -12.7738727 Jul 10 121.40351 5.72732529 9.5964940 Aug 10 126.89520 5.55318621 10.1048012 Sep 10 130.48728 4.10392357 4.5127237 Oct 10 128.62115 -0.30815725 -4.6211468 Nov 10 125.35960 -2.49104293 -7.3596000 Dec 10 123.47361 -2.04376606 -2.4736079 Jan 11 118.86661 -3.93856996 2.1333921 Feb 11 113.59924 -4.92038501 4.4007583 Mar 11 109.31149 -4.45319970 3.6885083 Apr 11 105.75469 -3.79121650 1.2453064 May 11 109.41079 1.71079866 -9.4107857 Jun 11 114.38592 4.12275261 -12.3859186 Jul 11 120.33643 5.47286652 9.6635734 Aug 11 125.67516 5.37381338 10.3248413 Sep 11 127.69511 2.89733834 5.3048905 Oct 11 125.05167 -1.19412988 -5.0516727 Nov 11 120.28285 -3.83414636 -8.2828531 Dec 11 111.67563 -7.35979162 -2.6756303 Jan 12 106.75952 -5.55485333 3.2404824 Feb 12 101.34769 -5.44925595 4.6523077 Mar 12 97.77322 -4.06573562 4.2267834 Apr 12 97.77019 -1.06754193 0.2298085 May 12 101.45153 2.43827244 -9.4515296 Jun 12 104.98781 3.24897207 -12.9878103 Jul 12 110.30592 4.77632564 9.6940819 Aug 12 116.02967 5.47546672 10.9703325 Sep 12 118.16239 3.00904687 5.8376056 Oct 12 118.55478 1.07827105 -4.5547758 Nov 12 115.68888 -1.83250198 -7.6888770 Dec 12 109.85710 -4.78440929 -3.8570994 Jan 13 107.05648 -3.32026127 3.9435249 Feb 13 105.27344 -2.18614399 4.7265577 Mar 13 101.01653 -3.71329957 2.9834737 Apr 13 100.23961 -1.54772477 -0.2396143 May 13 104.65541 2.85172005 -8.6554082 Jun 13 111.13057 5.52518899 -13.1305731 Jul 13 113.45435 3.16359399 8.5456497 Aug 13 121.44010 6.71969197 12.5598991 Sep 13 127.04804 5.89995718 5.9519603 > m$resid Jan Feb Mar Apr May 1 0.000000000 -0.922790445 -0.502711826 0.318982775 -0.154135093 2 1.135572991 -0.062065858 -0.447915951 -0.357251321 0.118386644 3 1.069816739 0.339953438 -0.108122120 -0.634848758 0.235286428 4 0.491704165 0.186591208 -0.353801347 -0.131549194 0.523542127 5 0.356474904 0.212189534 -0.488657586 -0.155924402 1.045166108 6 -0.027264008 0.461983751 0.307578025 0.198595227 0.906800820 7 -0.205093755 0.590236504 0.303175828 0.738937494 1.025309578 8 1.052732691 -0.459417967 0.324219456 0.860242572 0.503349949 9 -0.386359118 0.047525076 0.191204353 0.674716203 0.333774791 10 -0.420208701 0.144649711 0.505239990 0.144729501 0.394408789 11 -0.372185068 -0.192814572 0.091765187 0.130073322 1.081087275 12 0.354524287 0.020738244 0.271755040 0.589090578 0.688837545 13 0.287582081 0.222732667 -0.299969914 0.425482301 0.864403038 Jun Jul Aug Sep Oct 1 0.470763373 4.236493365 1.092140210 -3.717559662 -2.720492258 2 1.187622193 3.072204219 0.328180288 -1.598355827 -2.707189740 3 1.762869404 2.979992880 -0.705837413 -1.846655667 -2.079316685 4 2.167378220 1.827812932 -0.655839645 -1.533601972 -1.884404727 5 1.728147553 1.035419190 0.085555604 -1.286855648 -1.873062534 6 1.157816677 0.855637782 -0.045803544 -0.860784795 -1.249615013 7 0.840268949 0.101715177 -0.083712441 -0.574051649 -1.563483852 8 0.603484487 0.521330709 0.088061103 -1.248461482 -0.189605918 9 1.271554666 -0.014722631 -0.491495292 -0.176973913 -0.492016450 10 0.678603185 -0.054384129 -0.034201974 -0.284679628 -0.866750935 11 0.473775312 0.265151991 -0.019454404 -0.486454001 -0.803780657 12 0.159248053 0.299966608 0.137313035 -0.484477299 -0.379310751 13 0.525165572 -0.463816642 0.698424945 -0.161019539 Nov Dec 1 -0.617929817 1.090925520 2 -1.908219516 0.847714524 3 -0.605658831 0.048505663 4 -0.854539698 0.251673028 5 -0.810306497 -0.297583382 6 -1.166565394 -0.662857185 7 -0.890076129 -0.611791488 8 -0.777773364 0.007130318 9 -0.091458157 -0.288841398 10 -0.428864000 0.087874320 11 -0.518677740 -0.692644977 12 -0.571873484 -0.579915362 13 > 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 > postscript(file="/var/wessaorg/rcomp/tmp/17bpc1353674302.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') > grid() > dev.off() null device 1 > mylagmax <- nx/2 > postscript(file="/var/wessaorg/rcomp/tmp/2zgu21353674302.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/3x25z1353674302.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/4jx7s1353674302.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/58ar71353674302.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/6b47e1353674302.tab") > > try(system("convert tmp/17bpc1353674302.ps tmp/17bpc1353674302.png",intern=TRUE)) character(0) > try(system("convert tmp/2zgu21353674302.ps tmp/2zgu21353674302.png",intern=TRUE)) character(0) > try(system("convert tmp/3x25z1353674302.ps tmp/3x25z1353674302.png",intern=TRUE)) character(0) > try(system("convert tmp/4jx7s1353674302.ps tmp/4jx7s1353674302.png",intern=TRUE)) character(0) > try(system("convert tmp/58ar71353674302.ps tmp/58ar71353674302.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 3.649 0.514 4.149