The Mathematics Behind Time Series ML

The mathematics under the forecast: why a thousand correlated observations are worth fifty, what stationarity actually buys you, how two unrelated random walks look significant, and why shuffled folds lie.

Every technique in the preceding articles rested on an assumption so quiet that it is rarely stated aloud. Points are drawn independently from a fixed distribution. That assumption is what licenses the train-test split, what makes cross-validation an estimator of anything, what turns the variance of a mean into $\sigma^2/n$, and what allows a model fitted on some rows to say something about other rows. Remove it and the machinery does not degrade gracefully. It breaks, silently, while continuing to print numbers.

Time series removes it. Yesterday’s temperature is not independent of today’s; this quarter’s revenue is not a fresh draw unrelated to last quarter’s. The dependence is not a nuisance to be washed out, it is the entire signal, since a series with no dependence across time is white noise and cannot be forecast at all. So the mathematics of time series is the mathematics of what to do when the independence assumption is both false and load-bearing: how much information a correlated sample really contains, what conditions make estimation possible at all, how to convert a sequence into something a regressor can eat, and why the validation scheme from every other article will lie to you here.

The cost of dependence, quantified

Begin with the simplest possible estimate, the sample mean, and ask what correlation does to it. For independent observations the result is familiar:

Var(x)=σ2n\operatorname{Var}(\bar{x}) = \frac{\sigma^2}{n}

For a dependent sequence, the variance of a sum picks up every pairwise covariance, and the general expression is:

Var(x)=σ2n[1+2k=1n1(1kn)ρk]\operatorname{Var}(\bar{x}) = \frac{\sigma^2}{n}\left[1 + 2\sum_{k=1}^{n-1}\left(1 – \frac{k}{n}\right)\rho_k\right]

where ρk\rho_k is the correlation between observations k steps apart. The bracket is the damage. When all ρk\rho_k are zero it equals 1 and we recover the familiar formula; when the ρk\rho_k are positive, as they overwhelmingly are in real series, the bracket exceeds 1 and the variance of our estimate is larger than the independent calculation claims.

Make it concrete with the workhorse model, the first-order autoregression, whose correlations decay geometrically as ρk=ϕk\rho_k = \phi^k. Summing the geometric series and letting n grow:

Var(x)σ2n1+ϕ1ϕ\operatorname{Var}(\bar{x}) \approx \frac{\sigma^2}{n} \cdot \frac{1 + \phi}{1 – \phi}

Comparing this against σ2/neff\sigma^2/n_{\text{eff}} defines the effective sample size, the number of independent observations your correlated sample is actually worth:

neff=n1ϕ1+ϕn_{\text{eff}} = n \cdot \frac{1 – \phi}{1 + \phi}

Put a realistic number in. A daily series with ϕ=0.9\phi = 0.9, which is unremarkable for temperatures, prices, or traffic, gives:

neff=1000×0.11.953n_{\text{eff}} = 1000 \times \frac{0.1}{1.9} \approx 53

A thousand observations carrying the statistical weight of fifty-three. Every standard error computed as if the rows were independent is understated by a factor of roughly $

194.4\sqrt{19} \approx 4.4, and every confidence interval is four times too narrow. This is the first and most general lesson: a long time series is a small sample in disguise, and the more persistent the series, the smaller it is.

Autocovariance, and the condition that makes estimation possible

To go further we need to name the dependence structure. Define the autocovariance at lag k:

γ(k)=Cov(xt,xt+k)=𝔼[(xtμ)(xt+kμ)] \gamma(k) = \operatorname{Cov}(x_t, x{t+k}) = \mathbb{E}\left[(x_t – \mu)(x{t+k} – \mu)\right]

and normalise it by the variance to get the autocorrelation function:

ρ(k)=γ(k)γ(0)\rho(k) = \frac{\gamma(k)}{\gamma(0)}

