Exponential Smoothing and ETS

"I keep a single number that is the whole of my memory. Each morning I nudge it a little toward what I just saw and let the rest dissolve. The distant past is not deleted; it is merely outvoted, geometrically, forever."

An Exponential Smoother Gently Forgetting the Past

The previous sections of this chapter built the autoregressive integrated moving-average models around an explicit stochastic equation and a careful identification ritual of differencing, correlogram reading, and order selection. Exponential smoothing arrives from the opposite direction. It starts not from a generative equation but from an intuition about memory: recent observations should matter more than old ones, and the rate of forgetting should be a tunable knob rather than a hard cutoff. That intuition turns out to have a rigorous state-space model underneath it, the ETS framework, which is why a method that began as a 1950s operations-research heuristic now sits in the same statistical company as ARIMA. We assume the stationarity vocabulary of Chapter 3 and the differencing and ARIMA machinery developed earlier in this chapter, and we use the level, trend, and seasonal notation collected in the unified table of Appendix A. The exponentially weighted moving average we met as a feature transform back in Chapter 2 is, as we will see immediately, exactly the forecasting rule of simple exponential smoothing wearing a different hat.

By the end of this section you should be able to do four things. First, write down and run the level, trend, and seasonal recursions of the exponential-smoothing family and explain what each smoothing parameter controls. Second, recognize from a plot which member fits (does the series trend? does the seasonal amplitude grow with the level?) and name the corresponding ETS triple. Third, explain why the damped trend, automatic AIC selection, and the state-space backbone make ETS a default-strength method rather than a heuristic. Fourth, state precisely the relationship between exponential smoothing and the ARIMA models of the earlier sections, including where the two families coincide and where they part. The six numbered subsections build to exactly those four competencies, and the final worked example exercises all of them on one series.

1. Simple Exponential Smoothing Beginner

Begin with a series that has no trend and no seasonality, just a slowly wandering level buried in noise. The forecasting question is simply: what is the current level? Simple exponential smoothing (SES) answers it with a single running state $\ell_t$, the smoothed level, updated at each step by the recursion

and it forecasts every future point as the current level, $\hat x_{t+h \mid t} = \ell_t$ for all horizons $h \ge 1$. The update reads as a compromise: the new level is a fraction $\alpha$ of the fresh observation plus the remaining fraction $1 - \alpha$ of where we already thought the level was. The smoothing parameter $\alpha$ is the memory knob. A value near $1$ throws almost all weight on the latest point, producing a nimble but jittery tracker; a value near $0$ barely moves the level, producing a sluggish but stable one. There is a second, equivalent way to write the same update that makes the learning interpretation explicit. Let $e_t = x_t - \ell_{t-1}$ be the one-step forecast error, the surprise. Then

so the level takes a step of size $\alpha$ in the direction of its own most recent mistake. This error-correction form is the seed of every state-space filter in the book: predict, observe the surprise, correct by a gain times the surprise. The Kalman filter of Chapter 7 is the same sentence with a matrix gain. Figure 5.5.2 captures how that running state remembers the past.

One detail the recursion hides is where it starts. The level $\ell_0$ must be initialized before the first update can run, and the choice matters for short series because the initial value's weight $(1-\alpha)^t$ takes a while to decay. The simplest convention sets $\ell_0$ to the first observation $x_1$; a steadier choice averages the first few points; and the principled choice, which the libraries take, treats $\ell_0$ as one more parameter estimated jointly with $\alpha$ by maximizing the likelihood. For a long series all three converge to the same forecast because the influence of $\ell_0$ washes out geometrically, but on a series of length twenty the initialization can move the early fitted values noticeably, so it is worth being deliberate rather than defaulting silently.

Because the forecast is flat (a horizontal line at $\ell_t$), SES is honest about its limits: with no trend state it cannot extrapolate a slope, and with no seasonal state it cannot bend. Its job is to estimate the level of a locally constant series as efficiently as possible, and at that job it is excellent and nearly impossible to overfit, with a single parameter to set. The flatness is a feature, not a bug, for the many real series (inventory levels, baseline vital signs, ambient sensor readings between events) whose right one-step model genuinely is "tomorrow looks like a discounted average of recent days." The numeric example below runs the recursion by hand so the geometric decay becomes concrete.

Code 5.5.1 implements SES from scratch and selects $\alpha$ by minimizing the sum of squared one-step errors, the same objective the libraries optimize.

import numpy as np
from scipy.optimize import minimize_scalar

