Stationarity, autocorrelation, ACF and PACF

Before this lesson: lesson 01 (components, moving averages, lag plots). From the probability-and-statistics track: expectation, variance, covariance and correlation; p-values; the t-test of a regression slope. Time: about 80 minutes, plus the notebook.

---

Why this matters #

• Model choice. AR, MA, ARMA and ARIMA models (lessons 05 and 06) assume a stationary series after differencing. The ACF and PACF tell you which model to try.
• Residual checks. After every model in this track, you check that the residuals look like white noise: an ACF inside the band and a large Ljung–Box p-value. If not, the model missed a pattern.
• Features. Lags with large autocorrelation are the most useful inputs for machine-learning forecasters (lesson 08).
• Correct statistics. Autocorrelated data contain less information than the same number of independent points. This matters for confidence intervals, for the effective sample size of MCMC chains, and for regressions between trending series, which can be pure nonsense (spurious regression, section 3.9).

1. Predict first #

2. The picture #

White noise and the random walk #

White noise is the simplest series: values with mean 0, the same variance $\sigma^2$, and no correlation between any two times. It has no memory. A random walk adds up white noise: $y_t = y_{t-1} + \varepsilon_t$. It has perfect memory: every shock stays in the series forever. So the random walk wanders further and further away from where it started.

Interactive figure: loading…

Memory that fades: the AR(1) process #

Between "no memory" and "perfect memory" there is the AR(1) process $y_t = \phi\, y_{t-1} + \varepsilon_t$ with $|\phi| < 1$. Each value keeps a fraction $\phi$ of the previous value and adds a new shock. Old shocks shrink by the factor $\phi$ at every step, so the series always returns towards its mean. Its autocorrelation at lag $k$, the correlation between $y_t$ and $y_{t-k}$, is $\phi^k$: it decays geometrically.

Interactive figure: loading…

Stationary does not mean flat #

The AR(1) path with $\phi = 0.9$ makes long excursions above and below 0. It is still stationary: the process has the same mean, variance and autocorrelations at every time. A single path can look like it has a trend for a while. Stationarity is a property of the process that generates the data, not of one picture. A random walk is different: its variance really grows with $t$.

3. Definitions and main results #

3.1 White noise and the random walk #

Definition (white noise). $\varepsilon_t \sim \mathrm{WN}(0, \sigma^2)$ means $\mathbb E[\varepsilon_t] = 0$, $\operatorname{Var}(\varepsilon_t) = \sigma^2$ and $\operatorname{Cov}(\varepsilon_t, \varepsilon_s) = 0$ for $t \ne s$. Gaussian white noise is the special case $\varepsilon_t \overset{iid}{\sim} \mathcal N(0, \sigma^2)$.

Definition (random walk). $y_t = y_{t-1} + \varepsilon_t$ for $t \ge 1$, with $y_0 = 0$. With a constant $c$, $y_t = c + y_{t-1} + \varepsilon_t$ is a random walk with drift.

Result 1 (variance of a random walk). $\mathbb E[y_t] = 0$ and $\operatorname{Var}(y_t) = t\sigma^2$. With drift, $\mathbb E[y_t] = ct$ and the variance is the same.

Proof.

1. Unroll the recursion: $y_t = \varepsilon_1 + \varepsilon_2 + \dots + \varepsilon_t$ (plus $ct$ with drift).
2. The mean of a sum is the sum of the means: $\mathbb E[y_t] = 0$ (or $ct$).
3. The $\varepsilon_i$ are uncorrelated, so the variance of the sum is the sum of the variances: $t\sigma^2$. $\blacksquare$

The same steps give $\operatorname{Cov}(y_s, y_t) = \min(s, t)\,\sigma^2$: the covariance depends on the times themselves, not only on their distance.

3.2 Stationarity #

Definition. A process is strictly stationary if the joint distribution of $(y_{t_1}, \dots, y_{t_r})$ is the same as that of $(y_{t_1 + h}, \dots, y_{t_r + h})$ for every shift $h$ and all times. It is weakly stationary (or covariance stationary) if

1. $\mathbb E[y_t] = \mu$ is the same for all $t$;
2. $\operatorname{Var}(y_t) < \infty$;
3. $\operatorname{Cov}(y_t, y_{t+k}) = \gamma(k)$ depends only on the lag $k$, not on $t$.

Condition 3 with $k = 0$ makes the variance $\gamma(0)$ constant as well. A strictly stationary process with finite variance is weakly stationary. For a Gaussian process the two are the same, because means and covariances fix all its joint distributions. In this track "stationary" means weakly stationary.

