skforecast/skforecast-ai/autocorrelation-and-lag-selection

Autocorrelation and Lag Selection

When to Use

Use this skill when the user wants to:

• Inspect the ACF / PACF of a series before training a forecaster.
• Pick a sensible starting set of lags for ForecasterRecursive, ForecasterDirect, or any multi-series forecaster.
• Detect seasonality (peaks at multiples of the period) from autocorrelation.
• Replace statsmodels.tsa.stattools.acf / pacf calls with a faster equivalent on long series.

Related skills

• Next: feature-engineering (turn the candidate lags into a full feature pipeline with rolling and calendar features)
• Next: hyperparameter-optimization (refine the candidate set via cross-validation with bayesian_search_forecaster)
• Next: feature-selection (prune redundant lags or exog variables with select_features once the forecaster is configured)

Overview

The skforecast functions mirror the public API of statsmodels.tsa.stattools.acf / pacf but are reimplemented from scratch (FFT for the ACF, Levinson–Durbin for the PACF) to maximise speed on long series. Only the most common configuration options are exposed.

If the user needs the wider catalogue offered by statsmodels — alternative estimators ('ywunbiased', 'ols', 'ld', 'burg'), Ljung–Box statistics, etc. — recommend calling statsmodels.tsa directly.

For plotting, keep the statsmodels helpers with their fastest configuration: plot_acf(..., alpha=0.05, fft=True) and plot_pacf(..., alpha=0.05, method='burg').

Computing the ACF

from skforecast.stats import acf
# Just the values (lag 0 included; output length = nlags + 1)
acf_vals = acf(y, nlags=20)
# With 95% Bartlett confidence intervals
acf_vals, confint = acf(y, nlags=20, alpha=0.05)
# Unbiased estimator (denominator n - k)
acf_vals = acf(y, nlags=20, adjusted=True)

Notes:

• nlags must satisfy 0 < nlags < len(x). If None, defaults to min(int(10 * log10(n)), n - 1).
• acf_vals[0] is always 1.0.
• Default is the biased estimator (adjusted=False, denominator n). It guarantees a positive semi-definite Toeplitz matrix and is required for downstream PACF computation.

Computing the PACF

from skforecast.stats import pacf
pacf_vals = pacf(y, nlags=20)
pacf_vals, confint = pacf(y, nlags=20, alpha=0.05)

Notes:

• nlags must satisfy 0 < nlags < len(x) // 2 — Levinson–Durbin becomes unreliable when the AR order approaches half the sample size; the function raises explicitly.
• Confidence intervals are asymptotic under the white-noise null: ±z_{α/2} / sqrt(n) for all lags ≥ 1. Lag 0 has no uncertainty.
• The implementation uses the biased ACF internally, so values can differ slightly (~2e-2 for typical series) from statsmodels.pacf(method='yw'), which uses the unbiased estimator.

Reading ACF / PACF Patterns

import matplotlib.pyplot as plt
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
fig, axes = plt.subplots(1, 2, figsize=(11, 3))
plot_acf(y, lags=36, alpha=0.05, fft=True, ax=axes[0])
plot_pacf(y, lags=36, alpha=0.05, method='burg', ax=axes[1])
plt.tight_layout()

Ranking Lags with calculate_lag_autocorrelation

from skforecast.stats import calculate_lag_autocorrelation
# Default: rows ordered by |PACF| descending, lags 1..n_lags
lag_table = calculate_lag_autocorrelation(
    data=y,
    n_lags=24,
    sort_by='partial_autocorrelation_abs',
)
# Other valid sort keys
calculate_lag_autocorrelation(data=y, n_lags=24, sort_by='lag')
calculate_lag_autocorrelation(data=y, n_lags=24, sort_by='autocorrelation_abs')

Returns a DataFrame with columns ['lag', 'partial_autocorrelation_abs', 'partial_autocorrelation', 'autocorrelation_abs', 'autocorrelation'] covering lags 1..n_lags (lag 0 dropped because it is always 1.0).

last_n_samples keeps only the most recent observations — useful when the series is very long and only the recent dynamics matter:

calculate_lag_autocorrelation(data=y, n_lags=24, last_n_samples=2000)

data accepts a pd.Series or a single-column pd.DataFrame. n_lags must be strictly less than len(data) // 2.

End-to-End: From PACF to a Forecaster

import numpy as np
from skforecast.stats import calculate_lag_autocorrelation
from skforecast.recursive import ForecasterRecursive
from sklearn.ensemble import HistGradientBoostingRegressor
# 1. Rank lags by |PACF|
lag_table = calculate_lag_autocorrelation(data=y_train, n_lags=24)
# 2. Threshold against the asymptotic 95% white-noise band (1.96 / sqrt(n))
threshold = 1.96 / np.sqrt(len(y_train))
selected_lags = (
    lag_table
    .loc[lag_table['partial_autocorrelation_abs'] > threshold, 'lag']
    .astype(int)
    .sort_values()
    .tolist()
)
# 3. Feed directly into any forecaster
forecaster = ForecasterRecursive(
    estimator=HistGradientBoostingRegressor(random_state=963),
    lags=selected_lags,
)
forecaster.fit(y=y_train)
predictions = forecaster.predict(steps=12)

The selected lags are an informed starting point, not a final answer — refine with bayesian_search_forecaster (see hyperparameter-optimization) or prune with select_features (see feature-selection).

Common Mistakes

1. Using `nlags >= n // 2` for PACF. Levinson–Durbin loses stability as the AR order approaches half the sample size; the function raises explicitly. Reduce nlags or use a longer series.
2. Reading the ACF to choose AR order. ACF is contaminated by indirect dependencies through earlier lags and rarely cuts off cleanly for AR processes — use the PACF for AR order, the ACF for MA order.
3. Treating the confidence bands as joint tests. Bartlett (ACF) and white-noise (PACF) bands are pointwise. Several lags will cross the band by chance; use them to spot candidates, not to formally reject white noise.
4. Skipping detrending / differencing. A strong trend dominates the ACF and masks higher-lag structure. If the series shows visible trend, read ACF / PACF on y.diff().dropna() (or use differentiation=1 in the forecaster).
5. Switching to `adjusted=True` by default. The biased estimator is preferred for forecasting workflows; the unbiased one can produce |ACF| > 1 at high lags and breaks the positive semi-definite property required by Levinson–Durbin.
6. Expecting statsmodels-only options. method='ols', method='burg', alternative confidence intervals, etc. are not exposed by skforecast.stats.pacf. For those, call statsmodels.tsa.stattools.pacf directly — the API is the same.
7. Passing a multi-column DataFrame to `calculate_lag_autocorrelation`. It accepts a Series or a single-column DataFrame and raises otherwise. Select the target column first.