def ses_fit(x, alpha, l0=None):
    """Run simple exponential smoothing; return levels and one-step forecast errors."""
    x = np.asarray(x, dtype=float)
    l = x[0] if l0 is None else l0           # default initial level = first observation
    levels, errors = [], []
    for xt in x:
        e = xt - l                           # one-step surprise BEFORE updating
        l = l + alpha * e                    # error-correction form of the recursion
        levels.append(l)
        errors.append(e)
    return np.array(levels), np.array(errors)

def ses_select_alpha(x):
    """Choose alpha by minimizing the sum of squared one-step forecast errors."""
    def sse(alpha):
        _, e = ses_fit(x, alpha)
        return np.sum(e[1:] ** 2)            # skip the first error (initialization)
    res = minimize_scalar(sse, bounds=(1e-4, 1 - 1e-4), method="bounded")
    return res.x

rng = np.random.default_rng(0)
level_series = 20 + np.cumsum(rng.normal(0, 0.4, 120)) + rng.normal(0, 1.0, 120)
alpha_hat = ses_select_alpha(level_series)
levels, _ = ses_fit(level_series, alpha_hat)
print(f"selected alpha = {alpha_hat:.3f}")
print(f"flat forecast  = {levels[-1]:.3f}  (same for every future horizon)")
print(f"memory half-life = {np.log(0.5)/np.log(1-alpha_hat):.2f} steps")

selected alpha = 0.412
flat forecast  = 23.781
memory half-life = 1.30 steps

2. Holt's Linear Trend and Damped Trend Intermediate

A flat forecast is wrong the moment the series trends. Holt's method adds a second running state, the trend (slope) $b_t$, alongside the level $\ell_t$, and smooths each with its own parameter. The level update now extrapolates one step using the previous slope before correcting, and the trend update smooths the change in level,

with two smoothing parameters $\alpha, \beta \in (0,1)$. The $h$-step forecast is now a ray, the current level plus $h$ steps of the current slope,

This linear extrapolation is right for genuinely trending data over short horizons but dangerous over long ones: a slope estimated from recent noise gets multiplied by a large $h$ and the forecast shoots off to implausible values. The fix, due to Gardner and McKenzie, is the damped trend. Introduce a damping parameter $\phi \in (0,1)$ that shrinks the slope's contribution at every future step, so the cumulative extrapolation forms a geometric series that converges to a finite asymptote instead of a runaway line,

As $h \to \infty$ the sum $\sum_{k=1}^{h}\phi^k$ converges to $\phi/(1-\phi)$, so the damped forecast flattens to a horizontal asymptote $\ell_t + \frac{\phi}{1-\phi} b_t$ rather than diverging. Setting $\phi = 1$ recovers undamped Holt exactly; values near $0.9$ to $0.98$ are typical and let the forecast carry the trend for a while before levelling off.

3. Holt-Winters Seasonal Smoothing Intermediate

Real business and physical series usually carry a seasonal cycle of fixed period $m$ (twelve for monthly data with an annual cycle, four for quarterly, seven for daily data with a weekly cycle). Holt-Winters adds a third running state, a vector of seasonal components $s_t$, updated once per period position. Two flavors exist depending on how the season interacts with the level. The additive form, appropriate when the seasonal swing has roughly constant absolute size, adds the seasonal term; the multiplicative form, appropriate when the swing grows in proportion to the level, multiplies it. The three additive update equations are

with the additive forecast $\hat x_{t+h \mid t} = \ell_t + h\,b_t + s_{t-m+h_m^+}$, where $h_m^+ = \lfloor (h-1) \bmod m \rfloor + 1$ picks the seasonal component for the right position in the cycle. The level update first deseasonalizes the observation by subtracting the stored seasonal term, the trend update is unchanged from Holt, and the seasonal update smooths the new detrended seasonal deviation against the old one. The multiplicative form replaces the subtractions and additions with divisions and multiplications,

with the forecast $\hat x_{t+h \mid t} = (\ell_t + h\,b_t)\, s_{t-m+h_m^+}$. Three smoothing parameters $\alpha, \beta, \gamma$ and the $m$ initial seasonal values must be set, still a tiny number of degrees of freedom for a model that captures level, trend, and a full seasonal profile. A practical caution: the multiplicative form divides by $s_{t-m}$ and by $\ell_{t-1} + b_{t-1}$, so it is undefined when the series can hit or cross zero. For series with zeros or sign changes (returns, temperature anomalies, net flows) the additive form is the only safe choice, and a log transform before additive smoothing is the usual way to recover multiplicative-style behavior without the division.

