Trend, seasonality and noise

The classical decomposition model#

The additive model assumes

$$ y_t = T_t + S_t + R_t $$

and the multiplicative model assumes

$$ y_t = T_t \cdot S_t \cdot R_t $$

• $y_t$ — the observed value at time $t$.
• $T_t$ — the trend-cycle component, a smooth low-frequency term.
• $S_t$ — the seasonal component, with period $m$, constrained so that $\sum_{j=1}^{m} S_j = 0$ (additive) or $\prod_{j=1}^{m} S_j = 1$ (multiplicative).
• $R_t$ — the remainder.

The identification constraint on $S_t$ is what makes the split unique. Without it, any constant could be moved between $T_t$ and $S_t$ with no change to $y_t$. Taking logs maps the multiplicative model to the additive one exactly: $\log y_t = \log T_t + \log S_t + \log R_t$, valid only for strictly positive $y_t$.

The classical algorithm and its defects#

The seasonal_decompose procedure is:

1. Estimate $T_t$ by a centred moving average of order $m$. For even $m$ this is the $2 \times m$ MA, with weights $\frac{1}{2m}, \frac{1}{m}, \dots, \frac{1}{m}, \frac{1}{2m}$ over $m+1$ terms, which avoids the half-period phase shift a plain even-order MA introduces.
2. Detrend: $d_t = y_t - T_t$.
3. Average $d_t$ over each seasonal index $j$, then centre so the values sum to zero, giving $\hat{S}_j$.
4. $R_t = y_t - T_t - S_t$.

Four defects follow directly:

• $T_t$ is undefined for the first and last $\lfloor m/2 \rfloor$ observations, so $R_t$ is too. This is fatal for real-time monitoring at the current edge.
• $\hat{S}_j$ is constant across all years. A seasonal shape that evolves is forced into a single average.
• A single outlier propagates into $\hat{S}_j$ for that seasonal index in every year, because step 3 uses a mean.
• Rapid level shifts are smeared by the moving average across $m$ observations.

STL#

Cleveland, Cleveland, McRae and Terpenning (1990), STL: A Seasonal-Trend Decomposition Procedure Based on Loess, fixes the first three. STL runs an inner loop alternating two loess smooths — a cycle-subseries smooth producing $S_t$ and a low-pass smooth producing $T_t$ — inside an outer loop that computes robustness weights

$$ \rho_t = B!\left(\frac{|R_t|}{6\,\mathrm{median}(|R|)}\right), \qquad B(u) = \begin{cases} (1-u^2)^2 & 0 \le u < 1 \ 0 & u \ge 1 \end{cases} $$

where $B$ is Tukey's bisquare. Observations with large remainders are downweighted in the next iteration, so an outlier is absorbed by $R_t$ rather than bent into $S_t$.

The key tuning parameters are seasonal (the loess span for the cycle-subseries smooth, controlling how fast $S_t$ may evolve) and trend (the span for $T_t$). Setting seasonal to a large odd number recovers a near-constant seasonal shape; small values track a changing one at the cost of variance. Both must be odd.

STL is additive and univariate, handles any integer $m$, and has $O(n)$ cost per iteration. It does not handle multiple seasonalities or calendar effects. MSTL (Bandara, Hyndman and Bergmeir, 2021) extends it to multiple periods by iterative extraction, and is available as statsmodels.tsa.seasonal.MSTL.

Official seasonal adjustment#

National statistical agencies use X-13ARIMA-SEATS, maintained by the US Census Bureau, which combines the X-11 moving-average filters with a RegARIMA pre-adjustment step. The pre-adjustment handles trading-day effects, moving holidays, level shifts and additive outliers before decomposition, then forecast-extends the series so the filters are not truncated at the endpoints. The endpoint problem is the reason published seasonally adjusted figures get revised as new data arrives.

For a treatment of the filters and their frequency responses, see Ladiray and Quenneville (2001), Seasonal Adjustment with the X-11 Method.

Multiple and non-integer seasonality#

Half-hourly electricity demand carries a daily period of 48, a weekly period of 336 and an annual period of 365.25 × 48, which is not an integer. Estimating a separate index per seasonal position is infeasible for the annual term.

The standard remedy is a Fourier representation: replace the seasonal indices with

$$ S_t = \sum_{k=1}^{K} \left[ \alpha_k \sin!\left(\frac{2\pi k t}{m}\right) + \beta_k \cos!\left(\frac{2\pi k t}{m}\right) \right] $$

• $m$ — the period, which may be non-integer.
• $K$ — the number of harmonics, chosen by information criterion; $K \le m/2$.
• $\alpha_k, \beta_k$ — coefficients estimated by regression.

This costs $2K$ parameters instead of $m-1$, and handles several periods at once by summing terms with different $m$. It is the mechanism behind TBATS (De Livera, Hyndman and Snyder, 2011) and behind the seasonality module in Prophet.

Reading#

• Cleveland et al. (1990), STL: A Seasonal-Trend Decomposition Procedure Based on Loess, Journal of Official Statistics 6(1).
• Hyndman and Athanasopoulos, Forecasting: Principles and Practice, 3rd ed., chapter 3 — otexts.com/fpp3/decomposition.html.
• De Livera, Hyndman and Snyder (2011), Forecasting time series with complex seasonal patterns using exponential smoothing, JASA 106(496).
• Bandara, Hyndman and Bergmeir (2021), MSTL: A Seasonal-Trend Decomposition Algorithm for Time Series with Multiple Seasonal Patterns — arxiv.org/abs/2107.13462.

What to learn next#

• Stationarity — why most models want the trend and seasonality gone before they start.
• Moving averages — the smoothing step that produced the trend line above.
• Anomaly detection — turning the residual into an alarm.