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(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) > 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 0.000000 36.208110 4.650892 0.000000 > m$fitted level slope sea Jan 1 149.00000 0.00000000 0.00000000 Feb 1 140.64268 -8.49307815 -1.64268058 Mar 1 133.50007 -7.48479596 1.49992567 Apr 1 129.28328 -4.92518488 0.71672376 May 1 126.65138 -3.16631002 0.34861744 Jun 1 122.31185 -4.06655159 -0.31184994 Jul 1 117.09187 -4.95319532 -0.09186731 Aug 1 111.89855 -5.13780075 0.10145373 Sep 1 111.82496 -1.24511431 1.17504117 Oct 1 142.44672 23.24742381 6.55328075 Nov 1 160.10422 18.95106308 -3.10421870 Dec 1 160.51967 4.70448912 -3.51966687 Jan 2 149.42313 -7.41655389 -2.42313093 Feb 2 138.65136 -9.99854288 -1.65136447 Mar 2 130.53499 -8.61043185 1.46500647 Apr 2 124.48676 -6.70924648 0.51324031 May 2 121.87569 -3.66906383 1.12430542 Jun 2 117.15475 -4.44511553 -0.15474567 Jul 2 112.24192 -4.79063041 1.75807819 Aug 2 109.91471 -2.96936430 1.08529493 Sep 2 116.97624 4.44713747 -4.97624230 Oct 2 135.03977 14.51411010 8.96023493 Nov 2 150.23712 15.01915649 -0.23711739 Dec 2 152.11774 5.31475787 -3.11773550 Jan 3 139.24826 -8.11523536 -5.24825696 Feb 3 125.53757 -12.25286189 -2.53757411 Mar 3 114.62277 -11.27519939 1.37722509 Apr 3 114.93419 -2.80887045 2.06580996 May 3 110.12284 -4.27666918 0.87716482 Jun 3 104.31250 -5.39670563 0.68750260 Jul 3 99.42151 -5.02752635 2.57848610 Aug 3 95.72145 -4.05752352 -0.72145396 Sep 3 100.41436 2.33841009 -7.41436219 Oct 3 114.35518 10.81852104 9.64481842 Nov 3 127.55530 12.55850168 2.44469876 Dec 3 125.73458 2.06014740 -1.73457923 Jan 4 119.77631 -3.79779900 -4.77631147 Feb 4 110.05878 -8.12081464 -4.05878135 Mar 4 105.90411 -5.24063090 -0.90410869 Apr 4 102.15770 -4.15587499 2.84229540 May 4 99.07897 -3.37194083 1.92102594 Jun 4 93.93705 -4.65785530 1.06294777 Jul 4 89.62051 -4.41014267 3.37949359 Aug 4 86.58700 -3.41035570 -2.58700368 Sep 4 95.63734 5.64313017 -8.63734446 Oct 4 106.97588 9.78108235 9.02412267 Nov 4 114.69223 8.28158484 5.30776797 Dec 4 117.67598 4.43523333 -0.67597633 Jan 5 113.87490 -1.54813087 -4.87490408 Feb 5 110.47291 -2.89384824 -5.47291470 Mar 5 108.30025 -2.37204397 -1.30024548 Apr 5 106.30520 -2.09933142 2.69479731 May 5 105.74182 -0.98629899 3.25818327 Jun 5 105.87081 -0.17868212 2.12919346 Jul 5 103.82919 -1.52632853 3.17080531 Aug 5 104.54645 0.09714549 -5.54645453 Sep 5 112.13104 5.51737273 -9.13104231 Oct 5 121.36637 8.20885579 9.63362989 Nov 5 130.38094 8.79195570 6.61905689 Dec 5 134.49814 5.40894924 0.50186147 > m$resid Jan Feb Mar Apr May Jun 1 0.00000000 -1.53765266 0.17964659 0.42925615 0.29107935 -0.14978954 2 -2.01924576 -0.43543369 0.23049192 0.32002040 0.50384762 -0.12890454 3 -2.23786367 -0.68856124 0.16226182 1.41518045 -0.24393609 -0.18589315 4 -0.97525779 -0.71829141 0.47830145 0.18079808 0.13037576 -0.21347629 5 -0.99539323 -0.22354861 0.08668039 0.04539754 0.18511668 0.13412391 Jul Aug Sep Oct Nov Dec 1 -0.14735085 -0.03067898 0.64691914 4.07033839 -0.71399948 -2.36759559 2 -0.05744911 0.30266238 1.23251262 1.67305691 0.08393235 -1.61305751 3 0.06137145 0.16121987 1.06288313 1.40936625 0.28915845 -1.74557285 4 0.04116389 0.16617477 1.50454521 0.68768885 -0.24920418 -0.63958571 5 -0.22389993 0.26981961 0.90077728 0.44729139 0.09691246 -0.56250657 > 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/12n6r1259672224.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/2sj961259672224.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/3euda1259672224.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/4vr171259672224.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/5aqay1259672224.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/6yq391259672224.tab") > system("convert tmp/12n6r1259672224.ps tmp/12n6r1259672224.png") > system("convert tmp/2sj961259672224.ps tmp/2sj961259672224.png") > system("convert tmp/3euda1259672224.ps tmp/3euda1259672224.png") > system("convert tmp/4vr171259672224.ps tmp/4vr171259672224.png") > system("convert tmp/5aqay1259672224.ps tmp/5aqay1259672224.png") > > > proc.time() user system elapsed 1.259 0.810 1.520