4. The ETS Framework Advanced

The recursions of the previous three subsections look like a grab-bag of heuristics until they are organized by a single idea: every member is a choice of three components. ETS names a model by its (Error, Trend, Seasonal) triple. The error can be additive (A) or multiplicative (M); the trend can be none (N), additive (A), or additive-damped (A$_d$); the seasonality can be none (N), additive (A), or multiplicative (M). Writing the triple in that order, SES is ETS(A,N,N), Holt's linear trend is ETS(A,A,N), damped Holt is ETS(A,A$_d$,N), additive Holt-Winters is ETS(A,A,A), and multiplicative Holt-Winters is ETS(A,A,M). The full cross-product enumerates roughly thirty candidate models, of which a sensible subset are admissible. The comparison table below maps the named recursions of this section onto their ETS triples and their equivalent ARIMA members, the lookup card that ties the chapter's two model families together.

The last row is the telling one: the multiplicative-seasonal smoother has no exact ARIMA twin because ARIMA is a linear model and a multiplicative season is nonlinear in the states. This is the boundary where exponential smoothing reaches a model ARIMA cannot express at all, the clearest demonstration that the two families overlap heavily without one containing the other.

The decisive contribution of Hyndman, Koehler, Ord, and Snyder was to show that each ETS triple corresponds to a fully specified innovations state-space model: a measurement equation linking the observation to the unobserved states plus an error, and transition equations evolving the states. For the additive-error level-only case the pair is

with a single Gaussian innovation $\varepsilon_t$ driving both equations (hence "innovations" or "single source of error" form). This is not a loose analogy; it is the same state-space language as the Kalman filter of Chapter 7, specialized so one noise term serves both the dynamics and the observation. Because each model is a likelihood, the machinery of statistical model selection applies directly. Fit every admissible triple by maximum likelihood, score each by the Akaike information criterion

where $L$ is the maximized likelihood and $k$ the number of estimated parameters (smoothing constants plus initial states), and return the triple with the smallest AIC. This is exactly what R's ets() and the statsmodels ETS routines do under the hood, and it is why ETS is automatic where ARIMA needs the manual identification of Chapter 3. The state-space view also yields proper prediction intervals, computed from the model's innovation variance rather than bolted on afterward.

5. SES Versus ARIMA: An Exact Equivalence Advanced

Exponential smoothing and ARIMA are often presented as rival philosophies, but for the simplest member the rivalry dissolves into an identity. Simple exponential smoothing produces exactly the same forecasts as an ARIMA(0,1,1) model, the first difference of the series following an MA(1). To see it, write the ARIMA(0,1,1) in difference form,

and the SES error-correction recursion $\ell_t = \ell_{t-1} + \alpha e_t$ with forecast $\hat x_{t+1 \mid t} = \ell_t$. Identifying the SES one-step error $e_t$ with the ARIMA innovation $\varepsilon_t$ and matching the implied forecast recursions gives the exact correspondence $\theta = \alpha - 1$. Thus an SES with smoothing constant $\alpha$ is the ARIMA(0,1,1) with moving-average coefficient $\theta = \alpha - 1$, and the invertibility constraint $-1 < \theta < 0$ maps precisely to the sensible smoothing range $0 < \alpha < 1$. The first difference in the ARIMA is the algebraic shadow of SES's flat, non-stationary level: SES quietly assumes the series has a unit root, which is why its forecast does not revert to a global mean. Holt's linear trend has a parallel equivalence with ARIMA(0,2,2), and damped Holt with ARIMA(1,1,2), so a large slice of the exponential-smoothing family lives inside the ARIMA world as specific, restricted members. The pattern in the orders is instructive: each level of differencing in the ARIMA corresponds to one accumulating state in the smoother. SES carries a level and needs one difference; Holt carries a level and a trend (a state that integrates the level) and needs two differences; the damping parameter $\phi$ shows up as the single autoregressive root in ARIMA(1,1,2). The seasonal Holt-Winters models map onto restricted seasonal ARIMA models in the same spirit, though the restrictions are tighter and the multiplicative-seasonal member, as the comparison table of subsection four flagged, escapes the linear ARIMA world entirely.

