25

I have a time series that contains double seasonal components and I would like to decompose the series into the following time series components (trend, seasonal component 1, seasonal component 2 and irregular component). As far as I know, the STL procedure for decomposing a series in R only allows one seasonal component, so I have tried decomposing the series twice. First, by setting the frequency to be the first seasonal component using the following code:

ser = ts(data, freq=48)
dec_1 = stl(ser, s.window="per")

Then, I decomposed the irregular component of the decomposed series (dec_1) by setting the frequency to be the second seasonal component, such that:

ser2 = ts(dec_1$time.series[,3], freq=336)
dec_2 = stl(ser2, s.window="per")

I'm not very confident with this approach. And I would like to know if there are any other ways to decompose a series that has multiple seasonalities. Also,I have noticed that the tbats() function in the R forecast package allows one to fit a model to a series with multiple seasonalities however, it doesn't say how to decompose a series with it.

Stephan Kolassa
  • 95,027
  • 13
  • 197
  • 357
ace
  • 251
  • 1
  • 3
  • 3
  • Hi there and welcome to the site. For your two seasonal components, do they have different periodicity, e.g. is one weekly and another monthly? – Michelle Mar 24 '12 at 19:57
  • 1
    Chapter 14 of Rob Hyndman, Koehler, Ord & Snyder "Forecasting with Exponential Smoothing" covers this. Hyndman also has a forecasting package in R. I seem to recall Hyndman having posted on this site on this topic, but it might have been on his blog. – zbicyclist Mar 24 '12 at 20:30
  • @Michelle Hi thanks for the reply. Yeah the two seasonal components have different periodicity. The first one has a periodicity of 48 (daily seasonality), while the second has a periodicity of 336 (weekly seasonality). It is a half hourly time series. – ace Mar 24 '12 at 20:30
  • @zbicyclist I believe the forecasting package that you're on about is the 'forecast' package that I mentioned in the original post. I have had a look at the tbats function of this package but it doesnt say how to use it for decomposing. I will check out the book to see if I can find any further illustration. – ace Mar 24 '12 at 20:37
  • 2
    Here's what I was thinking of. It was on Hyndman's blog. http://robjhyndman.com/papers/complex-seasonality/ – zbicyclist Mar 24 '12 at 20:38
  • @zbicyclist Thanks for the link to the article. The components of Figure 5 on Page 28 of the article is exactly what I'm looking for. However, it doesn't say how to implement (the decomposition) in R. – ace Mar 25 '12 at 12:57
  • Does anyone know how to apply the tbats() function in R for decomposition? – ace Mar 26 '12 at 21:40

3 Answers3

12

R's forecast package bats() and tbats() functions can fit BATS and TBATS models to the data. The functions return lists with a class attribute either "bats" or "tbats". One of the elements on this list is a time series of state vectors, $x(t)$, for each time, $t$.

See http://robjhyndman.com/papers/complex-seasonality/ for the formula's and Hyndman et al (2008) for a better description of ETS models. BATS and TBATS are an extension of ETS.

For example:

fit <- bats(myTimeseries)
fit$x

In this case, each row of x will be on fourier-like harmonic.

There are also plot.tbats() and plot.bats() functions to automatically decompose and view the components.

power
  • 1,564
  • 1
  • 16
  • 29
2

R's forecast package now has a function mstl() to handle multiple seasonal time series decomposition.

This page has got more details how to use it: https://pkg.robjhyndman.com/forecast/reference/mstl.html

1

The facebook prophet package supports multiple seasonalities.

Yearly, weekly and daily seasonalities are built-in but custom seasonalities can be specified.

Here is a custom monthly seasonality:

df <- ...     # data to build model on or decompose
future <- ... # data to make forecasts on

m <- prophet(weekly.seasonality = FALSE)
m <- add_seasonality(m, name = 'monthly', period = 30.5, fourier.order = 5)
m <- fit.prophet(m, df)
forecast <- predict(m, future)
prophet_plot_components(m, forecast)

If predict() is called without passing in a data frame, then it will decompose the time series used to build the model.