Process	Stationary?	Why
white noise	yes	mean 0, variance $\sigma^2$, $\gamma(k) = 0$ for $k \ne 0$
random walk	no	$\operatorname{Var}(y_t) = t\sigma^2$ grows
linear trend + noise	no	the mean $a + bt$ moves
seasonal pattern + noise	no	the mean $S_t$ depends on the season
AR(1) with $\lvert\phi\rvert < 1$	yes	section 3.3

3.3 The AR(1) process #

Definition. $y_t = \phi\, y_{t-1} + \varepsilon_t$ with $\varepsilon_t \sim \mathrm{WN}(0, \sigma^2)$, where $\varepsilon_t$ is uncorrelated with all past values $y_{t-1}, y_{t-2}, \dots$.

Why it forgets. Insert the equation into itself again and again:

A shock from $j$ steps ago has weight $\phi^j$. For $|\phi| < 1$ the weights shrink to 0, and the sum has variance $\sigma^2 (1 + \phi^2 + \phi^4 + \dots) = \sigma^2 / (1 - \phi^2)$. For $\phi = 1$ the weights stay at 1: this is the random walk, and the sum has no finite variance.

Result 2 (moments of a stationary AR(1)). If $|\phi| < 1$ and the process is stationary, then

Proof.

1. Mean. Take expectations: $\mu = \phi\mu + 0$, so $(1 - \phi)\mu = 0$ and $\mu = 0$ (because $\phi \ne 1$).
2. Variance. $\varepsilon_t$ is uncorrelated with $y_{t-1}$, so $\operatorname{Var}(y_t) = \phi^2 \operatorname{Var}(y_{t-1}) + \sigma^2$. Stationarity gives $\gamma(0) = \phi^2 \gamma(0) + \sigma^2$, so $\gamma(0) = \sigma^2 / (1 - \phi^2)$.
3. Autocovariance. Multiply the equation by $y_{t-k}$ for $k \ge 1$ and take expectations: $\mathbb E[y_t y_{t-k}] = \phi\, \mathbb E[y_{t-1} y_{t-k}] + \mathbb E[\varepsilon_t y_{t-k}]$. The last term is 0 (the new shock is uncorrelated with the past). So $\gamma(k) = \phi\, \gamma(k - 1)$.
4. Repeat step 3: $\gamma(k) = \phi^k \gamma(0)$, and dividing by $\gamma(0)$ gives $\rho(k) = \phi^k$. $\blacksquare$

Numbers. For $\phi = 0.8$ and $\sigma = 1$: $\operatorname{Var}(y_t) = 1/0.36 = 2.78$ and $\rho(1), \rho(2), \rho(3) = 0.8, 0.64, 0.512$. A simulation of length $2 \times 10^5$ gives variance 2.77 and sample autocorrelations 0.799, 0.639, 0.511. For $\phi = -0.5$: $\rho(k) = (-0.5)^k = -0.5, 0.25, -0.125, \dots$, which alternates in sign.

The case $\phi = 1$ is called a unit root: the characteristic equation $1 - \phi z = 0$ has its root at $z = 1$. Lesson 05 generalises this to AR($p$): a process is stationary when all roots lie outside the unit circle.

3.4 The autocorrelation function and its sample version #

Definition. For a stationary process, the autocovariance function is $\gamma(k) = \operatorname{Cov}(y_t, y_{t+k})$ and the autocorrelation function (ACF) is $\rho(k) = \gamma(k)/\gamma(0)$. It has $\rho(0) = 1$, $|\rho(k)| \le 1$ and $\rho(-k) = \rho(k)$.

Definition (sample ACF). From data $y_1, \dots, y_n$ with mean $\bar y$:

The sum has only $n - k$ terms, but we divide by $n$. This choice has a reason.

Result 3 (dividing by $n$ gives a valid autocovariance). The $n \times n$ matrix $\hat\Gamma_n$ with entries $\hat\gamma(|i - j|)$ is positive semidefinite. In particular $|\hat\rho(k)| \le 1$.

Proof of Result 3

Let $d_t = y_t - \bar y$ for $1 \le t \le n$ and $d_t = 0$ for all other integers $t$ (zero padding). Then $\sum_{u} d_u d_{u+k} = \sum_{t=1}^{n-|k|} d_t d_{t+|k|} = n\,\hat\gamma(|k|)$. For any vector $a \in \R^n$:

$$a^\top \hat\Gamma_n a = \sum_{i,j} a_i a_j \hat\gamma(|i - j|) = \frac1n \sum_{i,j} a_i a_j \sum_u d_{u} d_{u + i - j} = \frac1n \sum_s \Bigl(\sum_i a_i d_{s - i}\Bigr)^2 \ \ge\ 0 ,$$