The ACF is the series’ fingerprint. It says how much of what happened k steps ago is still visible now, and its shape distinguishes model families: geometric decay points to autoregression, an abrupt cutoff after lag q points to a moving average, and periodic spikes announce seasonality.

But there is a problem hiding in the definition, and it is the deepest issue in the subject. The expectation 𝔼[(xtμ)(xt+kμ)]\mathbb{E}[(x_t – \mu)(x_{t+k} – \mu)] is over repeated realisations of the process. We have one realisation. History ran once. There is no ensemble of parallel timelines to average over, so the expectation is not directly estimable from what we hold.

Stationarity is the assumption that rescues us. A series is weakly stationary when three conditions hold:

𝔼[xt]=μfor all t\mathbb{E}[x_t] = \mu \quad \text{for all } t
Var(xt)=γ(0)for all t\operatorname{Var}(x_t) = \gamma(0) \quad \text{for all } t
Cov(xt,xt+k)=γ(k)depending on k only, not on t \operatorname{Cov}(x_t, x_{t+k}) = \gamma(k) \quad \text{depending on } k \text{ only, not on } t

Each condition says the same thing in a different place: the process’s statistical properties do not drift. And that is exactly what licenses the estimator everyone actually computes:

γ^(k)=1nt=1nk(xtx)(xt+kx) \hat{\gamma}(k) = \frac{1}{n}\sum_{t=1}^{n-k}\left(x_t – \bar{x}\right)\left(x_{t+k} – \bar{x}\right)

Look at what this does. It averages across time and calls the result an estimate of an average across realisations. That substitution is only legitimate if the quantity being averaged is the same at every t, which is precisely the stationarity condition, and the formal property that time averages converge to ensemble averages is ergodicity. Without stationarity, the sum above is averaging different quantities at different dates and estimates nothing.

This reframes stationarity from a box-ticking exercise into the enabling condition of the entire field. It is not a technicality demanded by a test. It is what makes one run of history informative about the process that generated it.

Random walks, and the regression that fools everyone

The canonical non-stationary process is the random walk:

xt=xt1+ϵtx_t = x_{t-1} + \epsilon_t

Iterate it back to the start and the whole series is a cumulative sum of shocks:

xt=x0+s=1tϵsx_t = x_0 + \sum_{s=1}^{t}\epsilon_s

so its variance is:

Var(xt)=t,σ2\operatorname{Var}(x_t) = t,\sigma^2

The variance grows without bound as time passes, violating the second stationarity condition outright. The process has no mean to revert to and no fixed scale; shocks are permanent rather than transient. Note the boundary case in the AR(1) family: xt=ϕxt1+ϵtx_t = \phi x_{t-1} + \epsilon_t is stationary when |ϕ|<1|\phi| < 1 and becomes a random walk at exactly ϕ=1\phi = 1, which is why the unit root is the dividing line and why tests for it matter.

The remedy is differencing, which recovers stationarity in one step:

Δxt=xtxt1=ϵt\Delta x_t = x_t – x_{t-1} = \epsilon_t

Now the danger. Take two completely independent random walks, generated from unrelated noise, and regress one on the other. The regression reports a large R2R^2, a coefficient far from zero, and a t-statistic that would be overwhelming evidence in any textbook. There is no relationship whatsoever. This is spurious regression, and its most disturbing feature is asymptotic: as n grows the t-statistic does not settle down, it diverges. More data makes the illusion stronger, not weaker.

The mechanism is the effective sample size argument taken to its limit. With ϕ=1\phi = 1 the factor (1ϕ)/(1+ϕ)(1-\phi)/(1+\phi) is zero, so neff=0n_{\text{eff}} = 0 regardless of n. The regression’s standard errors are computed as though there were n independent observations when there are effectively none, and everything downstream, the t-statistic, the p-value, the confidence interval, is meaningless. Two trending series will appear related simply because both are trending, and difference them before regressing or the model is measuring the calendar.

