Прогнозные интервалы для регрессии градиентного бустинга

Построение нелинейных регрессоров квантилей и наименьших квадратов

Обучите модели градиентного бустинга, обученные с потерей квантиля и alpha=0.05, 0.5, 0.95.

Полученные модели для alpha=0.05 и alpha=0.95 создают 90%-ный доверительный интервал (95% - 5% = 90%).

Модель, обученная с alpha=0.5, создаёт регрессию медианы: в среднем должно быть одинаковое количество целевых наблюдений выше и ниже предсказанных значений.

from sklearn.ensemble import GradientBoostingRegressor
from sklearn.metrics import mean_pinball_loss, mean_squared_error

all_models = {}
common_params = dict(
    learning_rate=0.05,
    n_estimators=200,
    max_depth=2,
    min_samples_leaf=9,
    min_samples_split=9,
)
for alpha in [0.05, 0.5, 0.95]:
    gbr = GradientBoostingRegressor(loss="quantile", alpha=alpha, **common_params)
    all_models["q %1.2f" % alpha] = gbr.fit(X_train, y_train)

Обратите внимание, что HistGradientBoostingRegressor намного быстрее, чем GradientBoostingRegressor начиная со средних наборов данных (n_samples >= 10_000), что не наблюдается в данном примере.

Для сравнения, мы также обучили базовая модель, обученная с помощью обычной (среднеквадратической) ошибки (MSE).

gbr_ls = GradientBoostingRegressor(loss="squared_error", **common_params)
all_models["mse"] = gbr_ls.fit(X_train, y_train)

Создайте равномерно распределённый набор оценочных значений входных данных, охватывающий диапазон [0, 10].

xx = np.atleast_2d(np.linspace(0, 10, 1000)).T

Постройте истинную условную среднюю функцию f, прогнозы условной средней (потеря равна квадратичной ошибке), условную медиану и условный 90%-ный интервал (от 5-го до 95-го условных перцентилей).

import matplotlib.pyplot as plt

y_pred = all_models["mse"].predict(xx)
y_lower = all_models["q 0.05"].predict(xx)
y_upper = all_models["q 0.95"].predict(xx)
y_med = all_models["q 0.50"].predict(xx)

fig = plt.figure(figsize=(10, 10))
plt.plot(xx, f(xx), "g:", linewidth=3, label=r"$f(x) = x\,\sin(x)$")
plt.plot(X_test, y_test, "b.", markersize=10, label="Test observations")
plt.plot(xx, y_med, "r-", label="Predicted median")
plt.plot(xx, y_pred, "r-", label="Predicted mean")
plt.plot(xx, y_upper, "k-")
plt.plot(xx, y_lower, "k-")
plt.fill_between(
    xx.ravel(), y_lower, y_upper, alpha=0.4, label="Predicted 90% interval"
)
plt.xlabel("$x$")
plt.ylabel("$f(x)$")
plt.ylim(-10, 25)
plt.legend(loc="upper left")
plt.show()

Сравнивая предсказанную медиану с предсказанной средней, мы отмечаем, что медиана в среднем ниже средней, поскольку шум смещён в сторону высоких значений (больших выбросов). Оценка медианы также кажется более гладкой из-за её естественной устойчивости к выбросам.

Также обратите внимание, что индуцированный предвзятость деревьев градиентного бустинга, к сожалению, мешает нашему 0.05 квантилю полностью захватить синусоидальную форму сигнала, особенно около x=8. Настройка гиперпараметров может уменьшить этот эффект, как показано в последней части этого блокнота.

Анализ метрик ошибок

Измерьте модели с помощью mean_squared_error и mean_pinball_loss метрик на обучающей выборке.

import pandas as pd


def highlight_min(x):
    x_min = x.min()
    return ["font-weight: bold" if v == x_min else "" for v in x]


results = []
for name, gbr in sorted(all_models.items()):
    metrics = {"model": name}
    y_pred = gbr.predict(X_train)
    for alpha in [0.05, 0.5, 0.95]:
        metrics["pbl=%1.2f" % alpha] = mean_pinball_loss(y_train, y_pred, alpha=alpha)
    metrics["MSE"] = mean_squared_error(y_train, y_pred)
    results.append(metrics)

pd.DataFrame(results).set_index("model").style.apply(highlight_min)

Один столбец показывает все модели, оцененные по одной метрике. Минимальное значение в столбце должно быть получено, когда модель обучена и измерена с помощью той же метрики. Это всегда должно происходить на обучающей выборке, если обучение сошлось.

Обратите внимание, что из-за асимметрии распределения целевой переменной ожидаемые условные среднее и условное медианное значения существенно отличаются, и поэтому модель с квадратичной ошибкой не может дать хорошую оценку условной медианы и наоборот.

Если распределение целевой переменной было бы симметричным и не имело выбросов (например, с нормальным шумом), то оценка медианы и оценка методом наименьших квадратов дали бы похожие предсказания.

Затем мы делаем то же самое на тестовой выборке.

results = []
for name, gbr in sorted(all_models.items()):
    metrics = {"model": name}
    y_pred = gbr.predict(X_test)
    for alpha in [0.05, 0.5, 0.95]:
        metrics["pbl=%1.2f" % alpha] = mean_pinball_loss(y_test, y_pred, alpha=alpha)
    metrics["MSE"] = mean_squared_error(y_test, y_pred)
    results.append(metrics)

pd.DataFrame(results).set_index("model").style.apply(highlight_min)

Ошибки выше, что означает, что модели немного переобучились. Это всё ещё показывает, что лучшая тестовая метрика получается, когда модель обучается с минимизацией этой же метрики.

Обратите внимание, что оценка условной медианы конкурирует с оценкой методом наименьших квадратов по метрике MSE на тестовой выборке: это можно объяснить тем, что оценка методом наименьших квадратов очень чувствительна к большим выбросам, что может привести к значительному переобучению. Это видно на правой стороне предыдущего графика. Оценка условной медианы смещена (недооценка для этого асимметричного шума), но также естественно устойчива к выбросам и переобучается меньше.

Калибровка доверительного интервала

Мы также можем оценить способность двух крайних оценок квантилей к созданию хорошо откалиброванного условного 90%-ного доверительного интервала.

Для этого мы можем вычислить долю наблюдений, которые попадают между прогнозами:

def coverage_fraction(y, y_low, y_high):
    return np.mean(np.logical_and(y >= y_low, y <= y_high))


coverage_fraction(
    y_train,
    all_models["q 0.05"].predict(X_train),
    all_models["q 0.95"].predict(X_train),
)

np.float64(0.9)

На обучающей выборке калибровка очень близка к ожидаемому значению охвата для 90%-ного доверительного интервала.

coverage_fraction(
    y_test, all_models["q 0.05"].predict(X_test), all_models["q 0.95"].predict(X_test)
)

np.float64(0.868)

На тестовой выборке оценочный доверительный интервал немного слишком узкий. Тем не менее, необходимо будет заключить эти метрики в цикл перекрестной проверки, чтобы оценить их изменчивость при повторном применении данных.