where the last step substitutes $s = u + i$ and multiplies out the square. A positive semidefinite matrix has $2 \times 2$ minors $\hat\gamma(0)^2 - \hat\gamma(k)^2 \ge 0$, which gives $|\hat\rho(k)| \le 1$. $\blacksquare$

With the divisor $n - k$ this can fail. For $y = (1, 2, 3)$: $d = (-1, 0, 1)$, the lag-2 sum is $-1$, and dividing by $n - k = 1$ gives $\hat\gamma(2) = -1$ while $\hat\gamma(0) = 2/3$. Then "$\hat\rho(2)$" $= -1.5$, which is impossible for a correlation, and the matrix has the negative eigenvalue $-1/3$. The matrix $\hat\Gamma_n$ is a Toeplitz matrix (constant along each diagonal; see linear-algebra lesson 19). Lesson 05 solves linear systems with it (the Yule–Walker equations), which needs it to be positive definite.

The white-noise band. If the data are i.i.d. with finite variance, then for large $n$ each $\hat\rho(k)$, $k \ge 1$, is approximately $\mathcal N(0, 1/n)$ (Bartlett's formula). So about 95 % of the sample autocorrelations lie inside $\pm 1.96/\sqrt n$. This band is for one lag at a time. With 20 lags of true white noise you expect about one value outside the band by chance; with 40 lags, about two. One bar just outside the band at an odd lag is not evidence of anything.

What the sample ACF looks like:

• white noise: all bars (except lag 0) inside the band;
• AR(1) with $\phi > 0$: geometric decay $\phi^k$; with $\phi < 0$: alternating signs;
• trend or random walk: very slow, almost linear decay from close to 1. The deviations from $\bar y$ come in long runs of one sign, so products at many lags are positive;
• seasonal series: peaks at the lags $m, 2m, 3m, \dots$.

The last two patterns say "not stationary yet: difference first" (section 3.7).

3.5 The Ljung–Box test #

Instead of looking at each bar, test them together.

Definition (Ljung–Box test). The null hypothesis is that the series is white noise up to lag $\ell$: $\rho(1) = \dots = \rho(\ell) = 0$. The statistic is

Under the null hypothesis, $Q$ is approximately $\chi^2$ with $\ell$ degrees of freedom. A large $Q$ (a small p-value) means that some autocorrelation is left.

Why $\chi^2$. Each $\sqrt n\,\hat\rho(k)$ is approximately $\mathcal N(0, 1)$, and for different lags they are approximately independent. A sum of $\ell$ squared independent standard normals is $\chi^2_\ell$. The factor $(n + 2)/(n - k)$ corrects for small samples; without it you get the older Box–Pierce statistic $n \sum_k \hat\rho(k)^2$.

In practice. Use $\ell = 10$ for non-seasonal data and $\ell = 2m$ for seasonal data, but not more than $n/5$ (FPP3 §5.4). When you apply the test to the residuals of a fitted ARMA($p, q$) model, use $\ell - p - q$ degrees of freedom, because fitting the model has already removed some autocorrelation (lesson 06).

3.6 The partial autocorrelation function #

In an AR(1) process, $y_t$ and $y_{t-2}$ are correlated ($\rho(2) = \phi^2$), but only through $y_{t-1}$. The PACF removes such indirect effects.

Definition (PACF). The partial autocorrelation at lag $k$, written $\phi_{kk}$, is the correlation between $y_t$ and $y_{t-k}$ after removing from both the best linear prediction from the values in between, $y_{t-1}, \dots, y_{t-k+1}$. For $k = 1$ there is nothing in between, so $\phi_{11} = \rho(1)$.

Result 4 (PACF as a regression coefficient). Let $\hat y_t = a_1 y_{t-1} + \dots + a_k y_{t-k}$ be the best linear predictor of $y_t$ (smallest mean squared error) from the last $k$ values of a zero-mean stationary process. Then the last coefficient is the PACF: $a_k = \phi_{kk}$ (Shumway and Stoffer, ch. 3; Brockwell and Davis, ch. 2–3).

So the sample PACF can be computed by regression: for $k = 1, 2, \dots$, fit an AR($k$) model by least squares and keep its last coefficient.

Lag 2 by hand. The best coefficients $(a_1, a_2)$ solve the normal equations (the Yule–Walker equations), written with correlations:

By Cramer's rule,

Result 5 (the PACF of an AR process cuts off). For AR(1), $\phi_{11} = \phi$ and $\phi_{kk} = 0$ for all $k \ge 2$. For AR($p$), $\phi_{kk} = 0$ for all $k > p$.

Proof for AR(1), lag 2. $\rho(1) = \phi$, $\rho(2) = \phi^2$, so $\phi_{22} = (\phi^2 - \phi^2)/(1 - \phi^2) = 0$. In general: for $k > p$, the AR($p$) equation $y_t = \phi_1 y_{t-1} + \dots + \phi_p y_{t-p} + \varepsilon_t$ is already the best linear predictor from the last $k$ values, because the error $\varepsilon_t$ is uncorrelated with all past values. Its coefficient on $y_{t-k}$ is 0. $\blacksquare$

Interactive figure: loading…

This gives the classic identification table (lesson 05 explains the MA column):

	ACF	PACF
AR($p$)	decays (geometric or damped waves)	cuts off after lag $p$
MA($q$)	cuts off after lag $q$	decays

The Durbin–Levinson recursion

Fitting $K$ separate regressions costs a lot. The Durbin–Levinson recursion computes all $\phi_{kk}$ from $\rho(1), \dots, \rho(K)$ in $O(K^2)$ operations. With $\phi_{11} = \rho(1)$, for $k \ge 2$:

$$\phi_{kk} = \frac{\rho(k) - \sum_{j=1}^{k-1} \phi_{k-1,j}\, \rho(k - j)}{1 - \sum_{j=1}^{k-1} \phi_{k-1,j}\, \rho(j)}, \qquad \phi_{k,j} = \phi_{k-1,j} - \phi_{kk}\, \phi_{k-1,k-j} \ \ (j < k).$$

For $k = 2$ this is exactly the formula above. It solves the Toeplitz Yule–Walker systems of increasing size, reusing the previous solution (Brockwell and Davis, ch. 2). Using the sample ACF in it gives the "Yule–Walker" sample PACF; it differs slightly from the regression version on short series.

3.7 Making a series stationary #

The backshift operator $B$ shifts back one step: $B y_t = y_{t-1}$. The difference is $\nabla y_t = (1 - B) y_t = y_t - y_{t-1}$ and the seasonal difference is $\nabla_m y_t = (1 - B^m) y_t = y_t - y_{t-m}$.

• A random walk becomes white noise: $\nabla y_t = \varepsilon_t$ (plus $c$ with drift).
• A linear trend becomes a constant: $\nabla(a + bt) = b$.
• A fixed seasonal pattern disappears: $\nabla_m S_t = S_t - S_{t-m} = 0$.
• Log, then difference ≈ growth rate: $\log y_t - \log y_{t-1} = \log(1 + r_t) \approx r_t$, where $r_t = (y_t - y_{t-1})/y_{t-1}$ is the relative change. Take logs first when the variation grows with the level (lesson 01).
• The two differences commute: $(1 - B)(1 - B^m) = (1 - B^m)(1 - B)$. Often one seasonal and one ordinary difference are enough.

Over-differencing. Differencing a series that is already stationary does harm.

Result 6. If $\varepsilon_t$ is white noise, then $x_t = \nabla \varepsilon_t = \varepsilon_t - \varepsilon_{t-1}$ has $\rho_x(1) = -1/2$ and $\rho_x(k) = 0$ for $k \ge 2$.

Proof. $\operatorname{Var}(x_t) = \sigma^2 + \sigma^2 = 2\sigma^2$. $\operatorname{Cov}(x_t, x_{t-1}) = \operatorname{Cov}(\varepsilon_t - \varepsilon_{t-1}, \varepsilon_{t-1} - \varepsilon_{t-2}) = -\operatorname{Var}(\varepsilon_{t-1}) = -\sigma^2$. Lags $k \ge 2$ share no $\varepsilon$, so their covariance is 0. Hence $\rho_x(1) = -\sigma^2/(2\sigma^2) = -1/2$. $\blacksquare$

So a sample ACF with $\hat\rho(1)$ near $-0.5$ after differencing is a sign that you differenced once too often. The same happens when you difference a trend-stationary series $a + bt + \varepsilon_t$: you get $b + \varepsilon_t - \varepsilon_{t-1}$. For such a series it is better to subtract the fitted line.

Interactive figure: loading…

3.8 Unit-root tests: ADF and KPSS #

Looking at plots is not always enough. Two tests are standard, and their null hypotheses are opposite.

• ADF (augmented Dickey–Fuller). Regress $\nabla y_t = \alpha + \beta t + \gamma\, y_{t-1} + \sum_{j=1}^{p} \delta_j \nabla y_{t-j} + e_t$ (the constant and the trend term are optional). Null hypothesis: a unit root ($\gamma = 0$, so differencing is needed). Alternative: $\gamma < 0$, stationary. The test statistic does not follow a t-distribution under the null; it has its own Dickey–Fuller tables. A small p-value means: reject the unit root.
• KPSS (Kwiatkowski–Phillips–Schmidt–Shin). Null hypothesis: stationary (around a constant, or around a line in the trend version). A small p-value means: not stationary, difference. FPP3 uses KPSS repeatedly to choose the number of differences.

ADF (null: unit root)	KPSS (null: stationary)	Reading
small p	large p	stationary
large p	small p	unit root: difference
small p	small p	conflict: maybe trend-stationary, or a structural break
large p	large p	conflict: the data are not informative (short or near a unit root)

A large p-value is never proof of the null hypothesis. It only says that the data do not contradict it. Short series give large p-values for almost everything. Unit-root tests also have low power against $\phi$ close to 1 (for example 0.95).

3.9 Spurious regression #

Spurious regression (Granger and Newbold, 1974): regress one random walk on another, independent random walk. The true slope is 0, so a 5 % t-test should reject in about 5 % of cases. It rejects far more often.

Interactive figure: loading…

Why it happens. The t-test assumes that the regression errors are stationary and uncorrelated. Here the error $y_t - a - b x_t$ is itself a random walk: it is strongly autocorrelated, so the usual standard error is far too small. Phillips (1986) proved that the t-statistic does not converge to a distribution at all; it grows like $\sqrt n$. So more data make the false result more significant. In simulations with $n = 100$, about 77 % of the experiments reject; with $n = 200$, about 83 %. Warning signs are a high $R^2$ together with strongly autocorrelated residuals.

Fixes. Difference both series (regress $\nabla y$ on $\nabla x$), add lagged values of both series to the regression, or test for cointegration.

Cointegration

Two non-stationary series $x_t$ and $y_t$ (each with one unit root) are cointegrated if some combination $y_t - \beta x_t$ is stationary. Then they share a common stochastic trend, and the levels regression is meaningful: its slope estimate converges even faster than usual. Prices of the same good in two markets, or short and long interest rates, are classic examples. The Engle–Granger test (Engle and Granger, 1987) regresses $y$ on $x$ and runs a unit-root test on the residuals, with special critical values. Differencing cointegrated series throws away the long-run relation; error correction models keep it.

The periodogram: the ACF in frequency space

For a stationary process with $\sum_k |\gamma(k)| < \infty$, the spectral density is $f(\omega) = \sum_{k=-\infty}^{\infty} \gamma(k)\, e^{-2\pi i \omega k}$ for $\omega \in [-\tfrac12, \tfrac12]$ (cycles per time step). It says how the variance is spread over frequencies. The ACF and the spectral density are a Fourier pair (the Wiener–Khinchin theorem). White noise has a flat spectrum $f(\omega) = \sigma^2$ ("white", like white light); an AR(1) process has $f(\omega) = \sigma^2 / |1 - \phi e^{-2\pi i\omega}|^2$, which for $\phi > 0$ is largest at low frequencies ("red noise", the standard null model for climate series).

The sample version is the periodogram $I(\omega_j) = \frac1n \bigl|\sum_{t=1}^n (y_t - \bar y)\, e^{-2\pi i \omega_j t}\bigr|^2$ at the Fourier frequencies $\omega_j = j/n$. For $j \ne 0$ it equals exactly $\sum_{|k| < n} \hat\gamma(k)\, e^{-2\pi i\omega_j k}$, with $\hat\gamma$ divided by $n$. A seasonal series with period $m$ has a peak at $\omega = 1/m$ and at its harmonics $2/m, 3/m, \dots$. The periodogram is not a consistent estimator of $f$ (its variance does not shrink with $n$); smoothing it fixes this (Shumway and Stoffer, ch. 4).

Answers to section 1. (1) $\operatorname{Var}(y_{100}) = 100\sigma^2 = 100$, so the standard deviation is 10 (Result 1). (2) Negative: $\rho(3) = (-0.5)^3 = -0.125$ (Result 2). (3) Far more than 5 %: about 77 % of the time for $n = 100$ (section 3.9).

4. Worked example #

A short series with $n = 8$ (too short for real use, but good for computing by hand):

Step 1: mean and deviations. $\bar y = 28/8 = 3.5$. The deviations $d_t = y_t - \bar y$ are

Step 2: autocovariances (divide by $n = 8$).

• Lag 0: $\sum d_t^2 = 0.25 + 2.25 + 6.25 + 0.25 + 2.25 + 6.25 + 0.25 + 0.25 = 18$, so $\hat\gamma(0) = 18/8 = 2.25$.
• Lag 1: the 7 products $d_t d_{t+1}$ are $-0.75,\ 3.75,\ 1.25,\ -0.75,\ 3.75,\ 1.25,\ -0.25$. Their sum is $8.25$, so $\hat\gamma(1) = 8.25/8 = 1.03125$.
• Lag 2: the 6 products $d_t d_{t+2}$ are $-1.25,\ 0.75,\ -3.75,\ -1.25,\ 0.75,\ -1.25$. Their sum is $-6$, so $\hat\gamma(2) = -6/8 = -0.75$.

Step 3: autocorrelations. $\hat\rho(1) = 1.03125/2.25 = 0.458$ and $\hat\rho(2) = -0.75/2.25 = -0.333$.

Step 4: the band. $1.96/\sqrt 8 = 0.693$. Both values are inside the band. With only 8 points, even $\hat\rho(1) = 0.46$ is not significant.

Step 5: PACF at lag 2.

The partial autocorrelation is larger in size than $\hat\rho(2)$. Once $y_{t-1}$ is known, a high $y_{t-2}$ points to a lower $y_t$.

Step 6: Ljung–Box with $\ell = 2$.

For $\chi^2$ with 2 degrees of freedom, $P(\chi^2_2 > q) = e^{-q/2}$, so the p-value is $e^{-1.94} = 0.144$. We cannot reject "white noise up to lag 2", which agrees with step 4.

Step 7: an AR(1) check. If this series came from an AR(1) process, we would expect $\rho(2) = \rho(1)^2 > 0$. The estimate $\hat\rho(2) = -0.33$ has the wrong sign, and $\hat\phi_{22}$ is far from 0. So AR(1) is not a good description (the series looks more like a short wave), although 8 points cannot prove it.

5. In code #

The first part checks the worked example and Results 2 and 6. The second part runs ADF and KPSS with statsmodels.

import numpy as np
from scipy import stats

def sample_acf(y, max_lag):
    d = np.asarray(y, float) - np.mean(y)
    n = len(d)
    gamma = np.array([d[: n - k] @ d[k:] / n for k in range(max_lag + 1)])   # divide by n
    return gamma / gamma[0]

y = np.array([3, 5, 6, 4, 2, 1, 3, 4.])
r = sample_acf(y, 2)
print(r.round(4))                                   # [ 1.      0.4583 -0.3333]
print(round((r[2] - r[1]**2) / (1 - r[1]**2), 4))   # -0.6879  (PACF at lag 2)

n, ell = len(y), 2
Q = n * (n + 2) * np.sum(r[1:]**2 / (n - np.arange(1, ell + 1)))
print(round(Q, 3), round(stats.chi2.sf(Q, df=ell), 4))   # 3.882 0.1435

# AR(1) with phi = 0.8: variance and ACF against theory
rng = np.random.default_rng(0)
phi, n = 0.8, 200_000
e = rng.normal(size=n)
x = np.empty(n); x[0] = e[0] / np.sqrt(1 - phi**2)      # start in the stationary distribution
for t in range(1, n):
    x[t] = phi * x[t - 1] + e[t]
print(round(x.var(), 2), 1 / (1 - phi**2))          # 2.77 2.7777777777777786
print(sample_acf(x, 3)[1:].round(3), phi ** np.arange(1, 4))   # [0.799 0.639 0.511] [0.8   0.64  0.512]

# over-differencing white noise gives rho(1) = -0.5
w = rng.normal(size=100_000)
print(sample_acf(np.diff(w), 1)[1].round(3))          # -0.5

import warnings
import numpy as np
from statsmodels.tsa.stattools import adfuller, kpss     # statsmodels 0.15

rng = np.random.default_rng(3)
walk = np.cumsum(rng.normal(size=300))           # unit root
ar = np.zeros(300)
for t in range(1, 300):
    ar[t] = 0.5 * ar[t - 1] + rng.normal()        # stationary AR(1)
with warnings.catch_warnings():
    warnings.simplefilter("ignore")               # KPSS warns when p is outside its table
    for name, s in [("walk", walk), ("diff(walk)", np.diff(walk)), ("AR(1)", ar)]:
        print(f"{name:11s} ADF p = {adfuller(s)[1]:.3f}   KPSS p = {kpss(s, nlags='auto')[1]:.3f}")
# walk        ADF p = 0.858   KPSS p = 0.010
# diff(walk)  ADF p = 0.000   KPSS p = 0.100
# AR(1)       ADF p = 0.000   KPSS p = 0.100

statsmodels reports KPSS p-values only between 0.01 and 0.1 (from a table): 0.010 means "0.01 or smaller" and 0.100 means "0.1 or larger". The random walk: ADF cannot reject a unit root and KPSS rejects stationarity, so difference it. After one difference both tests agree that it is stationary.

6. Where it shows up #

• Model identification and residual checks. The ACF/PACF table of section 3.6 is the first step of the Box–Jenkins method (lessons 05 and 06). Every forecasting model in this track ends with an ACF plot and a Ljung–Box test of its residuals.
• Feature selection for ML. Lags with large ACF or PACF values, and the seasonal lags $m, 2m$, are the standard lag features of gradient-boosting forecasters (lesson 08). Differencing or scaling the input window is the same idea as making the series stationary; deep forecasters such as RevIN-based models normalise each window for this reason (lesson 09).
• Finance. Daily returns (log-differences of prices) have almost no autocorrelation, but squared returns are strongly autocorrelated: volatility clusters. This is what GARCH models describe (lesson 06).
• MCMC and simulation. The samples of a Markov chain are autocorrelated. The effective sample size $n / (1 + 2\sum_{k \ge 1} \rho(k))$ uses the ACF of the chain; for an AR(1)-like chain it is about $n (1 - \phi)/(1 + \phi)$. With $\phi = 0.9$, 10 000 draws are worth only about 526 independent ones.
• ML for science. Climate scientists test whether a spectral peak is real against an AR(1) "red noise" background, and use the effective sample size when they compute the significance of trends in autocorrelated data.
• Drift and monitoring. Many drift tests assume independent samples; with autocorrelated metrics they raise false alarms (see the MLOps track).

7. Common mistakes #

• Reading one bar outside the band as a discovery. The band is for one lag at a time. Look for a pattern (decay, cut-off, seasonal peaks) or use the Ljung–Box test.
• Dividing the sample autocovariance by $n - k$. It can give "correlations" larger than 1 and a matrix that is not positive semidefinite. Divide by $n$.
• "Stationary means no visible pattern." An AR(1) path with $\phi = 0.95$ looks like it has trends, but it is stationary. A seasonal series looks regular, but it is not stationary.
• Treating a large ADF p-value as proof of a unit root (or a large KPSS p-value as proof of stationarity). Absence of evidence is not evidence. Use both tests and the plots.
• Differencing too often. $\hat\rho(1) \approx -0.5$ after differencing means the last difference was not needed.
• Regressing one trending series on another. Two unrelated random walks give a "significant" slope most of the time. Difference both, or test for cointegration.
• Mixing up ACF and PACF. The ACF of an AR(1) process decays; its PACF cuts off after lag 1.

8. Check yourself #

1. Random walk with drift: $y_t = 0.5 + y_{t-1} + \varepsilon_t$, $y_0 = 0$, $\sigma = 2$. Give the mean and the variance of $y_{25}$.
   Answer

   $\mathbb E[y_{25}] = 25 \cdot 0.5 = 12.5$ and $\operatorname{Var}(y_{25}) = 25 \cdot 2^2 = 100$ (standard deviation 10).
2. A stationary AR(1) process has $\operatorname{Var}(y_t) = 4$ and $\sigma^2 = 3$. What can $\phi$ be, and what is $\rho(2)$?
   Answer

   $4 = 3/(1 - \phi^2)$ gives $\phi^2 = 1/4$, so $\phi = 0.5$ or $\phi = -0.5$. In both cases $\rho(2) = \phi^2 = 0.25$.
3. You compute the sample ACF of 400 residuals at lags 1 to 20. Where is the band, and how many bars outside it are normal for white noise?
   Answer

   $\pm 1.96/\sqrt{400} = \pm 0.098$. About $20 \times 0.05 = 1$ bar outside the band is normal by chance.
4. For a stationary process $\rho(1) = 0.5$ and $\rho(2) = 0.4$. Compute $\phi_{22}$. Could the process be AR(1)?
   Answer

   $\phi_{22} = (0.4 - 0.25)/(1 - 0.25) = 0.15/0.75 = 0.2$. For AR(1) we would need $\rho(2) = \rho(1)^2 = 0.25$ and $\phi_{22} = 0$, so this is not AR(1).
5. After one difference, the sample ACF of a long series has $\hat\rho(1) = -0.48$ and nothing else outside the band. What happened, and what do you do?
   Answer

   This is the signature of over-differencing (Result 6: $\rho(1) = -0.5$). The original series was probably already stationary (or trend-stationary). Model the original series, or subtract a fitted trend instead of differencing.
6. ADF gives p = 0.01 and KPSS gives p = 0.01 on the same series. What could explain it?
   Answer

   Both nulls are rejected, so the tests conflict. Common causes: a trend-stationary series tested with a KPSS version without a trend, a structural break (a sudden level shift), or changing variance. Look at the plot, try the trend versions of the tests, and consider a model for the break.

Recall cardsAnswer in your head first, then reveal. Add the cards to your review deck to see them again just before you would forget.

1. Define weak (covariance) stationarity.

   Show answer

2. What is the variance of a random walk $y_t = y_{t-1} + \varepsilon_t$ with $y_0 = 0$ at time $t$?

   Show answer

3. What is the variance of a stationary AR(1) process $y_t = \phi y_{t-1} + \varepsilon_t$?

   Show answer

4. What is the ACF of a stationary AR(1) process? Give $\rho(3)$ for $\phi = -0.5$.

   Show answer

5. Why does the sample autocovariance divide by $n$ and not by $n - k$?

   Show answer

6. Where is the approximate 95 % band of a sample ACF for white noise, and how must it be read?

   Show answer

7. What does the Ljung–Box test test, and what is its statistic?

   Show answer

8. What is the PACF at lag $k$?

   Show answer

9. What are the null hypotheses of the ADF and KPSS tests?

   Show answer

10. What is the sign of over-differencing in the ACF?

    Show answer

Add these cards to my reviewOpen the review →

9. Interview angle #

• "What is stationarity, and why do we care?" — Constant mean, variance and lag-dependent covariance. Models estimated on the past can only describe the future if these properties stay the same; ARMA theory and most tests assume it.
• "How do you make a series stationary?" — Log or Box–Cox for growing variance, seasonal difference for seasons, ordinary difference for trends or unit roots; check with plots, ACF and KPSS/ADF; do not over-difference.
• "ACF vs PACF: how do you use them?" — The ACF of an AR($p$) decays and its PACF cuts off after $p$; the reverse for MA($q$).
• "Why can a regression between two time series be misleading?" — Spurious regression: with non-stationary series the t-test is invalid; the t-statistic grows with $n$. Difference, or test for cointegration.
• "Derive the ACF of an AR(1)." — Multiply by $y_{t-k}$, take expectations, use that the shock is uncorrelated with the past: $\gamma(k) = \phi\gamma(k-1)$.

10. Practice #

The notebook for this lesson runs in Colab. Each exercise has a check cell and a folded solution:

1. Write sample_acf(y, max_lag) (divide by $n$) and check it on the worked example.
2. Write simulate_ar1(phi, n, seed) that starts in the stationary distribution, and compare its sample variance and ACF with $\sigma^2/(1-\phi^2)$ and $\phi^k$ for several $\phi$.
3. Write pacf_by_regression(y, max_lag) with np.linalg.lstsq, and check it on AR(1) and AR(2) data.
4. Write ljung_box(x, lag) returning $Q$ and the p-value, and check that p-values of white noise are roughly uniform while AR(1) data give small p-values.
5. Spurious regression: write slope_t_stat(x, y) and estimate the rejection rate for 1000 pairs of independent random walks, in levels and in differences.

11. Watch and read #

Time Series Talk : Stationarityritvikmath · 10 min

The three conditions of weak stationarity, with pictures of series that break each one.

Time Series Talk : Autocorrelation and Partial Autocorrelationritvikmath · 13 min

The difference between the ACF and the PACF, explained with regressions.

Unit Roots : Time Series Talkritvikmath · 14 min

Why $\phi = 1$ is special, and what the Dickey–Fuller test checks.

Spurious regressionBen Lambert · 3 min

Two unrelated trending series that look strongly related.

8. Time Series Analysis IMIT OpenCourseWare · 18.S096 Topics in Mathematics with Applications in Finance · 76 min

A university lecture on stationarity, the Wold representation and ARMA models; watch it after lesson 05 too.

• Hyndman and Athanasopoulos, Forecasting: Principles and Practice, 3rd ed. (free online): §2.8 autocorrelation and §2.9 white noise, §5.4 residual diagnostics (Ljung–Box), §9.1 stationarity and differencing (KPSS).
• Shumway and Stoffer, Time Series Analysis and Its Applications, 5th ed. (Springer, 2025), chapters 1–3 (stationarity, ACF, PACF, AR models) and 4 (spectral analysis); the R code of the book is free on GitHub.
• Penn State STAT 510, Applied Time Series Analysis (free course notes): lessons 1 and 2 on the ACF, the PACF and AR models.
• Brockwell and Davis, Introduction to Time Series and Forecasting (Springer): Bartlett's formula for the sample ACF, the PACF and the Durbin–Levinson recursion.
• Granger and Newbold (1974), Spurious regressions in econometrics, Journal of Econometrics 2(2) — the original simulation study.