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(0.0589999409995141 + ,-1.41417450501306 + ,-18.3846475755861 + ,4.24266776527226 + ,-9.89940750475467 + ,-1.41416601970523 + ,-4.94968147844429 + ,8.48527118362683 + ,5.65684722004823 + ,2.12135156036077 + ,-12.7278137906086 + ,-1.41417238368611 + ,-5.65678640867546 + ,5.71546709301464 + ,11.430917448217 + ,4.0824915417843 + ,-5.71542749248443 + ,-4.89893277693375 + ,5.30722585896545 + ,-8.98140023237252 + ,-0.816480018948223 + ,-3.67419274460486 + ,2.44950292718234 + ,-7.34841447478369 + ,-12.2473676641251 + ,13.5676795115537 + ,-7.50550178552387 + ,-10.1035605510469 + ,4.61879420118808 + ,-12.9903024874634 + ,3.7527758629083 + ,-14.1450008210891 + ,-2.30938107887528 + ,-11.2582591729112 + ,-4.33008822361948 + ,1.29908651985677e-05 + ,16.454417696876 + ,7.82621147962468 + ,4.02491438422167 + ,-8.72060773063855 + ,-0.894411455568734 + ,-0.223594253653449 + ,-6.93176357788012 + ,-2.01244265593132 + ,-7.1553761680782 + ,0.223617958497356 + ,8.27342503551901 + ,15.2051986751483 + ,2.90688653528183 + ,-12.7801213366766 + ,-3.10373653589491 + ,-7.1203466912372 + ,-3.46888342810694 + ,-8.3983624583384 + ,2.55603701106735 + ,0.182580385499283 + ,10.5892565952768 + ,15.7013138196074 + ,6.75522353446613 + ,5.11206123957597 + ,8.76352852127385) > par1 = '12' > 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.000000000 0.002505017 0.000000000 61.193315047 > m$fitted level slope sea Jan 1 0.05899994 0.0000000000 0.0000000000 Feb 1 -0.12329239 -0.0607976827 -0.0607473259 Mar 1 -3.15995244 -0.8057354670 -0.8040358849 Apr 1 -2.25220638 -0.4622483748 -0.4619912279 May 1 -4.08133813 -0.6908963469 -0.6890517361 Jun 1 -3.93347767 -0.5704147145 -0.5698551715 Jul 1 -4.47841370 -0.5672044511 -0.5666950976 Aug 1 -2.22363437 -0.2501131098 -0.2566248062 Sep 1 -0.85242321 -0.0854535407 -0.0969880540 Oct 1 -0.35159149 -0.0310300680 -0.0447982076 Nov 1 -2.57324059 -0.2187025190 -0.2223236826 Dec 1 -2.51889230 -0.1969286683 -0.2020712083 Jan 2 -2.95596695 -0.2198945380 2.3211856036 Feb 2 -1.52010069 -0.1035597558 -0.1062336433 Mar 2 0.61875353 0.0459337951 0.0251075384 Apr 2 1.21003238 0.0806494991 0.0547095542 May 2 0.21351783 0.0147810103 0.0004934103 Jun 2 -0.51798541 -0.0293036033 -0.0343294705 Jul 2 0.27410940 0.0178185001 0.0011722013 Aug 2 -0.95302404 -0.0518866322 -0.0485924235 Sep 2 -0.97417583 -0.0501993691 -0.0474586206 Oct 2 -1.35149062 -0.0678812073 -0.0585642321 Nov 2 -0.93925407 -0.0422297963 -0.0436139098 Dec 2 -1.73554663 -0.0821867185 -0.0650636738 Jan 3 -2.78282910 -0.1486931684 0.6706615618 Feb 3 -0.93231905 -0.0458488178 -0.0101837785 Mar 3 -1.74513946 -0.0852622647 -0.0278552886 Apr 3 -2.77959072 -0.1341021589 -0.0475660289 May 3 -2.05743285 -0.0899040580 -0.0316561758 Jun 3 -3.34915916 -0.1522271387 -0.0514718278 Jul 3 -2.70012304 -0.1104382630 -0.0398602124 Aug 3 -4.03563125 -0.1747423218 -0.0552902796 Sep 3 -4.00008778 -0.1636285010 -0.0530187077 Oct 3 -4.91524639 -0.2036274290 -0.0598694111 Nov 3 -5.02879339 -0.1987994911 -0.0591904350 Dec 3 -4.66868763 -0.1686412275 -0.0557985547 Jan 4 -2.76825836 -0.0343022212 1.1617219368 Feb 4 -1.66317078 0.0267588098 -0.0983591329 Mar 4 -1.02749674 0.0596248730 -0.0952552230 Apr 4 -1.77519675 0.0157560754 -0.0984243884 May 4 -1.65803641 0.0213005646 -0.0981355555 Jun 4 -1.47780266 0.0300398559 -0.0978392393 Jul 4 -2.01414142 -0.0012621004 -0.0983894738 Aug 4 -2.00474104 -0.0006701709 -0.0983876965 Sep 4 -2.53734571 -0.0303161204 -0.0980932799 Oct 4 -2.26302996 -0.0132813258 -0.0984566763 Nov 4 -1.15199718 0.0497674503 -0.1004307409 Dec 4 0.63303749 0.1473052570 -0.1043291823 Jan 5 0.88582797 0.1537072970 1.1803065190 Feb 5 -0.40591451 0.0728595204 -0.1072798991 Mar 5 -0.61413521 0.0570986123 -0.1071440004 Apr 5 -1.23927483 0.0187626175 -0.1065306272 May 5 -1.44716188 0.0060043632 -0.1062441286 Jun 5 -2.16716794 -0.0349220400 -0.1050958742 Jul 5 -1.68598531 -0.0057950121 -0.1060531413 Aug 5 -1.48131717 0.0060932831 -0.1064923278 Sep 5 -0.18005960 0.0793007021 -0.1094460336 Oct 5 1.59457439 0.1751717233 -0.1135812675 Nov 5 2.31366702 0.2059396551 -0.1149767069 Dec 5 2.80871457 0.2222962586 -0.1157464117 Jan 6 3.51611392 0.2501822229 1.4169090591 > m$resid Jan Feb Mar Apr May Jun 1 0.00000000 -0.17208834 -2.06755436 1.01221982 -0.75066371 0.45200344 2 -0.65098223 1.04502332 1.52314349 0.39524986 -0.82714632 -0.60351227 3 -1.35865982 1.98160392 -0.78135640 -0.99016992 0.91154241 -1.30163441 4 2.47042814 1.29686483 0.69606948 -0.92563699 0.11648143 0.18273261 5 0.11597281 -1.65795195 -0.32201708 -0.78053809 -0.25892865 -0.82816285 6 0.52833545 Jul Aug Sep Oct Nov Dec 1 0.01390948 1.58940722 0.95155697 0.36036030 -1.41296477 0.18485189 2 0.69587456 -1.09973699 0.02818927 -0.31030091 0.46935326 -0.75714845 3 0.88042347 -1.36224177 0.23610374 -0.85029540 0.10252229 0.63888426 4 -0.65129335 0.01225550 -0.61084034 0.34937263 1.28747233 1.98374149 5 0.58783932 0.23936489 1.47095842 1.92290807 0.61618381 0.32715459 6 > 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/fisher/rcomp/tmp/1rowz1353967296.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/fisher/rcomp/tmp/2wjlm1353967296.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/fisher/rcomp/tmp/392ke1353967296.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/fisher/rcomp/tmp/4ga941353967296.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',type='b') > grid() > dev.off() null device 1 > postscript(file="/var/fisher/rcomp/tmp/5frfy1353967296.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/fisher/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/fisher/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/fisher/rcomp/tmp/6kdzw1353967296.tab") > > try(system("convert tmp/1rowz1353967296.ps tmp/1rowz1353967296.png",intern=TRUE)) character(0) > try(system("convert tmp/2wjlm1353967296.ps tmp/2wjlm1353967296.png",intern=TRUE)) character(0) > try(system("convert tmp/392ke1353967296.ps tmp/392ke1353967296.png",intern=TRUE)) character(0) > try(system("convert tmp/4ga941353967296.ps tmp/4ga941353967296.png",intern=TRUE)) character(0) > try(system("convert tmp/5frfy1353967296.ps tmp/5frfy1353967296.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 3.194 0.685 3.871