Графики частичной зависимости и индивидуального условного математического ожидания

Предварительная обработка для модели нейронной сети

Для масштабирования числовых признаков мы используем QuantileTransformer, а категориальные признаки кодируем с помощью OneHotEncoder.

from sklearn.compose import ColumnTransformer
from sklearn.preprocessing import OneHotEncoder, QuantileTransformer

mlp_preprocessor = ColumnTransformer(
    transformers=[
        ("num", QuantileTransformer(n_quantiles=100), numerical_features),
        ("cat", OneHotEncoder(handle_unknown="ignore"), categorical_features),
    ]
)
mlp_preprocessor

ColumnTransformer(transformers=[('num', QuantileTransformer(n_quantiles=100),
                                 ['temp', 'feel_temp', 'humidity',
                                  'windspeed']),
                                ('cat', OneHotEncoder(handle_unknown='ignore'),
                                 Index(['season', 'month', 'hour', 'holiday', 'weekday', 'workingday',
       'weather'],
      dtype='str'))])

['temp', 'feel_temp', 'humidity', 'windspeed']

Index(['season', 'month', 'hour', 'holiday', 'weekday', 'workingday',
       'weather'],
      dtype='str')

Предобработчик для модели градиентного бустинга

Для модели градиентного бустинга мы оставляем числовые признаки без изменений и кодируем только категориальные признаки с помощью OrdinalEncoder.

from sklearn.preprocessing import OrdinalEncoder

hgbdt_preprocessor = ColumnTransformer(
    transformers=[
        ("cat", OrdinalEncoder(), categorical_features),
        ("num", "passthrough", numerical_features),
    ],
    sparse_threshold=1,
    verbose_feature_names_out=False,
).set_output(transform="pandas")
hgbdt_preprocessor

ColumnTransformer(sparse_threshold=1,
                  transformers=[('cat', OrdinalEncoder(),
                                 Index(['season', 'month', 'hour', 'holiday', 'weekday', 'workingday',
       'weather'],
      dtype='str')),
                                ('num', 'passthrough',
                                 ['temp', 'feel_temp', 'humidity',
                                  'windspeed'])],
                  verbose_feature_names_out=False)

Index(['season', 'month', 'hour', 'holiday', 'weekday', 'workingday',
       'weather'],
      dtype='str')

['temp', 'feel_temp', 'humidity', 'windspeed']

Многослойный перцептрон

Обучим MLPRegressor и построим графики частичной зависимости для отдельных переменных.

from time import time

from sklearn.neural_network import MLPRegressor
from sklearn.pipeline import make_pipeline

print("Training MLPRegressor...")
tic = time()
mlp_model = make_pipeline(
    mlp_preprocessor,
    MLPRegressor(
        hidden_layer_sizes=(30, 15),
        learning_rate_init=0.01,
        early_stopping=True,
        random_state=0,
    ),
)
mlp_model.fit(X_train, y_train)
print(f"done in {time() - tic:.3f}s")
print(f"Test R2 score: {mlp_model.score(X_test, y_test):.2f}")

Training MLPRegressor...
done in 0.745s
Test R2 score: 0.61

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

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

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

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

Построим график усредненной частичной зависимости.

import matplotlib.pyplot as plt

from sklearn.inspection import PartialDependenceDisplay

common_params = {
    "subsample": 50,
    "n_jobs": 2,
    "grid_resolution": 20,
    "random_state": 0,
}

print("Computing partial dependence plots...")
features_info = {
    # features of interest
    "features": ["temp", "humidity", "windspeed", "season", "weather", "hour"],
    # type of partial dependence plot
    "kind": "average",
    # information regarding categorical features
    "categorical_features": categorical_features,
}
tic = time()
_, ax = plt.subplots(ncols=3, nrows=2, figsize=(9, 8), constrained_layout=True)
display = PartialDependenceDisplay.from_estimator(
    mlp_model,
    X_train,
    **features_info,
    ax=ax,
    **common_params,
)
print(f"done in {time() - tic:.3f}s")
_ = display.figure_.suptitle(
    (
        "Partial dependence of the number of bike rentals\n"
        "for the bike rental dataset with an MLPRegressor"
    ),
    fontsize=16,
)

