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(100,96.21064363,96.31280765,107.1793443,114.9066592,92.56060184,114.9995356,107.1236185,117.7765394,107.3650971,106.2970187,114.5072908,98.0031578,103.0649206,100.2879168,104.6066685,111.1544534,104.9874617,109.9284852,111.5352466,132.4974459,100.3436426,123.0983561,114.2379493,104.569518,109.0833101,106.9843039,133.6769759,124.8537197,122.5132349,116.8013374,116.0118882,129.7575926,125.1973623,143.7912139,127.9465032,130.2962757,108.4424631,129.3675118,143.6797622,131.8844618,117.6186496,118.9560695,104.8202842,134.624315,140.401226,143.8005015,153.4317823,153.2924677,127.3149438,153.5525216,136.9276493,131.7730101,144.3391845,107.4208229,113.6249652,124.2221603,102.0618557,96.36853348,111.6838488) > 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 32.72600 0.00000 13.61853 25.40493 > m$fitted level slope sea Jan 1 100.00000 0.000000000 0.00000000 Feb 1 97.53737 -0.143860738 -0.68171387 Mar 1 96.67605 -0.182276625 -0.14489012 Apr 1 102.12500 0.004317751 2.95980938 May 1 109.32658 0.156170055 2.75277520 Jun 1 101.43057 0.025257720 -5.66563470 Jul 1 107.56235 0.114662007 4.99959671 Aug 1 108.00754 0.119328934 -1.01595020 Sep 1 113.09646 0.188359029 2.69464568 Oct 1 110.67395 0.152550877 -2.26567071 Nov 1 108.11325 0.115808854 -0.73203695 Dec 1 111.11847 0.154429597 2.23410032 Jan 2 106.47289 0.414337131 -6.45891247 Feb 2 104.79329 0.384711891 -1.00189560 Mar 2 103.26542 0.332871415 -2.35600892 Apr 2 103.15876 0.322987967 1.60098536 May 2 104.52022 0.339372331 6.25675764 Jun 2 108.51833 0.381587090 -4.88095255 Jul 2 107.48391 0.367956484 2.96942263 Aug 2 110.59210 0.392634154 -0.07361318 Sep 2 119.83723 0.471586343 9.37517164 Oct 2 112.47096 0.402403526 -9.21849161 Nov 2 117.50901 0.438267831 3.86752661 Dec 2 114.72706 0.429886069 0.71136017 Jan 3 112.55167 0.470105711 -6.98151024 Feb 3 110.91517 0.457697890 -1.07565710 Mar 3 109.98319 0.434290852 -2.52689960 Apr 3 120.86125 0.614237465 9.19863485 May 3 121.64561 0.616593834 3.14749310 Jun 3 124.56934 0.641096412 -2.89062720 Jul 3 120.87324 0.603558852 -2.49373762 Aug 3 119.74201 0.590103457 -3.09761355 Sep 3 118.83764 0.579110777 11.46521076 Oct 3 126.91819 0.629649393 -4.45893863 Nov 3 133.44781 0.658628732 8.18777017 Dec 3 131.08258 0.656420040 -2.02734275 Jan 4 133.95655 0.645690061 -4.48348274 Feb 4 123.36096 0.598266385 -10.87297711 Mar 4 128.61796 0.652702168 -0.87033658 Apr 4 132.21424 0.692565429 10.43936042 May 4 131.33425 0.673483321 1.10783649 Jun 4 125.46760 0.608359722 -5.49990446 Jul 4 122.83429 0.581681911 -2.70625238 Aug 4 115.88128 0.527876482 -8.32920641 Sep 4 120.40572 0.553302666 12.76785286 Oct 4 133.50339 0.620345664 2.34143232 Nov 4 135.63130 0.625666354 7.62116718 Dec 4 145.50949 0.633616176 4.55045585 Jan 5 150.93520 0.629860458 0.60482782 Feb 5 147.01653 0.612005981 -18.06649347 Mar 5 150.95790 0.641620206 1.42289157 Apr 5 139.30023 0.506927433 1.94534277 May 5 132.74192 0.432088942 1.53592590 Jun 5 139.16948 0.487359612 3.02255580 Jul 5 125.18760 0.374528298 -12.55181920 Aug 5 123.79371 0.362771627 -9.52929085 Sep 5 120.50031 0.342135327 5.04593023 Oct 5 111.40868 0.299890063 -5.92800102 Nov 5 101.87239 0.271006811 -1.93685180 Dec 5 103.91625 0.273360174 7.12383290 > m$resid Jan Feb Mar Apr May Jun 1 0.00000000 -0.31017488 -0.10499447 0.94146683 1.24356834 -1.40104564 2 -0.98299615 -0.34189420 -0.30305673 -0.07360468 0.17900730 0.63638930 3 -0.47719005 -0.35989262 -0.22917580 1.75790143 0.02925581 0.40073672 4 0.39354021 -1.93772917 0.78432292 0.49877130 -0.27029152 -1.13463225 5 0.84066691 -0.78591930 0.56675337 -2.09631850 -1.21531375 1.03945457 Jul Aug Sep Oct Nov Dec 1 1.06408405 0.05761078 0.86622375 -0.45509077 -0.47294636 0.50366872 2 -0.24687062 0.47796701 1.54419185 -1.36727297 0.80853709 -0.56228221 3 -0.75600144 -0.30271167 -0.26085251 1.30918638 1.02927782 -0.52878895 4 -0.56478958 -1.31489170 0.69774665 2.18952253 0.26304938 1.61736831 5 -2.51970021 -0.30855298 -0.63826519 -1.64646165 -1.71630926 0.30960531 > 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/www/html/rcomp/tmp/1jua81259944577.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(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/www/html/rcomp/tmp/20ycl1259944577.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(mylevel,main='Level') > spectrum(myseas,main='Seasonal') > spectrum(myresid,main='Standardized Residals') > par(op) > dev.off() null device 1 > postscript(file="/var/www/html/rcomp/tmp/3k5qr1259944577.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(mylevel,main='Level') > cpgram(myseas,main='Seasonal') > cpgram(myresid,main='Standardized Residals') > par(op) > dev.off() null device 1 > postscript(file="/var/www/html/rcomp/tmp/4wvqf1259944577.ps",horizontal=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/www/html/rcomp/tmp/5k84k1259944577.ps",horizontal=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/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,'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/www/html/rcomp/tmp/66o831259944577.tab") > > system("convert tmp/1jua81259944577.ps tmp/1jua81259944577.png") > system("convert tmp/20ycl1259944577.ps tmp/20ycl1259944577.png") > system("convert tmp/3k5qr1259944577.ps tmp/3k5qr1259944577.png") > system("convert tmp/4wvqf1259944577.ps tmp/4wvqf1259944577.png") > system("convert tmp/5k84k1259944577.ps tmp/5k84k1259944577.png") > > > proc.time() user system elapsed 1.539 0.850 1.782