Lag features: how a series becomes a supervised problem

With stationarity secured, the bridge to ordinary machine learning is short. An autoregression of order p is:

xt=c+i=1pϕi,xti+ϵtx_t = c + \sum_{i=1}^{p}\phi_i , x_{t-i} + \epsilon_t

Read it as a regression, because that is what it is. Assemble a design matrix whose row at time t holds the p previous values, and a target vector holding the current one:

Xt=[xt1,;xt2,;,;xtp],yt=xtX_t = \left[x_{t-1}, ; x_{t-2}, ; \ldots, ; x_{t-p}\right], \qquad y_t = x_t

Fit that with ordinary least squares and you have recovered AR(p) exactly. Fit it with a gradient booster instead and you have a nonlinear autoregression. This is the whole conceptual content of “machine learning for time series”: the temporal structure lives in the feature engineering, and once the lags are built, any regressor from any other article applies unchanged. The sequence enters the model through the design matrix and nowhere else.

Which lags to include is answerable from the ACF rather than by guesswork, and for AR(1) the ACF is exactly the geometric decay we assumed earlier:

ρ(k)=ϕk\rho(k) = \phi^k

The general stationarity condition for AR(p) is that the roots of the characteristic polynomial lie outside the unit circle, which for p = 1 reduces to |ϕ|<1|\phi| < 1, agreeing with everything above.

Rolling statistics are convolutions, and direction matters

The other standard feature is the rolling mean over a window of width w:

mt=1wj=0w1xtjm_t = \frac{1}{w}\sum_{j=0}^{w-1} x_{t-j}

That expression should look familiar from the CNN article, and it is not an analogy. A rolling mean is a convolution of the series with the kernel [1/w,,1/w][1/w, \ldots, 1/w], the same operation a convolutional layer performs, applied along time rather than space, with a fixed kernel rather than a learned one. The reason a one-dimensional CNN works on sequences is that its layers are computing learned versions of exactly this.

The index in the sum carries the most important practical constraint in the field. It runs xtjx_{t-j}, backwards only. A rolling statistic centred on t would include xt+1x_{t+1} and later values, which do not exist at the moment of prediction. In pandas, .rolling(window=w, center=True) produces exactly that, a feature computed from the future, and any model using it will validate beautifully and fail in production, because at deployment the future is unavailable. Causality is a constraint on the summation limits, and it is violated by a keyword argument.

Exponential weighting replaces the flat kernel with a geometric one:

mt=αxt+(1α)mt1m_t = \alpha x_t + (1 – \alpha) m_{t-1}

which unrolls into an infinite weighted sum over the whole past:

mt=αj=0(1α)j,xtjm_t = \alpha \sum_{j=0}^{\infty}(1 – \alpha)^j , x_{t-j}

Every past value contributes, with weight decaying geometrically, so there is no arbitrary cutoff. The natural way to state the memory is the half-life, the lag at which a value’s weight has halved:

half-life=ln0.5ln(1α) \text{half-life} = \frac{\ln 0.5}{\ln(1 – \alpha)}

Distributed lags and dynamic causal effects

Forecasting asks what comes next; causal work asks what an intervention does, and when. The distributed lag model answers the second by letting an input act over multiple horizons:

yt=α+j=0qβj,xtj+ut y_t = \alpha + \sum_{j=0}^{q}\beta_j , x_{t-j} + u_t

Each βj\beta_j is the dynamic causal effect at horizon j: the change in y today attributable to a one-unit change in x j periods ago. The sequence β0,β1,,βq{\beta_0, \beta_1, \ldots, \beta_q} is the impulse response function, and it is the shape of the answer to “when does this work?”. A marketing spend with β0\beta_0 small and β2\beta_2 large is a campaign that takes two periods to bite.

Partial sums answer the cumulative question:

cumulative multiplier(h)=j=0hβj\text{cumulative multiplier}(h) = \sum_{j=0}^{h}\beta_j

