Thanks! Following Prof. Ripley's reply I inserted the following lines into StructTS (makeBSM)
if (nf < 3) { ind <- 2L } else { ind <- 3:nf T[cbind(ind + 1L, ind)] <- 1 } This worked and gave me the same variance estimates as the dlm package. I also see that the StructTS function was updated and can handle 2 seasons now. Thanks for the help! Birgit I have replreplied that StructI have replaced t to the help from Giovanni Petris and Brian Ripley From: Giovanni Petris <gpet...@uark.edu> To: Birgit Erni <birgit.e...@uct.ac.za> CC: <r-help@r-project.org> Date: 30 November 2010 05:40 PM Subject: Re: [R] StructTS with 2 seasons Hi Birgit, You are right, StructTS does not seem to be able to handle a time series with frequency 2 (and after a quick glance at the code I suspect a frequency 3 may not work either). There are contributed packages that you can use to fit (and smooth, filter, forecast) structural time series and more general state space models. Among them I would look at dlm, KFAS, dse, sspir. HTH, Giovanni (author of 'dlm') On Tue, 2010-11-30 at 10:40 +0200, Birgit Erni wrote: > Dear All, > > I am trying to fit a structural time series model using the StructTS > function (package stats) with only 2 seasons (summer and winter). More > than 2 seasons work fine but with 2 seasons I get this error: > > > fit <- StructTS(y.ts, type="BSM") > Error in T[cbind(ind + 1L, ind)] <- 1 : subscript out of bounds > > I have looked at Prof. Ripley's 2002 RNews article but cannot see why > theoretically it should not be possible to have only 2 seasons. I have > also looked at the StructTS code, and the problem is with the transition > matrix in makeBSM, but I am not sure whether I should overwrite this or > whether my thinking is wrong. > > I attach code that simulates data for the basic structural model (with > time-varying trend and season), the top line s determines the number of > seasons in a year (3 and 4 work fine, 2 gives the above error), and at > the end tries to fit the model. I am using R version 2.11.1, Windows XP. > > > # ----- Basic structural model ----- > > s <- 2 # frequency or period > > mu0 <- 6 > beta0 <- 1 > it0 <- 0 > it.states <- sample(1:6,s) > > mu <- numeric() # mean > beta <- numeric() # slope > it <- numeric() # seasonal state > > mu[1] <- mu0 + beta0 + rnorm(1,0,2) > beta[1] <- beta0 + rnorm(1,0,0.5) > it[1] <- - sum(it.states[2:s]) + rnorm(1,0,1) > > it.states <- c(it[1],it.states[2:s]) > > for (i in 2:30) > { mu[i] <- mu[i-1] + beta[i-1] + rnorm(1,0,2) > beta[i] <- beta[i-1] + rnorm(1,0,2) > it[i] <- - sum(it.states[1:(s-1)]) + rnorm(1,0,1) > it.states <- c(it[i],it.states[1:(s-1)]) > } > > y <- numeric() > > for (i in 1:30) > { y[i] <- mu[i] + it[i] + rnorm(1,0,3) } > > par(mfrow=c(2,2)) > plot(y,type="b") > plot(mu) > plot(beta) > plot(it) > > y.ts <- ts(y,start=c(1977,1),frequency=s) > y.ts > > fit <- StructTS(y.ts, type="BSM") > fit > > # ------------------------------------------------------ > > Thank you. > > Birgit Erni > > Department of Statistical Sciences > University of Cape Town, South Africa > > > > > > ### > UNIVERSITY OF CAPE TOWN > > This e-mail is subject to the UCT ICT policies and e-mai...{{dropped:13}} > > ______________________________________________ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Giovanni Petris <gpet...@uark.edu> Associate Professor Department of Mathematical Sciences University of Arkansas - Fayetteville, AR 72701 Ph: (479) 575-6324, 575-8630 (fax) http://definetti.uark.edu/~gpetris/ ### UNIVERSITY OF CAPE TOWN This e-mail is subject to the UCT ICT policies and e-mai...{{dropped:13}} ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.