Seasonal plot for time series with multiple seasonalities

Question

My data is like the following, half hourly multi-seasonal time series from 2011 to 2016.

HH_ID         Date_ID   Demand   Year_ID  MM_ID  day
201101010000  20110101  35090     2011   201101  Monday
201101010030  20110101  35612     2011   201101  Monday
201101010100  20110101  35245     2011   201101  Monday
201101010130  20110101  34161     2011   201101  Monday
...            ...       ...      ...    ....
demand.ts1<-msts(demand$Demand,seasonal.period=c(48,7*48,365*48),start=2011)

I would like to use season plot or whatever to plot it and explore its seasonality. The current one I got is quite ugly. Is there any recommendation how to get a better plot?

I also randomly manually chose a week from each year, plotted them.

The question is how to let xlim to show the day (Mon,Tues...), rather than the index of data? Thanks.

Answer 1

Let's simulate some data with your data dimensions:

library(forecast)
set.seed(1)
demand <- msts(rnorm(5*365*48),seasonal.period=c(48,7*48,365*48),start=2011)

I would suggest you plot intraday curves first. For this, simply recode your time series to have a single seasonality of 48, then call seasonplot():

seasonplot(ts(demand,frequency=48),col=rainbow(length(demand)/48))

Next, aggregate demands within days and plot their seasonality within weeks:

demand.daily <- ts(as.vector(by(demand,rep(1:(length(demand)/48),each=48),sum)),frequency=7)
seasonplot(demand.daily,col=rainbow(length(demand.daily)/7))

Finally, you can plot either daily or weekly aggregates and their seasonality over the year. The problem here is that a year has 365.25 days and 52.18 weeks, so you will need to drop a few data points here and there.

These graphics may already give you some ideas. However, you still have a lot of data (5 years are 1800 days, so the first plot has 1800 lines, and that's not really easy to look at). So one alternative might be to plot, e.g., beanplots of half-hourly data, then add lines connecting, say, the 10%, 25%, 50%, 75% and 90% quantiles.

library(beanplot)
beanplot(demand~rep(1:48,length(demand)/48),what=c(0,1,0,0),col="lightgray")
probs <- c(.1,.25,.5,.75,.9)
lines(rbind(matrix(1:48,nrow=48,ncol=length(probs)),NA),
    t(cbind(matrix(unlist(by(demand,rep(1:48,length(demand)/48),quantile,probs=probs)),
      nrow=length(probs)),NA)))

Finally, in your plot, you can suppress the x axis annotation by using the xaxt="n" parameter in your plotting command, then add a custom axis using the axis() command. (I wouldn't add all "Monday", "Tuesday" etc. x axis labels for a full year, though.)