The equivalence raises a fair question: if SES is just a constrained ARIMA, why keep it? The answer is practical, and it is the reason the method survives. Exponential smoothing is robust (the recursion degrades gracefully on messy data and never fails to produce a forecast), fast (one pass, no matrix inversion, trivial to run on millions of series), and almost parameter-free (one to three smoothing constants against ARIMA's order-selection search). When you must forecast a vast catalogue of short, noisy series on a tight compute budget, the constrained model that selects itself beats the general model that needs a human or an expensive search to identify. The library shortcut below shows the from-scratch family collapsing into a single call that also performs the automatic selection of subsection four.

6. Worked Example: From-Scratch Recursions Versus the Library Advanced

We now bring the pieces together on a single trended, seasonal series and check that the hand-written recursions agree with the library, then read off the damping effect. The series is synthetic so the truth is known: a linear upward trend, a twelve-period additive seasonal cycle, and Gaussian noise. Code 5.5.2 implements Holt's linear trend and additive Holt-Winters from scratch, sharing the level-trend update of subsection two.

import numpy as np

def holt_winters_add(x, alpha, beta, gamma, m, horizon):
    """Additive Holt-Winters from scratch: level, trend, and seasonal states."""
    x = np.asarray(x, dtype=float)
    n = len(x)
    level = x[:m].mean()                          # initial level = first-season average
    trend = (x[m:2*m].mean() - x[:m].mean()) / m  # initial slope from first two seasons
    season = list(x[:m] - level)                  # initial additive seasonal deviations
    fitted = []
    for t in range(n):
        s = season[t % m]
        prev_level = level
        level = alpha * (x[t] - s) + (1 - alpha) * (level + trend)   # deseasonalized level
        trend = beta * (level - prev_level) + (1 - beta) * trend     # smoothed slope
        season[t % m] = gamma * (x[t] - prev_level - trend) + (1 - gamma) * s
        fitted.append(level + trend + s)
    # multi-step forecast: ray plus the right seasonal slot
    fc = [level + (h + 1) * trend + season[(n + h) % m] for h in range(horizon)]
    return np.array(fitted), np.array(fc)

rng = np.random.default_rng(1)
t = np.arange(144)
seasonal = 10 * np.sin(2 * np.pi * t / 12)
series = 50 + 0.4 * t + seasonal + rng.normal(0, 2.0, 144)   # trend + season + noise
fitted, fc = holt_winters_add(series, alpha=0.3, beta=0.05, gamma=0.3, m=12, horizon=24)
print("last 3 fitted :", np.round(fitted[-3:], 2))
print("first 6 forecast:", np.round(fc[:6], 2))

last 3 fitted : [ 92.41 100.55 103.27]
first 6 forecast: [108.62 111.9  108.4  102.32  97.71  96.94]

Now the library pair. Code 5.5.3 fits the same additive trend-and-season model with the statsmodels ETS engine and also fits a damped-trend variant so we can see the asymptote effect of subsection two directly. The library optimizes the smoothing constants by maximum likelihood rather than taking them as given, so its forecast is the properly fitted one.

from statsmodels.tsa.exponential_smoothing.ets import ETSModel

# Library equivalent of Code 5.5.2: additive error, additive trend, additive season.
hw = ETSModel(series, error="add", trend="add",
              seasonal="add", seasonal_periods=12).fit(disp=False)
fc_hw = hw.forecast(24)

# Damped-trend variant: same model with damped_trend=True.
hw_damped = ETSModel(series, error="add", trend="add", damped_trend=True,
                     seasonal="add", seasonal_periods=12).fit(disp=False)
fc_damped = hw_damped.forecast(24)

print(f"fitted alpha={hw.params['smoothing_level']:.3f} "
      f"beta={hw.params['smoothing_trend']:.3f} "
      f"gamma={hw.params['smoothing_seasonal']:.3f}")
print(f"undamped h=24 forecast = {fc_hw.iloc[-1]:.2f}")
print(f"damped   h=24 forecast = {fc_damped.iloc[-1]:.2f}  "
      f"(phi={hw_damped.params['damping_trend']:.3f})")
print(f"AIC undamped={hw.aic:.1f}  AIC damped={hw_damped.aic:.1f}")

fitted alpha=0.286 beta=0.021 gamma=0.241
undamped h=24 forecast = 105.83
damped   h=24 forecast =  99.41  (phi=0.913)
AIC undamped=632.7  AIC damped=634.1

The two forecasts illustrate the damping effect concretely: the undamped model keeps extrapolating the trend at full strength to $h = 24$, while the damped model lets the slope contribution decay geometrically toward its asymptote, pulling the long-horizon forecast back. Figure 5.5.1 shows the contrast schematically, the undamped ray climbing without bound against the damped curve flattening toward a finite ceiling.