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(79.8,83.4,113.6,112.9,104,109.9,99,106.3,128.9,111.1,102.9,130,87,87.5,117.6,103.4,110.8,112.6,102.5,112.4,135.6,105.1,127.7,137,91,90.5,122.4,123.3,124.3,120,118.1,119,142.7,123.6,129.6,151.6,110.4,99.2,130.5,136.2,129.7,128,121.6,135.8,143.8,147.5,136.2,156.6,123.3,104.5,139.8,136.5,112.1,118.5,94.4,102.3,111.4,99.2,87.8,115.8) > 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 21.3152492 0.1826258 24.0940170 0.0000000 > m$fitted level slope sea Jan 1 79.80000 0.00000000 0.0000000 Feb 1 81.94475 0.07813984 1.4552460 Mar 1 99.63761 1.24462559 13.9623866 Apr 1 108.78271 1.71632885 4.1172925 May 1 108.38671 1.61183251 -4.3867119 Jun 1 108.80131 1.55663900 1.0986919 Jul 1 104.04778 1.25131090 -5.0477812 Aug 1 103.66862 1.16484715 2.6313791 Sep 1 115.29082 1.77028711 13.6091765 Oct 1 115.92369 1.69959716 -4.8236899 Nov 1 109.62656 1.17445025 -6.7265617 Dec 1 117.56464 1.63829357 12.4353602 Jan 2 110.84034 1.21105249 -23.8403397 Feb 2 104.85578 0.74077256 -17.3557751 Mar 2 104.97299 0.69215341 12.6270059 Apr 2 101.49345 0.35131699 1.9065515 May 2 106.16692 0.70202260 4.6330775 Jun 2 108.03282 0.79483181 4.5671784 Jul 2 108.79865 0.79254112 -6.2986504 Aug 2 112.05978 0.98755482 0.3402154 Sep 2 116.97842 1.29976012 18.6215821 Oct 2 113.72908 0.93673983 -8.6290806 Nov 2 124.01513 1.68269923 3.6848686 Dec 2 124.24566 1.56752910 12.7543375 Jan 3 118.59238 0.99728867 -27.5923843 Feb 3 112.40286 0.41927510 -21.9028584 Mar 3 109.36329 0.13507478 13.0367146 Apr 3 115.52713 0.63544563 7.7728699 May 3 118.96139 0.86798476 5.3386099 Jun 3 118.37018 0.74705785 1.6298224 Jul 3 122.17102 0.99936531 -4.0710167 Aug 3 122.46644 0.94134134 -3.4664416 Sep 3 123.42925 0.94310760 19.2707525 Oct 3 130.18026 1.42027008 -6.5802587 Nov 3 128.72034 1.18405550 0.8796624 Dec 3 130.96855 1.27123138 20.6314472 Jan 4 133.14560 1.34551904 -22.7455962 Feb 4 128.22272 0.82951924 -29.0227246 Mar 4 124.29759 0.43673773 6.2024144 Apr 4 126.25496 0.56253888 9.9450441 May 4 125.78590 0.47719727 3.9141026 Jun 4 126.80058 0.52164496 1.1994238 Jul 4 126.62045 0.46363595 -5.0204525 Aug 4 133.20259 0.96909004 2.5974052 Sep 4 132.02307 0.79174762 11.7769294 Oct 4 141.82816 1.53491238 5.6718425 Nov 4 141.37906 1.37147632 -5.1790566 Dec 4 138.70286 1.03814553 17.8971415 Jan 5 139.97316 1.05727488 -16.6731569 Feb 5 136.48626 0.68242495 -31.9862582 Mar 5 134.95968 0.50007001 4.8403202 Apr 5 130.67026 0.10469302 5.8297404 May 5 119.61393 -0.81637377 -7.5139290 Jun 5 116.69580 -0.98979948 1.8041992 Jul 5 109.50174 -1.50178539 -15.1017403 Aug 5 103.50392 -1.87280877 -1.2039238 Sep 5 102.06309 -1.83717128 9.3369139 Oct 5 96.38451 -2.15392586 2.8154935 Nov 5 91.96862 -2.34035476 -4.1686151 Dec 5 92.74415 -2.08358249 23.0558504 > m$resid Jan Feb Mar Apr May Jun 1 0.000000000 0.364584540 2.778881761 1.541026753 -0.439961061 -0.250021200 2 -1.760961880 -1.475893142 -0.117046266 -0.785454475 0.840338002 0.229176248 3 -1.424172118 -1.413494547 -0.665510284 1.150996429 0.539379029 -0.283726688 4 0.176453762 -1.219538149 -0.917659913 0.292105753 -0.198744328 0.104084290 5 0.045010802 -0.880387043 -0.426348128 -0.921941934 -2.150643153 -0.406165396 Jul Aug Sep Oct Nov Dec 1 -1.307667570 -0.335273848 2.135791851 -0.230976431 -1.616062012 1.361193119 2 -0.005706375 0.484216155 0.769396468 -0.888623308 1.822149600 -0.283278216 3 0.594783317 -0.136788746 0.004158775 1.122795473 -0.556609253 0.206305259 4 -0.136175318 1.185475298 -0.415106357 1.737840060 -0.382647759 -0.782735681 5 -1.200892640 -0.869533787 0.083367903 -0.740254823 -0.436019764 0.601698345 > 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/1gnzq1259940481.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/2rkcr1259940481.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/3jq651259940481.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/4xsfa1259940481.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/50qj81259940481.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/6sp8o1259940481.tab") > system("convert tmp/1gnzq1259940481.ps tmp/1gnzq1259940481.png") > system("convert tmp/2rkcr1259940481.ps tmp/2rkcr1259940481.png") > system("convert tmp/3jq651259940481.ps tmp/3jq651259940481.png") > system("convert tmp/4xsfa1259940481.ps tmp/4xsfa1259940481.png") > system("convert tmp/50qj81259940481.ps tmp/50qj81259940481.png") > > > proc.time() user system elapsed 1.388 0.808 1.579