Отношения правдоподобия классов для оценки качества классификации

Анализ до и после тестирования

Предположим, что у нас есть популяция людей с физиологическими измерениями X, которые, как мы надеемся, могут служить косвенными биомаркерами заболевания, и фактическими показателями заболевания y (эталонная разметка). Большинство людей в популяции не страдают этим заболеванием, но у меньшинства (в данном случае примерно у 10%) оно есть:

from sklearn.datasets import make_classification

X, y = make_classification(n_samples=10_000, weights=[0.9, 0.1], random_state=0)
print(f"Percentage of people carrying the disease: {100 * y.mean():.2f}%")

Percentage of people carrying the disease: 10.37%

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

from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)

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

from sklearn.linear_model import LogisticRegression
from sklearn.metrics import class_likelihood_ratios

estimator = LogisticRegression().fit(X_train, y_train)
y_pred = estimator.predict(X_test)
pos_LR, neg_LR = class_likelihood_ratios(y_test, y_pred, replace_undefined_by=1.0)
print(f"LR+: {pos_LR:.3f}")

LR+: 12.617

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

Кросс-валидация отношений правдоподобия

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

import pandas as pd


def scoring(estimator, X, y):
    y_pred = estimator.predict(X)
    pos_lr, neg_lr = class_likelihood_ratios(y, y_pred, replace_undefined_by=1.0)
    return {"positive_likelihood_ratio": pos_lr, "negative_likelihood_ratio": neg_lr}


def extract_score(cv_results):
    lr = pd.DataFrame(
        {
            "positive": cv_results["test_positive_likelihood_ratio"],
            "negative": cv_results["test_negative_likelihood_ratio"],
        }
    )
    return lr.aggregate(["mean", "std"])

Сначала мы проверим модель LogisticRegression с гиперпараметрами по умолчанию, как и в предыдущем разделе.

from sklearn.model_selection import cross_validate

estimator = LogisticRegression()
extract_score(cross_validate(estimator, X, y, scoring=scoring, cv=10))

Мы подтверждаем, что модель полезна: апостериорные шансы в 12–20 раз превышают априорные.

Для сравнения рассмотрим фиктивную модель, которая выдаёт случайные предсказания с шансами, близкими к средней распространённости заболевания в обучающем наборе:

from sklearn.dummy import DummyClassifier

estimator = DummyClassifier(strategy="stratified", random_state=1234)
extract_score(cross_validate(estimator, X, y, scoring=scoring, cv=10))

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

Другой вариант фиктивной модели — всегда предсказывать наиболее частый класс, которым в данном случае является «нет заболевания».

estimator = DummyClassifier(strategy="most_frequent")
extract_score(cross_validate(estimator, X, y, scoring=scoring, cv=10))

/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: No samples were predicted for the positive class and `positive_likelihood_ratio` is set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)

Отсутствие положительных предсказаний означает, что не будет ни истинно положительных, ни ложно положительных результатов, поэтому LR+ не определено. Это ни в коем случае не следует интерпретировать как бесконечное LR+ (то есть как идеальное выявление классификатором положительных случаев). В такой ситуации функция class_likelihood_ratios возвращает nan и по умолчанию выдаёт предупреждение. Действительно, значение LR- помогает нам отбросить эту модель.

Похожий сценарий может возникнуть при кросс-валидации сильно несбалансированных данных с небольшим числом образцов: в некоторых фолдах не будет образцов с заболеванием, поэтому при тестировании они не дадут ни истинно положительных, ни ложно отрицательных результатов. Математически это приводит к бесконечному LR+, что также не следует интерпретировать как идеальное выявление классификатором положительных случаев. Такое событие повышает дисперсию оценённых отношений правдоподобия, но его всё ещё можно интерпретировать как увеличение апостериорных шансов наличия заболевания.

estimator = LogisticRegression()
X, y = make_classification(n_samples=300, weights=[0.9, 0.1], random_state=0)
extract_score(cross_validate(estimator, X, y, scoring=scoring, cv=10))

/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: `positive_likelihood_ratio` is ill-defined and set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: `positive_likelihood_ratio` is ill-defined and set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: `positive_likelihood_ratio` is ill-defined and set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: `positive_likelihood_ratio` is ill-defined and set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)
/home/circleci/project/sklearn/utils/_param_validation.py:218: UndefinedMetricWarning: `positive_likelihood_ratio` is ill-defined and set to `np.nan`. Use the `replace_undefined_by` param to
  return func(*args, **kwargs)

Инвариантность относительно распространённости

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

Здесь мы обучаем базовую модель LogisticRegression на исследовании «случай — контроль» с распространённостью 50%. Затем модель оценивается на популяциях с различной распространённостью. Мы используем функцию make_classification, чтобы процесс генерации данных всегда оставался одинаковым, как показано на графиках ниже. Метка 1 соответствует положительному классу «заболевание», тогда как метка 0 обозначает класс «нет заболевания».

from collections import defaultdict

import matplotlib.pyplot as plt
import numpy as np

from sklearn.inspection import DecisionBoundaryDisplay

populations = defaultdict(list)
common_params = {
    "n_samples": 10_000,
    "n_features": 2,
    "n_informative": 2,
    "n_redundant": 0,
    "random_state": 0,
}
weights = np.linspace(0.1, 0.8, 6)
weights = weights[::-1]