and the total effect once everything has played out is the long-run multiplier:

long-run multiplier=j=0qβj\text{long-run multiplier} = \sum_{j=0}^{q}\beta_j

Here the opening variance formula returns to collect its debt. The error term utu_t in a distributed lag model is almost always autocorrelated, so the classical OLS standard errors, which assume independent errors, are wrong in exactly the way the effective sample size predicted: too small, by a factor governed by the persistence of utu_t. The standard repair is a heteroskedasticity and autocorrelation consistent estimator, of which Newey-West is the common choice, which corrects the variance by including the covariance terms the classical formula omits. The coefficients are fine; it is their standard errors that lie.

Validation, and the leak that k-fold guarantees

Every other article ends its validation section with cross-validation. Here, cross-validation is a bug.

Standard k-fold shuffles rows into folds, which means a fold used for training will contain observations from after the observations in the test fold. The model is fitted on the future and evaluated on the past. Since neighbouring points are highly correlated, the training set effectively contains near-copies of the test points, and the resulting score measures interpolation between known values rather than extrapolation into unknown ones. It will be excellent, and it will be meaningless, because the deployed model will never have that luxury.

The valid scheme respects the arrow of time. Forward chaining trains on a prefix and tests on what follows:

train:1,,ttest:t+1,,t+h\text{train}: {1, \ldots, t} \quad \longrightarrow \quad \text{test}: {t+1, \ldots, t+h}

repeated for increasing t, which is what TimeSeriesSplit implements. Every evaluation then answers the question deployment will ask: given everything up to now, what happens next?

The correlation penalty applies here too, and it is worth being blunt about the arithmetic. A forward-chaining scheme with five splits on a persistent series is not five independent assessments; the training sets are nested, so the fold scores are correlated for the same reason k-fold’s were in the validation article, and the spread across splits understates the true uncertainty. Combine that with neff=n(1ϕ)/(1+ϕ)n_{\text{eff}} = n(1-\phi)/(1+\phi) and the honest summary of most time series projects is that they have far less evidence than the row count suggests.

What the mathematics tells you to do

Test for stationarity before anything else, because it is the assumption that makes estimating γ(k)\gamma(k) from a single history legitimate, and without it the ACF and every model built on it estimate nothing. Difference trending series before regressing them on each other, because two independent random walks produce a significant relationship whose t-statistic grows with the sample. Build lag features and rolling statistics with backward-looking windows only, since the summation limits are the mathematical statement of causality and center=True breaks it. Read the ACF to choose lags rather than guessing, because ρ(k)=ϕk\rho(k) = \phi^k tells you where the information runs out. Use HAC standard errors on distributed lag models, because the coefficients are consistent but their classical variances are not. Validate by forward chaining, never by shuffling, since shuffled folds train on the future. And divide your confidence by roughly (1+ϕ)/(1ϕ)\sqrt{(1+\phi)/(1-\phi)} before believing any of it.

The single idea underneath all of this is that dependence is simultaneously the signal and the tax. It is the signal because a series with no autocorrelation is unforecastable, so everything a model can learn is stored in the ρ(k)\rho(k). It is the tax because those same correlations shrink the sample from n to neffn_{\text{eff}}, invalidate the classical standard errors, and make the ordinary validation toolkit hazardous. Practitioners who arrive from tabular machine learning tend to bring the tools and leave the assumption behind, and the resulting models validate at 0.95 and forecast at nothing. The mathematics is telling you the same thing in six different voices: in time series, the rows are not evidence in the way you are used to, and the arrow of time is not a formality.

See you soon.

View Comments (2)

Leave a Reply

Subscribe to My Newsletter

Subscribe to my email newsletter to get the latest posts delivered right to your email. Pure inspiration, zero spam.

Discover more from Discuss Data Science, Machine Learning and Analytics

Subscribe now to keep reading and get access to the full archive.

Continue reading