## Thursday, September 1, 2016

### Seasonal Time Series and Linear Regression

My last post explored an idea to automatically detect outliers on seasonal time series. These time series would have both weekly and daily cycles along with a small trend, such as the ones you will find on some data sources such as server performance meters or even in some business cycles. Data for these time series is typically sampled with a given frequency, like every five minutes or every hour. To detect the existence of unexpected data on such time series, the post followed an idea that instantiates a predictive model for every sampling period of the weekly period. The original idea suggested the use of an exponential moving average model that proved to work very well. I actually tested it with a four-year data set from one of our servers only to find that outlier detection was spot-on. Now, what else can we try with this idea?

What if we used a different model, say a linear regression? In this case we would have a linear regression model for each one of the sampling periods. As you may recall, the exponential moving average model was implemented as an on-line process, being able to read and discard one sampled measure at a time. Can we do this with a linear regression model? Let's see how to estimate the parameters of a simple linear regression model:

$\hat{\alpha}=\bar{y}-\hat{\beta}\bar{x}$

and:

$\hat{\beta}=\frac{\sum_{i=1}^{n}\left&space;(&space;x_{i}-\bar{x}&space;\right&space;)\left&space;(&space;y_{i}-\bar{y}&space;\right&space;)}{\sum_{i=1}^{n}\left&space;(&space;x_{i}-\bar{x}&space;\right&space;)^{2}}$

Looking at the formula to calculate the slope estimate, it might seem that the whole data set would be needed. After all, the above formula seems to imply that we must keep track of all the x and y values, so the whole data set would have to be kept in memory like a spreadsheet.

Fortunately this is not so. After some algebraic manipulations (see here), we can safely transform the equation to this

$\hat{\beta}=\frac{\overline{xy}-\bar{x}\bar{y}}{\overline{x^{2}}-\bar{x}^2}$

where the variables above are very simple averages:

$\bar{x}=\frac{\sum_{i=1}^{n}x_{i}}{n},\bar{y}=\frac{\sum_{i=1}^{n}y_{i}}{n},\overline{x^{2}}=\frac{\sum_{i=1}^{n}x_{i}^{2}}{n},\overline{xy}=\frac{\sum_{i=1}^{n}x_{i}y_{i}}{n}$

This means that we only need a few variables to store all the necessary information to calculate the regression coefficients. Note that instead of the averages themselves, only the sums will be stored that, along with the sample count (n), will allow to perform the average calculation when needed. The stored variables are then:

$n,S_{x}=\sum_{i=1}^{n}x_{i},S_{y}=\sum_{i=1}^{n}y_{i},S_{xx}=\sum_{i=1}^{n}x_{i}^{2},S_{xy}=\sum_{i=1}^{n}x_{i}y_{i}$

For completeness, we will also store Syy (see below why). To determine if a given value is unexpected, we will use the same approach as before: first the expected value is calculated around which an acceptability range is calculated using the expected value's standard error or other similar measure. If the observed value falls outside this range, it is signaled as a potential outlier.

In the case of linear regression there is a standard approach to solve the outlier detection issue. For each estimated value $\small&space;\hat{y}$, calculate its residual - the difference between the expected and the observed value - and divide it by the residuals standard error (RSE): the studentized residual in statistics lingo. The residual is pretty easy to calculate:

$\epsilon=\hat{y}-\hat{\alpha}-\hat{\beta}x$

So how do we calculate RSE? Easy:

$RSE=\sqrt{\frac{RSS}{n-2}}$

Where RSS, the residual sum of squares is defined as:

$RSS=\sum_{i=1}^{n}\left&space;(&space;y_{i}-\hat{y}_{i}&space;\right&space;)^{2}$

We can also rearrange this formula using

$\hat{y}_{i}=\hat{\alpha}+\hat{\beta}x_{i}$

so that we get

$RSS=\sum_{i=1}^{n}\left&space;(&space;y_{i}-\hat{\alpha}-\hat{\beta}x_{i}&space;\right&space;)^{2}$

Expanding the square, we get

$RSS=\hat{\beta}^{2}\sum_{i=1}^{n}x_{i}^{2}-2\hat{\beta}\sum_{i=1}^{n}x_{i}y_{i}+2\hat{\alpha}\hat{\beta}\sum_{i=1}^{n}x_{i}+\sum_{i=1}^{n}y_{i}^{2}-2\hat{\alpha}\sum_{i=1}^{n}y_{i}+n\hat{\alpha}^{2}$
or

$RSS=\hat{\beta}^{2}S_{xx}-2\hat{\beta}S_{xy}+2\hat{\alpha}\hat{\beta}S_{x}+S_{yy}-2\hat{\alpha}S_{y}+n\hat{\alpha}^{2}$

which is already expressed in terms of the variables we store. Now, it is very easy to implement the new LinearRegressionModel class. The update method is simply:

Outlier detection is a bit more convoluted:

As usual, the code can be found on GitHub.

Equations in this post were created using a very nice online equation editor. Besides the Wikipedia article, I also used some knowledge from this wonderful book: An Introduction to Statistical Learning: with Applications in R (Springer Texts in Statistics).