Computing partial dependence plots...
done in 0.956s

Градиентный бустинг

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

from sklearn.ensemble import HistGradientBoostingRegressor

print("Training HistGradientBoostingRegressor...")
tic = time()
hgbdt_model = make_pipeline(
    hgbdt_preprocessor,
    HistGradientBoostingRegressor(
        categorical_features=categorical_features,
        random_state=0,
        max_iter=50,
    ),
)
hgbdt_model.fit(X_train, y_train)
print(f"done in {time() - tic:.3f}s")
print(f"Test R2 score: {hgbdt_model.score(X_test, y_test):.2f}")

Training HistGradientBoostingRegressor...
done in 0.118s
Test R2 score: 0.62

Здесь мы использовали гиперпараметры модели градиентного бустинга по умолчанию и не выполняли предобработку, поскольку модели на основе деревьев по своей природе устойчивы к монотонным преобразованиям числовых признаков.

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

Построим графики частичной зависимости для некоторых числовых и категориальных признаков.

print("Computing partial dependence plots...")
tic = time()
_, ax = plt.subplots(ncols=3, nrows=2, figsize=(9, 8), constrained_layout=True)
display = PartialDependenceDisplay.from_estimator(
    hgbdt_model,
    X_train,
    **features_info,
    ax=ax,
    **common_params,
)
print(f"done in {time() - tic:.3f}s")
_ = display.figure_.suptitle(
    (
        "Partial dependence of the number of bike rentals\n"
        "for the bike rental dataset with a gradient boosting"
    ),
    fontsize=16,
)

Computing partial dependence plots...
done in 1.421s

Анализ графиков

Сначала рассмотрим графики PDP для числовых признаков. Для обеих моделей общий тренд PDP для температуры показывает, что число аренд велосипедов растет с повышением температуры. Аналогичный анализ для влажности показывает противоположный тренд: число аренд велосипедов уменьшается с повышением влажности. Наконец, мы видим такую же тенденцию для скорости ветра: в обеих моделях число аренд велосипедов снижается при увеличении скорости ветра. Мы также видим, что предсказания MLPRegressor значительно более плавные, чем предсказания HistGradientBoostingRegressor.

Теперь рассмотрим графики частичной зависимости для категориальных признаков.

Для признака «сезон» мы видим, что самый низкий столбец соответствует весне. Для признака «погода» самый низкий столбец соответствует дождю. Для признака «час» видны два пика — около 7 утра и 6 вечера. Эти результаты согласуются с выводами, сделанными нами ранее при анализе набора данных.

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

ICE и PDP

PDP представляет собой усреднение маргинальных эффектов признаков. Мы усредняем отклики всех наблюдений из предоставленного набора. В результате некоторые эффекты могут оставаться незаметными. Для этого можно построить график отклика каждого отдельного наблюдения. Такое представление называется графиком индивидуальных эффектов (ICE). На графике ниже показаны ICE для 50 случайно выбранных наблюдений по признакам температуры и влажности.

print("Computing partial dependence plots and individual conditional expectation...")
tic = time()
_, ax = plt.subplots(ncols=2, figsize=(6, 4), sharey=True, constrained_layout=True)

features_info = {
    "features": ["temp", "humidity"],
    "kind": "both",
    "centered": True,
}

display = PartialDependenceDisplay.from_estimator(
    hgbdt_model,
    X_train,
    **features_info,
    ax=ax,
    **common_params,
)
print(f"done in {time() - tic:.3f}s")
_ = display.figure_.suptitle("ICE and PDP representations", fontsize=16)

Computing partial dependence plots and individual conditional expectation...
done in 0.607s

ICE для признака температуры дает нам дополнительную информацию: некоторые линии ICE горизонтальны, а на других видно снижение зависимости от температуры выше 35 градусов Цельсия. Для признака влажности наблюдается похожая картина: некоторые линии ICE резко снижаются, когда влажность превышает 80%.

Не все линии ICE параллельны — это указывает на то, что модель обнаруживает взаимодействия между признаками. Повторим эксперимент, ограничив модель градиентного бустинга так, чтобы она не учитывала взаимодействия между признаками, с помощью параметра interaction_cst:

from sklearn.base import clone

interaction_cst = [[i] for i in range(X_train.shape[1])]
hgbdt_model_without_interactions = (
    clone(hgbdt_model)
    .set_params(histgradientboostingregressor__interaction_cst=interaction_cst)
    .fit(X_train, y_train)
)
print(f"Test R2 score: {hgbdt_model_without_interactions.score(X_test, y_test):.2f}")

Test R2 score: 0.38

_, ax = plt.subplots(ncols=2, figsize=(6, 4), sharey=True, constrained_layout=True)

features_info["centered"] = False
display = PartialDependenceDisplay.from_estimator(
    hgbdt_model_without_interactions,
    X_train,
    **features_info,
    ax=ax,
    **common_params,
)
_ = display.figure_.suptitle("ICE and PDP representations", fontsize=16)

Трехмерное представление

Построим такой же график частичной зависимости для взаимодействия двух признаков, но на этот раз в трех измерениях.

# unused but required import for doing 3d projections with matplotlib < 3.2
import mpl_toolkits.mplot3d  # noqa: F401
import numpy as np

from sklearn.inspection import partial_dependence

fig = plt.figure(figsize=(5.5, 5))

features = ("temp", "humidity")
pdp = partial_dependence(
    hgbdt_model, X_train, features=features, kind="average", grid_resolution=10
)
XX, YY = np.meshgrid(pdp["grid_values"][0], pdp["grid_values"][1])
Z = pdp.average[0].T
ax = fig.add_subplot(projection="3d")
fig.add_axes(ax)

surf = ax.plot_surface(XX, YY, Z, rstride=1, cstride=1, cmap=plt.cm.BuPu, edgecolor="k")
ax.set_xlabel(features[0])
ax.set_ylabel(features[1])
fig.suptitle(
    "PD of number of bike rentals on\nthe temperature and humidity GBDT model",
    fontsize=16,
)
# pretty init view
ax.view_init(elev=22, azim=122)
clb = plt.colorbar(surf, pad=0.08, shrink=0.6, aspect=10)
clb.ax.set_title("Partial\ndependence")
plt.show()

Пользовательские точки для анализа

Ни в одном из предыдущих примеров не указывалось, _какие именно_ точки используются для построения графиков частичной зависимости. По умолчанию мы используем процентили, определенные по входному набору данных. В некоторых случаях полезно задать точные точки, для которых нужно получить предсказания модели. Например, если пользователь хочет проверить поведение модели на данных, отличающихся от обучающих, или сравнить две модели, обученные на немного разных данных. Параметр custom_values позволяет передать значения, для которых нужно получить предсказания модели. Он переопределяет параметры grid_resolution и percentiles. Вернемся к рассмотренному выше примеру с градиентным бустингом, но зададим пользовательские значения.

print("Computing partial dependence plots with custom evaluation values...")
tic = time()
_, ax = plt.subplots(ncols=2, figsize=(6, 4), sharey=True, constrained_layout=True)

features_info = {
    "features": ["temp", "humidity"],
    "kind": "both",
}

display = PartialDependenceDisplay.from_estimator(
    hgbdt_model,
    X_train,
    **features_info,
    ax=ax,
    **common_params,
    # we set custom values for temp feature -
    # all other features are evaluated based on the data
    custom_values={"temp": np.linspace(0, 40, 10)},
)
print(f"done in {time() - tic:.3f}s")
_ = display.figure_.suptitle(
    (
        "Partial dependence of the number of bike rentals\n"
        "for the bike rental dataset with a gradient boosting"
    ),
    fontsize=16,
)

Computing partial dependence plots with custom evaluation values...
done in 0.728s

Общее время выполнения скрипта: (0 минут 30.470 секунд)

Download Jupyter notebook: plot_partial_dependence.ipynb

Download Python source code: plot_partial_dependence.py

Download zipped: plot_partial_dependence.zip

Похожие примеры

Расширенное построение графиков с частичной зависимостью

Создание признаков, связанных со временем

Основные нововведения scikit-learn 0.24

Основные нововведения scikit-learn 1.2