Dual seasonal cycles in ts object

Question

I want to strip out seasonality from a ts. This particular ts is daily, and has both yearly and weekly seasonal cycles (frequency 365 and 7).

In order to remove both, I have tried conducting stl() on the ts with frequency set to 365, before extracting trend and remainders, and setting the frequency of the new ts to 7, and repeat.

This doesn't seem to be working very well and I am wondering whether it's my approach, or something inherent to the ts which is causing me problems. Can anyone critique my methodology, and perhaps recommend an alternate approach?

Answer 1

There is a very easy way to do it using a TBATS model implemented in the forecast package. Here is an example assuming your data are stored as x:

library(forecast)
x2 <- msts(x, seasonal.periods=c(7,365))
fit <- tbats(x2)
x.sa <- seasadj(fit)

Details of the model are described in De Livera, Hyndman and Snyder (JASA, 2011).

Answer 2

Check if this is useful:
Start and End Values depends on your Data - Change the Frequency values accordingly

splot <- ts(Data1, start=c(2010, 2), end=c(2013, 9), frequency=12)

additive trend, seasonal, and irregular components can be decomposed using the stl() Function

fit <- stl(splot, s.window="period")
monthplot(splot) 
library(forecast)
vi <-seasonplot(splot)

vi should give seperate values for a seasonal indices

Also Check the below one:

splot.stl <- stl(splot,s.window="periodic",na.action=na.contiguous)
    trend <- splot.stl$time.series[,2]
    season <- splot.stl$time.series[,1]
    remainder <- splot - trend - season

Answer 3

An approach that can handle not only seasonal components (cyclically reoccurring events) but also trends (slow shifts in the norm) admirably is stl(), specifically as implemented by Rob J Hyndman.

The decomp function Hyndman gives there (reproduced below) is very helpful for checking for seasonality and then decomposing a time series into seasonal (if one exists), trend, and residual components.

decomp <- function(x,transform=TRUE)
{
  #decomposes time series into seasonal and trend components
  #from http://robjhyndman.com/researchtips/tscharacteristics/
  require(forecast)
  # Transform series
  if(transform & min(x,na.rm=TRUE) >= 0)
  {
    lambda <- BoxCox.lambda(na.contiguous(x))
    x <- BoxCox(x,lambda)
  }
  else
  {
    lambda <- NULL
    transform <- FALSE
  }
  # Seasonal data
  if(frequency(x)>1)
  {
    x.stl <- stl(x,s.window="periodic",na.action=na.contiguous)
    trend <- x.stl$time.series[,2]
    season <- x.stl$time.series[,1]
    remainder <- x - trend - season
  }
  else #Nonseasonal data
  {
    require(mgcv)
    tt <- 1:length(x)
    trend <- rep(NA,length(x))
    trend[!is.na(x)] <- fitted(gam(x ~ s(tt)))
    season <- NULL
    remainder <- x - trend
  }
  return(list(x=x,trend=trend,season=season,remainder=remainder,
    transform=transform,lambda=lambda))
}

As you can see it uses stl() (which uses loess) if there is seasonality and penal­ized regres­sion splines if there is no seasonality.