Significant autocorrelation in time series decomposition random component

I'm very new to time series analysis. The data below represents about 8 years of aggregate daily visitors to some tourist attractions. I'm trying to examine the random component of some time series data to see if there's anything meaningful in there - i.e. once the trend and seasonal components are removed, is the daily visitor count influenced by factors such as deviations from the seasonal mean weather, promotions, major sporting events (olympics, world cups) etc? However, separating out this random component is proving difficult.

Here's some R code and visualisation of the data:

plot(visitors, col="lightgray", ylim=c(0,700000), ylab="Visitors")
par(new=TRUE);
plot(SMA(visitors, 30), ylim=c(0,700000), type='l', col="black", ylab="Visitors")


#set up timeseries, ignore first year when visitor counts weren't reliable at all sites
visitors = ts(all$visitor[365:length(all$visitor)], start=1, frequency=365.25);
#add a fictional visitor to each zero day so we can take logarithms
visitors[which(visitors==0)]=1;
#plot decomposition
attribs = decompose(visitors);
plot(attribs);


At this point, I observe the periodic spikes in the random component and think "That decomposition wasn't great." So I look at the correlograms of the random component to confirm my suspicion:

layout(c(2,1))
acf(attribs$random, na.action=na.pass); pacf(attribs$random, na.action=na.pass);


So I have significant autocorrelation everywhere. I wonder if maybe the decompose algorithm doesn't work well for these particular data, so I try stl

plot(stl(visitors, s.window="periodic"))


This random component looks even worse. And so, at this point I don't know what to think. Why are the decomposition algorithms putting so much periodocity in the random component? What is the root cause of this issue, and what analytical approach should I adopt to decompose to a truly random remainder?

• Welcome to our site! That is a well-researched question. – whuber Oct 24 '14 at 15:47

So after some research into the algorithms, it seems that neither decompose() nor stl() can handle complex, multiple periods (which otherwise make it through to the remainder component). Like most time series analysis questions, it seems Rob Hyndman has authored the answer to this problem in bats() and tbats() in his forecast package, now I have to learn how to use them properly by wading through his paper that he wrote a few years ago on the algorithm.
• Yes, makes sense, have you had a chance to run bats() or tbats()? I think it would be a very helpful resource if you don't mind sharing your findings on how the forecast package fared. – Zhubarb Oct 27 '14 at 8:19