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(15836.8,17570.4,18252.1,16196.7,16643,17729,16446.1,15993.8,16373.5,17842.2,22321.5,22786.7,18274.1,22392.9,23899.3,21343.5,22952.3,21374.4,21164.1,20906.5,17877.4,20664.3,22160,19813.6,17735.4,19640.2,20844.4,19823.1,18594.6,21350.6,18574.1,18924.2,17343.4,19961.2,19932.1,19464.6,16165.4,17574.9,19795.4,19439.5,17170,21072.4,17751.8,17515.5,18040.3,19090.1,17746.5,19202.1,15141.6,16258.1,18586.5,17209.4,17838.7,19123.5,16583.6,15991.2,16704.4,17420.4,17872,17823.2) > 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 537160.7035 975.8086 161850.9246 425256.0460 > m$fitted level slope sea Jan 1 15836.80 0.000000 0.00000 Feb 1 17006.42 71.528251 254.21075 Mar 1 17782.38 107.589920 231.49985 Apr 1 16953.49 77.115768 -379.66699 May 1 16698.49 69.163753 83.64256 Jun 1 17223.34 79.071102 313.45383 Jul 1 16831.88 68.715461 -187.39396 Aug 1 16324.08 55.451177 -87.67801 Sep 1 16289.86 53.285742 121.28248 Oct 1 17120.92 72.942649 395.49076 Nov 1 20032.39 147.711073 1102.71244 Dec 1 21799.78 192.004812 311.36707 Jan 2 21130.94 219.577047 -2486.57854 Feb 2 21812.18 235.448452 418.71205 Mar 2 22790.21 266.767287 855.45435 Apr 2 22334.11 239.788299 -724.70389 May 2 22655.56 242.468708 265.74817 Jun 2 21824.73 209.288741 -40.06862 Jul 2 21414.00 190.259146 -12.81826 Aug 2 21121.48 175.198196 -30.59490 Sep 2 19468.26 116.834418 -893.59510 Oct 2 20085.18 133.137682 388.67469 Nov 2 20861.06 154.167157 1054.44550 Dec 2 20169.63 128.590813 -34.62011 Jan 3 20267.21 128.035614 -2519.78783 Feb 3 19879.62 109.927828 -52.11983 Mar 3 19836.12 103.615804 1061.73534 Apr 3 20243.15 115.893053 -528.72884 May 3 19184.35 70.740896 -159.74370 Jun 3 20131.07 103.233012 896.65240 Jul 3 19345.49 70.693441 -443.45848 Aug 3 18821.01 48.874678 322.69554 Sep 3 18594.45 38.710470 -1149.54658 Oct 3 19169.89 58.588429 593.65481 Nov 3 18960.74 48.784687 1069.83783 Dec 3 19213.77 55.933949 175.75195 Jan 4 18890.65 43.675302 -2584.25393 Feb 4 18213.86 16.505646 -377.25616 Mar 4 18438.52 25.000395 1283.34882 Apr 4 19073.77 50.076032 148.24663 May 4 18380.35 20.099941 -941.68529 Jun 4 19060.63 46.223961 1771.80605 Jul 4 18585.75 25.816792 -644.33928 Aug 4 17822.17 -4.994246 -19.49579 Sep 4 18524.91 22.614899 -741.94207 Oct 4 18536.50 22.185933 557.60826 Nov 4 17622.10 -13.897860 464.21120 Dec 4 18154.37 6.746303 849.50895 Jan 5 17873.80 -3.972932 -2627.16786 Feb 5 17283.33 -26.962625 -812.93019 Mar 5 17311.60 -24.718104 1255.24115 Apr 5 17067.72 -33.696943 219.82170 May 5 17985.13 5.024083 -488.44179 Jun 5 17626.26 -9.657137 1628.72523 Jul 5 17291.53 -22.684783 -590.35154 Aug 5 16724.63 -44.413468 -536.65123 Sep 5 16984.39 -32.298780 -389.89400 Oct 5 16785.75 -38.901106 694.68807 Nov 5 17145.65 -23.167642 582.56373 Dec 5 17004.02 -27.806550 861.94025 > m$resid Jan Feb Mar Apr May Jun 1 0.00000000 1.12373987 0.82710658 -1.23079517 -0.44660989 0.61469939 2 -1.35672553 0.56517049 0.90844002 -0.93411351 0.10762370 -1.41948514 3 -0.04296323 -0.66136795 -0.19211490 0.38923055 -1.52970107 1.14498227 4 -0.50471518 -0.93219610 0.26422630 0.78232345 -0.96315078 0.85822592 5 -0.37687259 -0.75978081 0.07065043 -0.28133887 1.22956501 -0.47191178 Jul Aug Sep Oct Nov Dec 1 -0.63420820 -0.77583546 -0.12045673 1.04310613 3.80073472 2.16536427 2 -0.81968960 -0.63757692 -2.41203728 0.65906199 0.84597489 -1.11140040 3 -1.16177959 -0.77752586 -0.35960418 0.70023753 -0.34880053 0.26629471 4 -0.67753094 -1.02584771 0.91918723 -0.01430626 -1.21396646 0.70940583 5 -0.42165851 -0.70546695 0.39397392 -0.21522878 0.51572264 -0.15357425 > 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/1ihl81259952362.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/2gc041259952362.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/3qziq1259952362.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/4e5sp1259952362.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/58yum1259952362.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/618ue1259952362.tab") > system("convert tmp/1ihl81259952362.ps tmp/1ihl81259952362.png") > system("convert tmp/2gc041259952362.ps tmp/2gc041259952362.png") > system("convert tmp/3qziq1259952362.ps tmp/3qziq1259952362.png") > system("convert tmp/4e5sp1259952362.ps tmp/4e5sp1259952362.png") > system("convert tmp/58yum1259952362.ps tmp/58yum1259952362.png") > > > proc.time() user system elapsed 1.637 0.807 2.064