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(104.31,103.88,103.88,103.86,103.89,103.98,103.98,104.29,104.29,104.24,103.98,103.54,103.44,103.32,103.3,103.26,103.14,103.11,102.91,103.23,103.23,103.14,102.91,102.42,102.1,102.07,102.06,101.98,101.83,101.75,101.56,101.66,101.65,101.61,101.52,101.31,101.19,101.11,101.1,101.07,100.98,100.93,100.92,101.02,101.01,100.97,100.89,100.62,100.53,100.48,100.48,100.47,100.52,100.49,100.47,100.44)
> 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.0210641663 0.0001903905 0.0000000000 0.0000000000 
> m$fitted
         level        slope          sea
Jan 1 104.3100  0.000000000  0.000000000
Feb 1 103.9023 -0.023550291 -0.022294706
Mar 1 103.9022 -0.023062792 -0.022212705
Apr 1 103.8822 -0.022974338 -0.022202348
May 1 103.9120 -0.021039007 -0.022029756
Jun 1 104.0017 -0.016199296 -0.021683756
Jul 1 104.0016 -0.015389438 -0.021635803
Aug 1 104.3107  0.002748251 -0.020726263
Sep 1 104.3107  0.002581050 -0.020733478
Oct 1 104.2609 -0.000853206 -0.020862497
Nov 1 104.0015 -0.018787852 -0.021454363
Dec 1 103.5623 -0.049353992 -0.022346560
Jan 2 103.2932 -0.060565364  0.146757989
Feb 2 103.3310 -0.051980109 -0.011042252
Mar 2 103.3110 -0.049438050 -0.011022703
Apr 2 103.2710 -0.048670476 -0.011017401
May 2 103.1511 -0.054582048 -0.011054149
Jun 2 103.1210 -0.052512774 -0.011042550
Jul 2 102.9211 -0.065088771 -0.011106204
Aug 2 103.2410 -0.031901961 -0.010954328
Sep 2 103.2409 -0.029128406 -0.010942840
Oct 2 103.1510 -0.034459138 -0.010962841
Nov 2 102.9210 -0.051686372 -0.011021429
Dec 2 102.4311 -0.090493997 -0.011141129
Jan 3 102.0419 -0.115736025  0.058129519
Feb 3 102.0731 -0.102383470 -0.003119677
Mar 3 102.0631 -0.094121432 -0.003117251
Apr 3 101.9831 -0.092855540 -0.003116912
May 3 101.8331 -0.097988169 -0.003118160
Jun 3 101.7531 -0.096369889 -0.003117803
Jul 3 101.5631 -0.104804433 -0.003119496
Aug 3 101.6631 -0.086334532 -0.003116126
Sep 3 101.6531 -0.079444164 -0.003114984
Oct 3 101.6131 -0.075881033 -0.003114447
Nov 3 101.5231 -0.077157245 -0.003114622
Dec 3 101.3131 -0.089171082 -0.003116118
Jan 4 101.1643 -0.094494764  0.025744861
Feb 4 101.1119 -0.090682321 -0.001885695
Mar 4 101.1019 -0.083377385 -0.001891691
Apr 4 101.0719 -0.078543515 -0.001895291
May 4 100.9819 -0.079581214 -0.001894588
Jun 4 100.9319 -0.076901403 -0.001896238
Jul 4 100.9219 -0.070839908 -0.001899632
Aug 4 101.0219 -0.055359568 -0.001907512
Sep 4 101.0119 -0.051249031 -0.001909415
Oct 4 100.9719 -0.050229555 -0.001909844
Nov 4 100.8919 -0.052927750 -0.001908811
Dec 4 100.6219 -0.072602716 -0.001901964
Jan 5 100.5130 -0.075875666  0.017038242
Feb 5 100.4812 -0.071893471 -0.001172370
Mar 5 100.4812 -0.065375873 -0.001178342
Apr 5 100.4712 -0.060355670 -0.001182514
May 5 100.5212 -0.050351118 -0.001190074
Jun 5 100.4912 -0.048506132 -0.001191342
Jul 5 100.4712 -0.045921822 -0.001192957
Aug 5 100.4412 -0.044478376 -0.001193777
> m$resid
          Jan         Feb         Mar         Apr         May         Jun
1  0.00000000 -1.68396802  0.16057652  0.02079609  0.35827107  0.74821509
2 -1.73779671  0.55514160  0.21140850  0.06232905 -0.47066402  0.16208820
3 -2.14835894  0.88181935  0.60740502  0.09283514 -0.37563496  0.11823450
4 -0.41661527  0.25949167  0.53014937  0.35072936 -0.07527696  0.19436736
5 -0.24992154  0.27551719  0.47236907  0.36384159  0.72508062  0.13371491
          Jul         Aug         Sep         Oct         Nov         Dec
1  0.10878958  2.17859634 -0.01835077 -0.35025963 -1.72266318 -2.79492399
2 -0.97191686  2.53640812  0.21003657 -0.40062820 -1.28658622 -2.88324803
3 -0.61538107  1.34599737  0.50165629  0.25920998 -0.09278086 -0.87293568
4  0.43958185  1.12251339  0.29803632  0.07391180 -0.19560624 -1.42626624
5  0.18729655  0.10461264                                                
> 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/1d7rg1292675817.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/www/html/rcomp/tmp/25h8j1292675817.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/www/html/rcomp/tmp/3gqp41292675817.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/www/html/rcomp/tmp/4gqp41292675817.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/www/html/rcomp/tmp/5rz7p1292675817.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/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/6uinc1292675817.tab") 
> 
> try(system("convert tmp/1d7rg1292675817.ps tmp/1d7rg1292675817.png",intern=TRUE))
character(0)
> try(system("convert tmp/25h8j1292675817.ps tmp/25h8j1292675817.png",intern=TRUE))
character(0)
> try(system("convert tmp/3gqp41292675817.ps tmp/3gqp41292675817.png",intern=TRUE))
character(0)
> try(system("convert tmp/4gqp41292675817.ps tmp/4gqp41292675817.png",intern=TRUE))
character(0)
> try(system("convert tmp/5rz7p1292675817.ps tmp/5rz7p1292675817.png",intern=TRUE))
character(0)
> 
> 
> proc.time()
   user  system elapsed 
  1.347   0.843   8.327