# fit and evaluate base model on balanced classes
X, y = make_classification(**common_params, weights=[0.5, 0.5])
estimator = LogisticRegression().fit(X, y)
lr_base = extract_score(cross_validate(estimator, X, y, scoring=scoring, cv=10))
pos_lr_base, pos_lr_base_std = lr_base["positive"].values
neg_lr_base, neg_lr_base_std = lr_base["negative"].values

Теперь мы покажем границу принятия решения для каждого уровня распространённости. Обратите внимание, что мы строим график лишь для части исходных данных, чтобы лучше рассмотреть линейную границу принятия решения модели.

fig, axs = plt.subplots(nrows=3, ncols=2, figsize=(15, 12))

for ax, (n, weight) in zip(axs.ravel(), enumerate(weights)):
    X, y = make_classification(
        **common_params,
        weights=[weight, 1 - weight],
    )
    prevalence = y.mean()
    populations["prevalence"].append(prevalence)
    populations["X"].append(X)
    populations["y"].append(y)

    # down-sample for plotting
    rng = np.random.RandomState(1)
    plot_indices = rng.choice(np.arange(X.shape[0]), size=500, replace=True)
    X_plot, y_plot = X[plot_indices], y[plot_indices]

    # plot fixed decision boundary of base model with varying prevalence
    disp = DecisionBoundaryDisplay.from_estimator(
        estimator,
        X_plot,
        response_method="predict",
        alpha=0.5,
        ax=ax,
    )
    scatter = disp.ax_.scatter(X_plot[:, 0], X_plot[:, 1], c=y_plot, edgecolor="k")
    disp.ax_.set_title(f"prevalence = {y_plot.mean():.2f}")
    disp.ax_.legend(*scatter.legend_elements())

Определим функцию для бутстрэп-оценивания.

def scoring_on_bootstrap(estimator, X, y, rng, n_bootstrap=100):
    results_for_prevalence = defaultdict(list)
    for _ in range(n_bootstrap):
        bootstrap_indices = rng.choice(
            np.arange(X.shape[0]), size=X.shape[0], replace=True
        )
        for key, value in scoring(
            estimator, X[bootstrap_indices], y[bootstrap_indices]
        ).items():
            results_for_prevalence[key].append(value)
    return pd.DataFrame(results_for_prevalence)

Оценим базовую модель для каждого уровня распространённости с помощью бутстрэп-оценивания.

results = defaultdict(list)
n_bootstrap = 100
rng = np.random.default_rng(seed=0)

for prevalence, X, y in zip(
    populations["prevalence"], populations["X"], populations["y"]
):
    results_for_prevalence = scoring_on_bootstrap(
        estimator, X, y, rng, n_bootstrap=n_bootstrap
    )
    results["prevalence"].append(prevalence)
    results["metrics"].append(
        results_for_prevalence.aggregate(["mean", "std"]).unstack()
    )

results = pd.DataFrame(results["metrics"], index=results["prevalence"])
results.index.name = "prevalence"
results

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

fig, (ax1, ax2) = plt.subplots(nrows=1, ncols=2, figsize=(15, 6))
results["positive_likelihood_ratio"]["mean"].plot(
    ax=ax1, color="r", label="extrapolation through populations"
)
ax1.axhline(y=pos_lr_base + pos_lr_base_std, color="r", linestyle="--")
ax1.axhline(
    y=pos_lr_base - pos_lr_base_std,
    color="r",
    linestyle="--",
    label="base model confidence band",
)
ax1.fill_between(
    results.index,
    results["positive_likelihood_ratio"]["mean"]
    - results["positive_likelihood_ratio"]["std"],
    results["positive_likelihood_ratio"]["mean"]
    + results["positive_likelihood_ratio"]["std"],
    color="r",
    alpha=0.3,
)
ax1.set(
    title="Positive likelihood ratio",
    ylabel="LR+",
    ylim=[0, 5],
)
ax1.legend(loc="lower right")

ax2 = results["negative_likelihood_ratio"]["mean"].plot(
    ax=ax2, color="b", label="extrapolation through populations"
)
ax2.axhline(y=neg_lr_base + neg_lr_base_std, color="b", linestyle="--")
ax2.axhline(
    y=neg_lr_base - neg_lr_base_std,
    color="b",
    linestyle="--",
    label="base model confidence band",
)
ax2.fill_between(
    results.index,
    results["negative_likelihood_ratio"]["mean"]
    - results["negative_likelihood_ratio"]["std"],
    results["negative_likelihood_ratio"]["mean"]
    + results["negative_likelihood_ratio"]["std"],
    color="b",
    alpha=0.3,
)
ax2.set(
    title="Negative likelihood ratio",
    ylabel="LR-",
    ylim=[0, 0.5],
)
ax2.legend(loc="lower right")

plt.show()

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

Download Jupyter notebook: plot_likelihood_ratios.ipynb

Download Python source code: plot_likelihood_ratios.py

Download zipped: plot_likelihood_ratios.zip

Связанные примеры

SVM: взвешенные образцы

Настройка порога функции принятия решений после обучения

Важность масштабирования признаков

Объединение предикторов с помощью стекинга