# Can you suggest a novice method for detrending and deseasoning time series to find relationship between two variables?

I have a time series with two variables: discount and sales. My goal is to figure out how the discount affects the sales. I want to perform regression. However, the sales data has a trend and seasonality. What could you suggest for removing the trend and the seasonality given that I'm a complete novice? First and foremost, my goal is to learn.

For clarity, I plotted sales and discount (normalized) in time:

• You might find this useful en.wikipedia.org/wiki/Cointegration Oct 21 at 8:11
• A combination of long-term trend modeling with a global regression spline and cyclic variation modeling with a harmonic function can be quite handy. See for example hbiostat.org/R/rms/gtrans.html Oct 21 at 11:43
• I second splines with harmonics and feed the exogenous to it as well. Since your end goal is to run a regression on the series with trend and seasonality accounted for you may find that something like Holt-Winters will eat up a lot of that signal as 'trend'. So your regression results with the other series might not be as useful. Since your goal is to learn you could do this with Holt-Winters and see the potential issues! I would then hand everything over to the package prophet and see what you get in time series cross validation results. Oct 21 at 13:18
• Thanks to all of you! Many unfamilliar words that I will have to learn :) Oct 24 at 5:17

I wanted to expand on my comment:

So let's simulate a series which exhibits a basic trend, seasonality, and reacts to an exogenous variable. (I am working in python)

import numpy as np
import matplotlib.pyplot as plt

seasonality = ((np.cos(np.arange(1, 101))*10 + 50))
np.random.seed(100)
true = np.linspace(-1, 100, 100)
noise = np.random.normal(0, 1, 100)
exo = np.random.randint(low=0, high=2, size=len(true))
y = true + seasonality + exo*100 + noise
plt.plot(y)
plt.show()


Here we have a basic trend line with a seasonal pulse of a cosine bouncing between -10 and 10 along with an exogenous factor which is either 1 or 0. If it is 1 then it impacts our series by +100 else 0. There is noise but it is fairly negligible which makes the further outcome even more concerning!

Now, one way to handle what you want would be to do time series analysis to de-trend and de-seasonalize your series and hand off the residuals to a linear regression. However, special care needs to be taken if you do that and here is why. I will be using a dev package ThymeBoost just because it is easy to switch between different trend/seasonal methods and I am developing it so I am biased.

from ThymeBoost import ThymeBoost as tb

boosted_model = ThymeBoost()

output = boosted_model.fit(y,
trend_estimator='linear',
seasonal_estimator='fourier',
seasonal_period=25,
split_cost='mse',
global_cost='maicc',
fit_type='global')


This will do simple decomposition with a simple trend line and harmonic seasonality. It also will do gradient boosting on the decomposition procedure but that isn't important now.

Let's take a look at the outcome:

and the components:

You can see that the seasonality signal is not quite right, in fact, it is eating up some of the exogenous.

If we try a smoother, let's use double exponential smoothing, we will see even more of the signal being eaten up.

boosted_model = ThymeBoost()

output = boosted_model.fit(y,
trend_estimator='des',
seasonal_estimator='fourier',
seasonal_period=25,
split_cost='mse',
global_cost='maicc',
fit_type='global',
)
boosted_model.plot_components(output)


So now that we are using a more complex trend/level estimator it is overfitting.

Now, let's try adding exogenous, normally in standard decomposition procedure we would just take the previous output so our exogenous regression would be quite wrong. But if the model will try to fit for all components at once such as the prophet package or what ThymeBoost emulates with boosting we can get a better approximation.

boosted_model = ThymeBoost()

output = boosted_model.fit(y,
trend_estimator='linear',
exogenous_estimator='ols',
seasonal_estimator='fourier',
seasonal_period=25,
split_cost='mse',
global_cost='maicc',
fit_type='global',
exogenous=exo,
)

boosted_model.plot_results(output)


And the components:

Alright much better.

This was a little long winded and a very naive example but hopefully it illustrated how traditional decomposition of approximating trend then seasonality and finally exogenous can be quite flawed.

Like I said in the comment, I would try prophet or some other GAM like model that will fit everything at once and is well maintained.

I would recommend looking at Holt-Winters seasonality (explained here).

It decomposes the time series in to three terms: level, trend and seasonality.

• And if the OP is using R, this is easy to implement (cannot say the same for Python, etc as I don't use those). Oct 22 at 1:21