Классификация#
Рассмотрим задачу классификации точек на плоскости \(Ox_1x_2\). Любую точку на плоскости можно охарактеризовать её координатами \(x_1\) и \(x_2\). Эти координаты и примем за их признаковое описание.
Бинарная классификация#
Для простоты будем считать, что все точки принадлежат одному из двух классов, а также что имеется выборка из \(n\) точек, для которых известен этот самый класс. Назовем первый класс отрицательным, второй класс — положительным. Поставим в соответствие всем точкам отрицательного класса значение \(y_\mathrm{true} = -1\), а всем точкам положительного класса \(y_\mathrm{true}= + 1\). Тогда известную выборку можно представить в виде следующей таблицы.
\(x_1\) |
\(x_2\) |
\(y_\mathrm{true}\) |
---|---|---|
\(x_1^1\) |
\(x_2^1\) |
\(y_\mathrm{true}^1\) |
\(\cdots\) |
\(\cdots\) |
\(\cdots\) |
\(x_1^n\) |
\(x_2^n\) |
\(y_\mathrm{true}^n\) |
Здесь верхний индекс соответствует номеру точки.
Необходимо на основе этих данных научиться классифицировать произвольную точку на плоскости, т.е. построить отображение \(f\colon \mathbb{R}^2\to\{-1, 1\}\), которое бы хорошо классифицировало не только точки из известной выборки, но и остальные точки из \(\mathbb{R}\).
Пока не будем уточнять, какой смысл вкладывается в “хорошо классифицировать”.
import os import numpy as np import pandas as pd from matplotlib import pyplot as plt import seaborn as sns folder = os.path.join("..", "_static", "lecture_specific", "ml") plt.rcParams.update({"font.size": 16})
def plot_decision_boundary(classifier, data, ax, title=None, n_points=300): x_min, x_max = data["x1"].min(), data["x1"].max() y_min, y_max = data["x2"].min(), data["x2"].max() delta_x = (x_max - x_min) / 10 delta_y = (x_max - x_min) / 10 x = np.linspace(x_min - delta_x, x_max + delta_x, n_points) y = np.linspace(y_min - delta_y, y_max + delta_y, n_points) xx, yy = np.meshgrid(x, y) df = pd.DataFrame({ "x1": xx.reshape(-1), "x2": yy.reshape(-1) }) zz = classifier.predict(df) ax.contourf(xx, yy, zz.reshape(n_points, n_points), levels=1, cmap="Set2") sns.scatterplot( data=data, x='x1', y='x2', hue='y', s=100, palette="bwr", hue_norm=(-1, 1), ax=ax ) ax.set_title(title)
Случай линейно разделимой выборки#
Один из самых простых случаев — когда можно провести прямую (гиперплоскость в общем случае) таким образом, что все точки первого класса окажутся по одну сторону от неё, а все точки второго класса — по другую.
Ниже приводится пример такой выборки.
path = os.path.join(folder, "linearly_separable.csv") data = pd.read_csv(path) fig, ax = plt.subplots(figsize=(10, 7)) sns.scatterplot( data=data, x='x1', y='x2', hue='y', s=100, palette="bwr", hue_norm=(-1, 1), ax=ax ) data.head()
x1 | x2 | y | |
---|---|---|---|
0 | 1.228312 | -0.757178 | -1 |
1 | 0.698409 | -1.380295 | -1 |
2 | 2.548817 | 2.502258 | 1 |
3 | 0.573579 | -1.352979 | -1 |
4 | 0.585900 | -1.337457 | -1 |
Из построенной диаграммы рассеяния отчетливо видно, что классы можно отделить друг от друга прямой.
В общем случае прямую можно задать уравнением
\[
\omega_1 x_1 + \omega_2 x_2 + b = 0,
\]
а определить с какой стороны от неё находится точка с координатами \((x_1^*, x_2^*)\) можно по знаку величины \(a(x_1^*, x_2^*) = \omega_1 x_1^* + \omega_2 x_2^* + b\).
Значит можно задать классификатор
\[
y_\mathrm{pred}(x_1, x_2) = \mathrm{sgn}\,(a(x_1, x_2)),
\]
где \(\mathrm{sgn}\) — функция сигнум.
Чисто теоретически можно подобрать коэффициенты \(\vec{\omega} = (\omega_1, \omega_2)\) и \(b\) в ручную. В ячейке ниже демонстрируется такой классификатор с \(\omega_1=-1\), \(\omega_2=1\) и \(b=0.3\).
class ManualLinearClassifier: def __init__(self, weights, bias): self.weights = np.array(weights) self.bias = bias def predict(self, x): y = (self.weights @ x.T + self.bias).to_numpy() return np.sign(y) classifier = ManualLinearClassifier([-1, 1], 0.3) fig, ax = plt.subplots(figsize=(10, 7)) plot_decision_boundary(classifier, data, ax)
Даже такой простой классификатор успешно отделил все точки. Однако если бы это была многомерная задача, то подбирать коэффициенты вручную оказалось бы гораздо сложнее. Хочется каким-то образом обучить эти коэффициенты из данных.
Это можно сделать, опираясь на уже рассмотренную задачу линейной регрессии. Подберем коэффициенты \((\vec\omega, b)\) таким образом, чтобы минимизировать среднее квадратичное отклонение
$\(
L(X_\mathrm{train}, \vec\omega, b) = \dfrac{1}{n}\sum_{i=1}^{n}(\omega_1 x_1^i + \omega_2 x_2^i + b — y^i)^2 \sim \min_{\substack{\vec\omega\in\mathbb{R}^2 \\ b\in\mathbb{R}}}.
\)$
Код в ячейке ниже расширяет класс LinearRegression
таким образом, чтобы метод predict
возвращал не действительное число из \(\mathbb{R}\), а только его знак.
from sklearn.linear_model import LinearRegression class RegressionClassifier(LinearRegression): def predict(self, x): y = super().predict(x) return np.sign(y) regressor = RegressionClassifier() regressor.fit(data[["x1", "x2"]], data["y"]) fig, ax = plt.subplots(figsize=(10, 7)) plot_decision_boundary(regressor, data, ax)
Видим, что удалось разумно разделить классы на основе линейной регрессии. Тем не менее обычно для решения задач классификации используют другие функции потерь. Одна из самых распространенных функций потерь — перекрестная энтропия, которая в нашей задаче имеет вид
\[
L(p, y) = — (y \ln p + (1-y)\ln(1-p)).
\]
Здесь \(p = \sigma(\omega_1 x_1 + \omega_2 x_2 + b) \in[0, 1]\) интерпретируется в качестве вероятности того, что точка \((x_1, x_2)\) принадлежит положительному классу, а \((1-p)\) — в качестве вероятности того, что эта точка принадлежит отрицательному классу. Функция \(\sigma\) — сигмоида, значение которой в точке \(x\in\mathbb{R}\) определяется формулой
\[
\sigma(x) = \dfrac{1}{1 + e^{-x}},
\]
которая принимает значения в интервале \((0, 1)\), монотонно растет и при этом
\[\begin{split}
\begin{gathered}
\lim\limits_{x\to-\infty}\sigma(x)=0, \\
\lim\limits_{x\to+\infty}\sigma(x)=1.
\end{gathered}
\end{split}\]
Это и позволяет интерпретировать значения этой функции в качестве вероятности.
Код в ячейке ниже строит графики сигмоиды, перекрестной энтропии и среднеквадратичной функции потерь.
from scipy.special import expit as sigmoid from sklearn.metrics import log_loss as cross_entropy_loss from sklearn.metrics import mean_squared_error fig, axs = plt.subplots(figsize=(15, 5), ncols=3) # sigmoid x = np.linspace(-5, 5) y = sigmoid(x) axs[0].plot(x, y) axs[0].set_title("Sigmoid") axs[0].set_xlabel("$x$") axs[0].set_ylabel(r"$\sigma(x)$") # cross entropy loss probabilities = np.linspace(0.01, 0.99) l0 = [cross_entropy_loss([0], [p], labels=[0, 1]) for p in probabilities] l1 = [cross_entropy_loss([1], [p], labels=[0, 1]) for p in probabilities] axs[1].plot(probabilities, l0, probabilities, l1) axs[1].set_title("Cross entropy loss") axs[1].set_xlabel("$p$") axs[1].set_ylabel(r"$L(y, p)$") axs[1].legend(["$y=0$", "$y=1$"]) # mean squared error y_predictions = np.linspace(-3, 3) l1 = [mean_squared_error([-1], [y_p]) for y_p in y_predictions] l2 = [mean_squared_error([1], [y_p]) for y_p in y_predictions] axs[2].plot(y_predictions, l1, y_predictions, l2) axs[2].set_title("Mean squared error") axs[2].set_xlabel(r"$y_\mathrm{pred}$") axs[2].set_ylabel(r"$L(y_\mathrm{true}, y_\mathrm{pred})$") axs[2].legend([r"$y_\mathrm{true}=-1$", r"$y_\mathrm{true}=+1$"]) _ = axs[2].set_xticks([-1, 1])
Функция перекрестной энтропии \(L(y, p)\) определена на интервале \((0, 1)\) и в зависимости от значения \(y\) монотонно растёт на этом интервале от 0 до бесконечности (\(y=0\)) или монотонно убывает от бесконечности до 0 (\(y=1\)).
Значение \(p\) интерпретируется в виде предсказанной вероятности какого-то события (в задаче классификации это вероятность того, что классифицируемый объект принадлежит классу с \(y=1\)), а \(y\) считается истинной вероятностью. В связи с этим и рядом других факторов обученные с помощью перекрестной энтропии модели считают вероятностными.
Саму модель, использующую такую функцию потерь, называют логистической регрессией, а метод linear_model.LogisticRegression из библиотеки scikit learn
её реализует.
from sklearn.linear_model import LogisticRegression log_classifier = LogisticRegression() log_classifier.fit(data[["x1", "x2"]], data["y"]) fig, ax = plt.subplots(figsize=(10, 7)) plot_decision_boundary(log_classifier, data, ax)
Получили другую разделяющую прямую, а значит и другой классификатор. В общем случае нельзя однозначно сказать, какая из функций потерь лучше подойдет для конкретной задачи.
Тем не менее продемонстрируем недостаток функции потерь MSE
на чуть более сложных данных. Добавим к данным ещё одну точку положительного класса в глубину этого класса и посмотрим, как отреагирует классификаторы с разными функциями потерь.
outlier = pd.DataFrame({"x1": -3, "x2": 7, "y": 1}, index=[100]) data_with_outlier = pd.concat([data, outlier]) X = data_with_outlier[["x1", "x2"]] y = data_with_outlier["y"] log_classifier.fit(X, y) regressor.fit(X, y) fig, axs = plt.subplots(figsize=(15, 7), ncols=2) plot_decision_boundary(regressor, data_with_outlier, axs[0], title="MSE loss") plot_decision_boundary(log_classifier, data_with_outlier, axs[1], title="Cross entropy loss") for ax in axs: ax.scatter(outlier.x1, outlier.x2, s=200, c="r", label="outlier", marker="*") ax.legend()
Классификатор с MSE
функцией потерь резко отреагировал на добавление новой точки, хотя она была бы верно классифицирована и до этого.
Чтобы объяснить этот эффект, рассмотрим величину \(a(x_1^*, x_2^*) = \omega_1 x_1^* + \omega_2 x_2^* + b\), по знаку которой мы оцениваем класс точки \((x_1^*, x_2^*)\). Вектор \(\vec\omega = (\omega_1, \omega_2)\) задаёт нормаль к прямой \(\omega_1 x_1 + \omega_2 x_2+ b = 0\), а значит величина \(a\) определяет насколько далеко точка располагается от разделяющей прямой: чем больше абсолютное значение этой величины, тем дальше эта точка от прямой.
Для заданной прямой и для любой точки \((x_1, x_2)\) из выборки можно ввести отступ \(M(x_1, x_2, y_\mathrm{true}) = y_\mathrm{true} a(x_1, x_2)\). Величина \(M\) положительна для верно классифицируемых точек и отрицательная для ошибочно классифицируемых точек, а абсолютное значение этой величины можно интерпретировать в качестве ступени уверенности классификатора, т.к. чем выше это значение, тем дальше точка находится от разделяющей прямой, а значит тем устойчивее она классифицируется: если разделяющую прямую немного изменить, то точка все равно останется с прежней стороны. Исходя из этого, можно ввести четыре категории точек:
-
точки с большим положительным значением \(M\), которые устойчиво классифицируются верно;
-
точки с небольшим положительным значением \(M\), которые неустойчиво классифицируются верно;
-
точки с небольшим отрицательным значением \(M\), которые неустойчиво классифицируются неверно;
-
точки с большим отрицательным значением \(M\), которые устойчиво классифицируются неверно;
Возвращаясь к функцией потерь MSE
, заметим, что она штрафует за любое отклонения отступа \(M\) от единицы, т.е. точки первой категории тоже штрафуются. Если же обратить внимание на функцию перекрестной энтропии, то она равна 0, только если \(\sigma(M)=0\), т.е. \(M=\infty\). Таким образом с точки зрения этой функции потерь чем выше отступ, тем лучше.
Назовем зазором классификатора расстояние от разделяющей поверхности до ближайшей точки из выборки \(\min\limits_{i=1}^n M(x_1^i, x_2^i, y^i)\). Метод опорных векторов ищет классификатор с максимальным зазором. Функция sklearn.svm.LinearSVC позволяет найти такой линейный классификатор.
from sklearn.svm import LinearSVC svm_classifier = LinearSVC() svm_classifier.fit(data[["x1", "x2"]], data["y"]) fig, axs = plt.subplots(figsize=(15, 7), ncols=2) plot_decision_boundary(svm_classifier, data, axs[0], title="Support Vector Machine") plot_decision_boundary(log_classifier, data, axs[1], title="Logistic Regression")
Случай не разделимой линейно выборки#
Рассмотрим другую выборку, которую линейно разделить не выйдет.
path = os.path.join(folder, "circles.csv") data = pd.read_csv(path) fig, ax = plt.subplots(figsize=(10, 7)) sns.scatterplot( data=data, x='x1', y='x2', hue='y', s=100, palette="bwr", hue_norm=(-1, 1), ax=ax )
<AxesSubplot:xlabel='x1', ylabel='x2'>
Исходя из диаграммы рассеяния очевидно, что разделить эту выборку прямой не выйдет.
Конструирование признаков.#
При решении задачи мы решили, что признаками точки будем считать координаты точки в декартовой системе координат. Такой выбор признакового описания точки можно считать произвольным: почему, например, не в полярной системе координат; почему только \(x_1\) и \(x_2\).
В ячейке ниже определяется классификатор RadiusLogisticClassifier
, который при обучении и оценивании класса вместо исходных признаков x1
и x2
использует признак радиуса r
.
class RadiusLogisticClassifier(LogisticRegression): def get_radii(self, data): new_data = pd.DataFrame({ "r": np.sqrt(data["x1"] ** 2 + data["x2"] ** 2) }) return new_data def fit(self, data, targets): new_data = self.get_radii(data) super().fit(new_data, targets) def predict(self, data): new_data = self.get_radii(data) return super().predict(new_data) classifier = RadiusLogisticClassifier() classifier.fit(data[["x1", "x2"]], data["y"]) fig, ax = plt.subplots(figsize=(10, 7)) plot_decision_boundary(classifier, data, ax) classifier.coef_
Практика показывает, что удачный набор признаков значительно упрощает задачу обучения алгоритма. Процесс составления признакового описания объектов называют конструированием признаков.
Нелинейные классификаторы#
Поиск подходящего признакового описания требует экспертизы в решаемой задаче. Далеко не всегда очевидно, какие признаки будут полезны для классификации, а какие нет. В качестве альтернативы можно пробовать строить нелинейные классификаторы, которые по своей природе позволяет строить нелинейную разделяющую поверхность.
В ячейке ниже на четырех выборках, три из которых линейно неразделимы, тестируются обучаются три нелинейных классификатора:
-
метод ближайших соседей;
-
деревья решений;
-
нелинейный метод опорных векторов.
from sklearn.neighbors import KNeighborsClassifier from sklearn.tree import DecisionTreeClassifier from sklearn.svm import SVC datasets = [ pd.read_csv(os.path.join(folder, "linearly_separable.csv")), pd.read_csv(os.path.join(folder, "circles.csv")), pd.read_csv(os.path.join(folder, "moons.csv")), pd.read_csv(os.path.join(folder, "spirals.csv")), pd.read_csv(os.path.join(folder, "quadrant.csv")) ] models = { "K nearest neighbors": KNeighborsClassifier(), "Decision tree": DecisionTreeClassifier(), "Support vector machine": SVC() } n_models = len(models) n_datasets = len(datasets) fig, axs = plt.subplots(figsize=(16, 16), nrows=n_datasets, ncols=n_models, layout="tight") for i, (model_name, classifier) in enumerate(models.items()): for j, data in enumerate(datasets): classifier.fit(data[["x1", "x2"]], data["y"]) plot_decision_boundary(classifier, data, axs[j, i]) axs[j, i].xaxis.set_visible(False) axs[j, i].yaxis.set_visible(False) axs[0, i].set_title(model_name)
Выбор метода тоже является весьма не тривиальный задачей. Нередко пробуют как можно больше разных и выбирают тот, который демонстрирует лучшие показатели.
Метрики качества#
До сих пор мы обучали модели на всей выборке, а оценивали качество обученной модели визуально. В реальных ситуациях количество признаков обычно превышает два и визуально оценить качество модели невозможно.
path = os.path.join(folder, "balanced.csv") data = pd.read_csv(path) print(data.info())
<class 'pandas.core.frame.DataFrame'> RangeIndex: 1000 entries, 0 to 999 Data columns (total 11 columns): # Column Non-Null Count Dtype --- ------ -------------- ----- 0 f1 1000 non-null float64 1 f2 1000 non-null float64 2 f3 1000 non-null float64 3 f4 1000 non-null float64 4 f5 1000 non-null float64 5 f6 1000 non-null float64 6 f7 1000 non-null float64 7 f8 1000 non-null float64 8 f9 1000 non-null float64 9 f10 1000 non-null float64 10 y 1000 non-null int64 dtypes: float64(10), int64(1) memory usage: 86.1 KB None
Например, в примере выше у каждого объекта 10 признаков. В таких случаях используются различные метрики, значение которых вычисляется на выборках.
В ячейке ниже загруженная выборка делится сначала на две таблицы X
и Y
, где X
— таблица признаков, на основе которых классифицируются точки, а Y
— номера классов этих точек.
Затем каждая из этих таблиц делится на тренировочную и валидационную части.
from sklearn.model_selection import train_test_split Y = data["y"] X = data.drop("y", axis=1) X_train, X_test, Y_train, Y_test = train_test_split(X, Y, train_size=0.7)
В ячейке ниже инициализируются 9 классификаторов, среди которых несколько классификаторов из библиотеки scikit-learn
, а несколько наивных классификаторов:
-
ConstantClassifier
— классификатор, который всегда выдает один и тот же класс, вне зависимости от входных данных; -
RandomClassifier
— классификатор, который пытается угадать класс, случайным образом выбирая номер класса.
Далее качество этих классификаторов будет оценивать с помощью разных метрик.
from abc import ABC, abstractmethod from sklearn.ensemble import RandomForestClassifier rng = np.random.default_rng() class CustomClassifier: """ abc for custom classifiers all of them should have methods fit and predict so they can be used in the same context as sklearn classifiers """ def fit(self, *args, **kwargs): pass @abstractmethod def predict(self, x): raise NotImplementedError() class ConstantClassifier(CustomClassifier): def __init__(self, constant=0): self.constant = constant def predict(self, x): return self.constant * np.ones(shape=len(x)) class RandomClassifier(CustomClassifier): def predict(self, x): return rng.choice([0, 1], size=len(x), replace=True) classifiers = { "always 0": ConstantClassifier(0), "always 1": ConstantClassifier(1), "random": RandomClassifier(), "logistic regression": LogisticRegression(), "linear SVM": LinearSVC(max_iter=10000), "nonlinear SVM": SVC(), "kNN": KNeighborsClassifier(), "tree": DecisionTreeClassifier(), "random forest": RandomForestClassifier(), } for classifier in classifiers.values(): classifier.fit(X_train, Y_train)
Ошибки первого и второго рода#
Для произвольного бинарного классификатора на любой выборке можно вычислить следующие значения:
-
TP
(true positive
) — количество корректно классифицированных объектов положительного класса; -
TN
(true negative
) — количество корректно классифицированных объектов отрицательного класса; -
FP
(false positive
) — количество объектов отрицательного класса, для которых классификатор выдал положительный класс; -
FN
(false negative
) — количество объектов положительного класса, для которых классификатор выдал отрицательный класс.
\(y_\mathrm{true} = -1\) |
\(y_\mathrm{true} = +1\) |
|
---|---|---|
\(y_\mathrm{pred} = -1\) |
|
|
\(y_\mathrm{pred} = +1\) |
|
|
На главную диагональ этой таблицы попадают корректно классифицируемые объекты. Все ошибки сконцентрированы на побочной диагонали. При этом ложноположительные ошибки называют ошибками первого рода, а ложно отрицательные ошибки — ошибками второго рода.
Метод sklearn.metric.confusion_matrix принимает на вход массивы истинных классов и выхода классификатора, а возвращает матрицу, схожую с приведенной таблицей. Методом seaborn.heatmap удобно эти данные визуализировать.
from sklearn.metrics import confusion_matrix def plot_confusion_matrix(ax, y_true, y_pred, title): matrix = confusion_matrix(y_true, y_pred) sns.heatmap(matrix, annot=True, ax=ax, cbar=False, cmap="coolwarm", fmt="d") ax.set_xlabel(r"$y_\mathrm{pred}$") ax.set_ylabel(r"$y_\mathrm{true}$") ax.set_title(title) fig, axs = plt.subplots(figsize=(12, 12), nrows=3, ncols=3, layout="tight") for ax, (name, clf) in zip(axs.flat, classifiers.items()): y_pred = clf.predict(X_test) plot_confusion_matrix(ax, Y_test, y_pred, name)
Такого рода таблицы очень полезны для анализа того, на каких классах чаще возникают ошибки.
В качестве недостатка можно отметить, что в одной такой таблице минимум 4 числа, и не совсем понятно, как сравнивать качество моделей между собой на основе них.
Точность#
Самая простая метрика — точность классификации, которая соответствует доли корректно классифицированных точек.
Используя ранее введенные обозначения формулу для точности можно записать в виде
$\(
\mathrm{accuracy} = \dfrac{\mathrm{TP} + \mathrm{TN}}{\mathrm{TP} + \mathrm{TN} + \mathrm{FP} + \mathrm{FN}}.
\)$
Её преимущество относительно предыдущих таблиц заключается в том, что это одно всего одно число. Вычислить точность можно методом sklearn.metrics.accuracy_score.
from sklearn.metrics import accuracy_score def evaluate_accuracy(classifier, X, Y_true): Y_pred = classifier.predict(X) return accuracy_score(Y_true, Y_pred) train_accuracies = {} test_accuracies = {} for model_name, classifier in classifiers.items(): train_accuracies[model_name] = evaluate_accuracy(classifier, X_train, Y_train) test_accuracies[model_name] = evaluate_accuracy(classifier, X_test, Y_test) metrics = pd.DataFrame({ "train accuracy": pd.Series(train_accuracies), "test accuracy": pd.Series(test_accuracies) }) metrics
train accuracy | test accuracy | |
---|---|---|
always 0 | 0.478571 | 0.553333 |
always 1 | 0.521429 | 0.446667 |
random | 0.510000 | 0.496667 |
logistic regression | 0.834286 | 0.833333 |
linear SVM | 0.834286 | 0.826667 |
nonlinear SVM | 0.940000 | 0.946667 |
kNN | 0.958571 | 0.916667 |
tree | 1.000000 | 0.800000 |
random forest | 1.000000 | 0.896667 |
Недостаток такой метрики тот же, что и её преимущество: это всего лишь одно число. Например, вычислив точность классификации, мы потеряли информацию о том, какого рода ошибки чаще встречаются.
Несбалансированные выборки#
Если классы представлены в выборке неравномерно, то метрика точности может оказаться практически бесполезной. В качестве демонстрации рассмотрим следующую ситуацию. Пусть 99% выборки — объекты первого отрицательного класса, а оставшийся 1% — объекты положительного класса. Тогда выдающий всегда отрицательный класс классификатор будет иметь точность 99%, т.к. он будет ошибаться только на объектах положительного класса, коих всего 1%.
Несбалансированные выборки часто возникают при поиска аномалий, медицинских данных, обнаружении атак и т.д. В таких ситуациях принято рассматривать другие метрики. Часто смотрят на метрики precision
(точность) и recall
(полнота), которые определяются формулами
\[\begin{split}
\begin{aligned}
\mathrm{precision}& = \dfrac{\mathrm{TP}}{\mathrm{TP} + \mathrm{FP}}, \\
\mathrm{recall}& = \dfrac{\mathrm{TP}}{\mathrm{TP} + \mathrm{FN}}.
\end{aligned}
\end{split}\]
Смысл этих метрик легко понять, если провести аналогию с медицинским тестированием на редкое заболевание. Пусть положительный класс — множество людей действительно больных заболеванием, а отрицательный класс — остальные люди. Тест может давать положительный (болен) и отрицательный результат (не болен). Тогда
-
полнота (
recall
или чувствительность в медицине) — доля больных, которым тест выдал корректный результат. Чем выше чувствительность (полнота), тем вероятнее этот тест обнаружит больного. -
точность (
precision
) — доля действительно больных пациентов среди пациентов с положительным тестом. Чем выше точность теста, тем выше вероятность того, что пациент действительно болен, если его тест дал положительный результат. -
вместо точности в медицине используют метрику специфичности теста (\(\mathrm{specificity} = \frac{\mathrm{TN}}{\mathrm{FP} + \mathrm{TN}}\)), которая характеризует в какой степени тест специфичен именно к данному заболеванию.
В ячейке ниже загружается выборка, в которой соотношение положительного класса к отрицательному 1:9.
path = os.path.join(folder, "unbalanced.csv") data = pd.read_csv(path) classes_representation = data["y"].value_counts() fig, ax = plt.subplots(figsize=(5, 5)) labels = ["Negative", "Positive"] ax.pie(classes_representation, labels=labels, autopct="%.0f%%") Y = data["y"] X = data.drop("y", axis=1) X_train, X_test, Y_train, Y_test = train_test_split(X, Y, train_size=0.7)
В следующей ячейке вычисляются accuracy
, recall
и precision
для того же списка моделей, что фигурировал при обучении на сбалансированной выборке.
-
для вычисления
precision
используется sklearn.metrics.precision_score; -
для вычисления
recall
используется sklearn.metrics.recall_score.
from sklearn.metrics import precision_score, recall_score def evaluate_precision(classifier, X, Y_true): Y_pred = classifier.predict(X) return precision_score(Y_true, Y_pred, zero_division=1) def evaluate_recall(classifier, X, Y_true): Y_pred = classifier.predict(X) return recall_score(Y_true, Y_pred) accuracies = {} precisions = {} recalls = {} for model_name, classifier in classifiers.items(): classifier.fit(X_train, Y_train) accuracies[model_name] = evaluate_accuracy(classifier, X_test, Y_test) precisions[model_name] = evaluate_precision(classifier, X_test, Y_test) recalls[model_name] = evaluate_recall(classifier, X_test, Y_test) metrics = pd.DataFrame({ "accuracy": pd.Series(accuracies), "precision": pd.Series(precisions), "recall": pd.Series(recalls), }) metrics
accuracy | precision | recall | |
---|---|---|---|
always 0 | 0.891667 | 1.000000 | 0.000000 |
always 1 | 0.108333 | 0.108333 | 1.000000 |
random | 0.523333 | 0.111486 | 0.507692 |
logistic regression | 0.940000 | 0.853659 | 0.538462 |
linear SVM | 0.938333 | 0.868421 | 0.507692 |
nonlinear SVM | 0.951667 | 0.928571 | 0.600000 |
kNN | 0.948333 | 0.840000 | 0.646154 |
tree | 0.930000 | 0.676923 | 0.676923 |
random forest | 0.953333 | 0.893617 | 0.646154 |
Из этой таблицы видно, что выдающий всегда отрицательный класс классификатор always 0
имеет accuracy
равную 89.5%, но нулевую полноту. Построим диаграмму рассеивания для показателей precision
и recall
по типам моделей.
metrics["model"] = metrics.index fig, ax = plt.subplots(figsize=(10, 10)) sns.scatterplot(data=metrics, x="precision", y="recall", hue="model", style="model", s=200, ax=ax)
<AxesSubplot:xlabel='precision', ylabel='recall'>
Предыдущий пост см. здесь.
Проверка статистических гипотез
Для статистиков и исследователей данных проверка статистической гипотезы представляет собой формальную процедуру. Стандартный подход к проверке статистической гипотезы подразумевает определение области исследования, принятие решения в отношении того, какие переменные необходимы для измерения предмета изучения, и затем выдвижение двух конкурирующих гипотез. Во избежание рассмотрения только тех данных, которые подтверждают наши субъективные оценки, исследователи четко констатируют свою гипотезу заранее. Затем, основываясь на данных, они применяют выборочные статистики с целью подтвердить либо отклонить эту гипотезу.
Проверка статистической гипотезы подразумевает использование тестовой статистики, т.е. выборочной величины, как функции от результатов наблюдений. Тестовая статистика (test statistic) — это вычисленная из выборочных данных величина, которая используется для оценивания прочности данных, подтверждающих нулевую статистическую гипотезу и служит для выявления меры расхождения между эмпирическими и гипотетическими значениями. Конкретные методы проверки называются тестами, например, z-тест, t-тест (соответственно z-тест Фишера, t-тест Студента) и т.д. в зависимости от применяемых в них тестовых статистик.
Примечание. В отечественной статистической науке используется «туманный» термин «статистика критерия». Туманный потому здесь мы снова наблюдаем мягкую подмену: вместо теста возникает критерий. Если уж на то пошло, то критерий — это принцип или правило. Например, выполняя z-тест, t-тест и т.д., мы соответственно используем z-статистику, t-статистику и т.д. в правиле отклонения гипотезы. Это хорошо резюмируется следующей ниже таблицей:
Тестирование гипотезы |
Тестовая статистика |
Правило отклонения гипотезы |
z-тесты |
z-статистика |
Если тестовая статистика ≥ z или ≤ -z, то отклонить нулевую гипотезу H0. |
t-тесты |
t-статистика |
Если тестовая статистика ≥ t или ≤ -t, то отклонить нулевую гипотезу H0. |
Анализ дисперсии (ANOVA) |
F-статистика |
Если тестовая статистика ≥ F, то отклонить нулевую гипотезу H0. |
Тесты хи-квадрат |
Статистика хи-квадрат |
Если тестовая статистика ≥ χ, то отклонить нулевую гипотезу H0. |
Для того, чтобы помочь сохранить поток посетителей веб-сайта, дизайнеры приступают к работе над вариантом веб-сайта с использованием всех новейших методов по поддержанию внимания аудитории. Мы хотели бы удостовериться, что наши усилия не напрасны, и поэтому стараемся увеличить время пребывания посетителей на обновленном веб-сайте.
Отсюда главный вопрос нашего исследования состоит в том, «приводит ли обновленный вид веб-сайта к увеличению времени пребывания на нем посетителей»? Мы принимаем решение проверить его относительно среднего значения времени пребывания. Теперь, мы должны изложить две наши гипотезы. По традиции считается, что изучаемые данные не содержат того, что исследователь ищет. Таким образом, консервативное мнение заключается в том, что данные не покажут ничего необычного. Все это называется нулевой гипотезой и обычно обозначается как H0.
При тестировании статистической гипотезы исходят из того, что нулевая гипотеза является истинной до тех пор, пока вес представленных данных, подтверждающих обратное, не сделает ее неправдоподобной. Этот подход к поиску доказательств «в обратную сторону» частично вытекает из простого психологического факта, что, когда люди пускаются на поиски чего-либо, они, как правило, это находят.
Затем исследователь формулирует альтернативную гипотезу, обозначаемую как H1. Она может попросту заключаться в том, что популяционное среднее отличается от базового уровня. Или же, что популяционное среднее больше или меньше базового уровня, либо больше или меньше на некоторую указанную величину. Мы хотели бы проверить, не увеличивает ли обновленный дизайн веб-сайта время пребывания, и поэтому нашей нулевой и альтернативной гипотезами будут следующие:
-
H0: Время пребывания для обновленного веб-сайта не отличается от времени пребывания для существующего веб-сайта
-
H1: Время пребывания для обновленного веб-сайта больше по сравнению с временем пребывания для существующего веб-сайта
Наше консервативное допущение состоит в том, что обновленный веб-сайт никак не влияет на время пребывания посетителей на веб-сайте. Нулевая гипотеза не обязательно должна быть «нулевой» (т.е. эффект отсутствует), но в данном случае, у нас нет никакого разумного оправдания, чтобы считать иначе. Если выборочные данные не поддержат нулевую гипотезу (т.е. если данные расходятся с ее допущением на слишком большую величину, чтобы носить случайный характер), то мы отклоним нулевую гипотезу и предложим альтернативную в качестве наилучшего альтернативного объяснения.
Указав нулевую и альтернативную гипотезы, мы должны установить уровень значимости, на котором мы ищем эффект.
Статистическая значимость
Проверка статистической значимости изначально разрабатывалась независимо от проверки статистических гипотез, однако сегодня оба подхода очень часто используются во взаимодействии друг с другом. Задача проверки статистической значимости состоит в том, чтобы установить порог, за пределами которого мы решаем, что наблюдаемые данные больше не поддерживают нулевую гипотезу.
Следовательно, существует два риска:
-
Мы можем принять расхождение как значимое, когда на самом деле оно возникло случайным образом
-
Мы можем приписать расхождение случайности, когда на самом деле оно показывает истинное расхождение с популяцией
Эти две возможности обозначаются соответственно, как ошибки 1-го и 2-го рода:
H0 ложная |
H0 истинная |
|
Отклонить H0 |
Истинноотрицательный исход |
Ошибка 1-го рода (ложноположительный исход) |
Принять H0 |
Ошибка 2-го рода (ложноотрицательный исход) |
Истинноположительный исход |
Чем больше мы уменьшаем риск совершения ошибок 1-го рода, тем больше мы увеличиваем риск совершения ошибок 2-го рода. Другими словами, с чем большей уверенностью мы хотим не заявлять о наличии расхождения, когда его нет, тем большее расхождение между выборками нам потребуется, чтобы заявить о статистической значимости. Эта ситуация увеличивает вероятность того, что мы проигнорируем подлинное расхождение, когда мы с ним столкнемся.
В статистической науке обычно используются два порога значимости. Это уровни в 5% и 1%. Расхождение в 5% обычно называют значимым, а расхождение в 1% — крайне значимым. В формулах этот порог часто обозначается греческой буквой α (альфа) и называется уровнем значимости. Поскольку, отсутствие эффекта по результатам эксперимента может рассматриваться как неуспех (эксперимента либо обновленного веб-сайта, как в нашем случае), то может возникнуть желание корректировать уровень значимости до тех пор, пока эффект не будет найден. По этой причине классический подход к проверке статистической значимости требует, чтобы мы устанавливали уровень значимости до того, как обратимся к нашим данным. Часто выбирается уровень в 5%, и поэтому мы на нем и остановимся.
Проверка обновленного дизайна веб-сайта
Веб-команда в AcmeContent была поглощена работой, конструируя обновленный веб-сайт, который будет стимулировать посетителей оставаться на нем в течение более длительного времени. Она употребила все новейшие методы и, в результате мы вполне уверены, что веб-сайт покажет заметное улучшение показателя времени пребывания.
Вместо того, чтобы запустить его для всех пользователей сразу, в AcmeContent хотели бы сначала проверить веб-сайт на небольшой выборке посетителей. Мы познакомили веб-команду с понятием искаженности выборки, и в результате там решили в течение одного дня перенаправлять случайные 5% трафика на обновленный веб-сайт. Результат с дневным трафиком был нам предоставлен одним текстовым файлом. Каждая строка показывает время пребывания посетителей. При этом, если посетитель пользовался исходным дизайном, ему присваивалось значение «0», и если он пользовался обновленным (и надеемся, улучшенным) дизайном, то ему присваивалось значение «1».
Выполнение z-теста
Ранее при тестировании с интервалами уверенности мы располагали лишь одним популяционным средним, с которым и выполнялось сравнение.
При тестировании нулевой гипотезы с помощью z-теста мы имеем возможность сравнивать две выборки. Посетители, которые видели обновленный веб-сайт, отбирались случайно, и данные для обеих групп были собраны в тот же день, чтобы исключить другие факторы с временной зависимостью.
Поскольку в нашем распоряжении имеется две выборки, то и стандартных ошибок у нас тоже две. Z-тест выполняется относительно объединенной стандартной ошибки, т.е. квадратного корня суммы дисперсий (вариансов), деленных на размеры выборок. Она будет такой же, что и результат, который мы получим, если взять стандартную ошибку обеих выборок вместе:
Здесь σ2a — это дисперсия выборки a, σ2b — дисперсия выборки b и соответственно na и nb — размеры выборок a и b. На Python объединенная стандартная ошибка вычисляется следующим образом:
def pooled_standard_error(a, b, unbias=False):
'''Объединенная стандартная ошибка'''
std1 = a.std(ddof=0) if unbias==False else a.std()
std2 = b.std(ddof=0) if unbias==False else b.std()
x = std1 ** 2 / a.count()
y = std2 ** 2 / b.count()
return sp.sqrt(x + y)
С целью выявления того, является ли видимое нами расхождение неожиданно большим, можно взять наблюдавшиеся расхождения между средними значениями на объединенной стандартной ошибке. Эту статистическую величину принято обозначать переменной z:
Используя функции pooled_standard_error
, которая вычисляет объединенную стандартную ошибку, z-статистику можно получить следующим образом:
def z_stat(a, b, unbias=False):
return (a.mean() - b.mean()) / pooled_standard_error(a, b, unbias)
Соотношение z объясняет, насколько средние значения отличаются относительно величины, которую мы ожидаем при заданной стандартной ошибке. Следовательно, z-статистика сообщает нам о том, на какое количество стандартных ошибок расходятся средние значения. Поскольку стандартная ошибка имеет нормальное распределение вероятностей, мы можем связать это расхождение с вероятностью, отыскав z-статистику в нормальной ИФР:
def z_test(a, b):
return stats.norm.cdf([ z_stat(a, b) ])
В следующем ниже примере z-тест используется для сравнения результативность двух веб-сайтов. Это делается путем группировки строк по номеру веб-сайта, в результате чего возвращается коллекция, в которой конкретному веб-сайту соответствует набор строк. Мы вызываем groupby('site')['dwell-time']
для конвертирования набора строк в набор значений времени пребывания. Затем вызываем функцию get_group
с номером группы, соответствующей номеру веб-сайта:
def ex_2_14():
'''Сравнение результативности двух вариантов
дизайна веб-сайта на основе z-теста'''
groups = load_data('new-site.tsv').groupby('site')['dwell-time']
a = groups.get_group(0)
b = groups.get_group(1)
print('a n: ', a.count())
print('b n: ', b.count())
print('z-статистика:', z_stat(a, b))
print('p-значение: ', z_test(a, b))
a n: 284
b n: 16
z-статистика: -1.6467438180091214
p-значение: [0.04980536]
Установление уровня значимости в размере 5% во многом аналогично установлению интервала уверенности шириной 95%. В сущности, мы надеемся убедиться, что наблюдавшееся расхождение попадает за пределы 95%-го интервала уверенности. Если это так, то мы можем утверждать, что нашли результат с 5%-ым уровнем значимости.
P-значение — это вероятность совершения ошибки 1-го рода в результате неправильного отклонения нулевой гипотезы, которая в действительности является истинной. Чем меньше p-значение, тем больше определенность в том, что нулевая гипотеза является ложной, и что мы нашли подлинный эффект.
Этот пример возвращает значение 0.0498, или 4.98%. Поскольку оно немногим меньше нашего 5% порога значимости, мы можем утверждать, что нашли нечто значимое.
Приведем еще раз нулевую и альтернативную гипотезы:
-
H0: Время пребывания на обновленном веб-сайте не отличается от времени пребывания на существующем веб-сайте
-
H1: Время пребывания на обновленном веб-сайте превышает время пребывания на существующем веб-сайте.
Наша альтернативная гипотеза состоит в том, что время пребывания на обновленном веб-сайте больше.
Мы готовы заявить о статистической значимости, и что время пребывания на обновленном веб-сайте больше по сравнению с существующим веб-сайтом, но тут есть одна трудность — чем меньше размер выборки, тем больше неопределенность в том, что выборочное стандартное отклонение совпадет с популяционным. Как показано в результатах предыдущего примера, наша выборка из обновленного веб-сайта содержит всего 16 посетителей. Столь малые выборки делают невалидным допущение о том, что стандартная ошибка нормально распределена.
К счастью, существует тест и связанное с ним распределение, которое моделирует увеличенную неопределенность стандартных ошибок для выборок меньших размеров.
t-распределение Студента
Популяризатором t-распределения был химик, работавший на пивоварню Гиннес в Ирландии, Уилльям Госсетт, который включил его в свой анализ темного пива Стаут.
В 1908 Уильям Госсет опубликовал статью об этой проверке в журнале Биометрика, но при этом по распоряжению своего работодателя, который рассматривал использованную Госсеттом статистику как коммерческую тайну, был вынужден использовать псевдоним. Госсет выбрал псевдоним «Студент».
В то время как нормальное распределение полностью описывается двумя параметрами — средним значением и стандартным отклонением, t-распределение описывается лишь одним параметром, так называемыми степенями свободы. Чем больше степеней свободы, тем больше t-распределение похоже на нормальное распределение с нулевым средним и стандартным отклонением, равным 1. По мере уменьшения степеней свободы, это распределение становится более широким с более толстыми чем у нормального распределения, хвостами.
Приведенный выше рисунок показывает, как t-распределение изменяется относительно нормального распределения при наличии разных степеней свободы. Более толстые хвосты для выборок меньших размеров соответствуют увеличенной возможности наблюдать более крупные отклонения от среднего значения.
Степени свободы
Степени свободы, часто обозначаемые сокращенно df от англ. degrees of freedom, тесно связаны с размером выборки. Это полезная статистика и интуитивно понятное свойство числового ряда, которое можно легко продемонстрировать на примере.
Если бы вам сказали, что среднее, состоящее из двух значений, равно 10 и что одно из значений равно 8, то Вам бы не потребовалась никакая дополнительная информация для того, чтобы суметь заключить, что другое значение равно 12. Другими словами, для размера выборки, равного двум, и заданного среднего значения одно из значений ограничивается, если другое известно.
Если напротив вам говорят, что среднее, состоящее из трех значений, равно 10, и первое значение тоже равно 10, то Вы были бы не в состоянии вывести оставшиеся два значения. Поскольку число множеств из трех чисел, начинающихся с 10, и чье среднее равно 10, является бесконечным, то прежде чем вы сможете вывести значение третьего, второе тоже должно быть указано.
Для любого множества из трех чисел ограничение простое: вы можете свободно выбрать первые два числа, но заключительное число ограничено. Степени свободы могут таким образом быть обобщены следующим образом: количество степеней свободы любой отдельной выборки на единицу меньше размера выборки.
При сопоставлении двух выборок степени свободы на две единицы меньше суммы размеров этих выборок, что равно сумме их индивидуальных степеней свободы.
t-статистика
При использовании t-распределения мы обращаемся к t-статистике. Как и z-статистика, эта величина количественно выражает степень маловероятности отдельно взятого наблюдавшегося отклонения. Для двухвыборочного t-теста соответствующая t-статистика вычисляется следующим образом:
Здесь Sa̅b̅ — это объединенная стандартная ошибка. Объединенная стандартная ошибка вычисляется таким же образом, как и раньше:
Однако это уравнение допускает наличие информации о популяционных параметрах σa и σb, которые можно аппроксимировать только на основе крупных выборок. t-тест предназначен для малых выборок и не требует от нас принимать допущения о поплуляционной дисперсии (вариансе).
Как следствие, объединенная стандартная ошибка для t-теста записывается как квадратный корень суммы стандартных ошибок:
На практике оба приведенных выше уравнения для объединенной стандартной ошибки дают идентичные результаты при заданных одинаковых входных последовательностях. Разница в математической записи всего лишь служит для иллюстрации того, что в условиях t-теста мы на входе зависим только от выборочных статистик. Объединенная стандартная ошибка может быть вычислена следующим образом:
def pooled_standard_error_t(a, b):
'''Объединенная стандартная ошибка для t-теста'''
return sp.sqrt(standard_error(a) ** 2 +
standard_error(b) ** 2)
Хотя в математическом плане t-статистика и z-статистика представлены по-разному, на практике процедура вычисления обоих идентичная:
t_stat = z_stat
def ex_2_15():
'''Вычисление t-статистики
двух вариантов дизайна веб-сайта'''
groups = load_data('new-site.tsv').groupby('site')['dwell-time']
a = groups.get_group(0)
b = groups.get_group(1)
return t_stat(a, b)
-1.6467438180091214
Различие между двумя выборочными показателями является не алгоритмическим, а концептуальным — z-статистика применима только тогда, когда выборки подчинены нормальному распределению.
t-тест
Разница в характере работы t-теста вытекает из распределения вероятностей, из которого вычисляется наше p-значение. Вычислив t-статистику, мы должны отыскать ее значение в t-распределении, параметризованном степенями свободы наших данных:
def t_test(a, b):
df = len(a) + len(b) - 2
return stats.t.sf([ abs(t_stat(a, b)) ], df)
Значение степени свободы обеих выборок на две единицы меньше их размеров, и для наших выборок составляет 298.
Напомним, что мы выполняем проверку статистической гипотезы. Поэтому выдвинем нашу нулевую и альтернативную гипотезы:
-
H0: Эта выборка взята из популяции с предоставленным средним значением
-
H1: Эта выборка взята из популяции со средним значением большего размера
Выполним следующий ниже пример:
def ex_2_16():
'''Сравнение результативности двух вариантов
дизайна веб-сайта на основе t-теста'''
groups = load_data('new-site.tsv').groupby('site')['dwell-time']
a = groups.get_group(0)
b = groups.get_group(1)
return t_test(a, b)
array([ 0.05033241])
Этот пример вернет p-значение, составляющее более 0.05. Поскольку оно больше α, равного 5%, который мы установили для проверки нулевой гипотезы, то мы не можем ее отклонить. Наша проверка с использованием t-теста значимого расхождения между средними значениями не обнаружила. Следовательно, наш едва значимый результат z-теста отчасти объясняется наличием слишком малой выборки.
Двухсторонние тесты
В нашей альтернативной гипотезе было принято неявное допущение, что обновленный веб-сайт будет работать лучше существующего. В процедуре проверки нулевой статистической гипотезы предпринимаются особые усилия для обеспечения того, чтобы при поиске статистической значимости мы не делали никаких скрытых допущений.
Проверки, при выполнении которых мы ищем только значимое количественное увеличение или уменьшение, называются односторонними и обычно не приветствуются, кроме случая, когда изменение в противоположном направлении было бы невозможным. Название термина «односторонний» обусловлено тем, что односторонняя проверка размещает всю α в одном хвосте распределения. Не делая проверок в другом направлении, проверка имеет больше мощности отклонить нулевую гипотезу в отдельно взятом направлении и, в сущности, понижает порог, по которому мы судим о результате как значимом.
Статистическая мощность — это вероятность правильного принятия альтернативной гипотезы. Она может рассматриваться как способность проверки обнаруживать эффект там, где имеется искомый эффект.
Хотя более высокая статистическая мощность выглядит желательной, она получается за счет наличия большей вероятности совершить ошибку 1-го рода. Правильнее было бы допустить возможность того, что обновленный веб-сайт может в действительности оказаться хуже существующего. Этот подход распределяет нашу α одинаково по обоим хвостам распределения и обеспечивает значимый результат, не искаженный под воздействием априорного допущения об улучшении работы обновленного веб-сайта.
В действительности в модуле stats библиотеки scipy уже предусмотрены функции для выполнения двухвыборочных t-проверок. Это функция stats.ttest_ind
. В качестве первого аргумента мы предоставляем выборку данных и в качестве второго — выборку для сопоставления. Если именованный аргумент equal_var
равен True
, то выполняется стандартная независимая проверка двух выборок, которая предполагает равные популяционные дисперсии, в противном случае выполняется проверка Уэлша (обратите внимание на служебную функцию t_test_verbose
, (которую можно найти среди примеров исходного кода в репо):
def ex_2_17():
'''Двухсторонний t-тест'''
groups = load_data('new-site.tsv').groupby('site')['dwell-time']
a = groups.get_group(0)
b = groups.get_group(1)
return t_test_verbose(a, sample2=b, fn=stats.ttest_ind) #t-тест Уэлша
{'p-значение': 0.12756432502462475,
'степени свободы ': 17.761382349686098,
'интервал уверенности': (76.00263198799597, 99.89877646270826),
'n1 ': 284,
'n2 ': 16,
'среднее x ': 87.95070422535211,
'среднее y ': 122.0,
'дисперсия x ': 10463.941024237296,
'дисперсия y ': 6669.866666666667,
't-статистика': -1.5985205593851322}
По результатам t-теста служебная функция t_test_verbose
возвращает много информации и в том числе p-значение. P-значение примерно в 2 раза больше того, которое мы вычислили для односторонней проверки. На деле, единственная причина, почему оно не совсем в два раза больше, состоит в том, что в модуле stats имплементирован легкий вариант t-теста, именуемый t-тестом Уэлша, который немного более робастен, когда две выборки имеют разные стандартные отклонения. Поскольку мы знаем, что для экспоненциальных распределений среднее значение и дисперсия тесно связаны, то этот тест немного более строг в применении и даже возвращает более низкую значимость.
Одновыборочный t-тест
Независимые выборки в рамках t-тестов являются наиболее распространенным видом статистического анализа, который обеспечивает очень гибкий и обобщенный способ установления, что две выборки представляют одинаковую либо разную популяцию. Однако в случаях, когда популяционное среднее уже известно, существует еще более простая проверка, представленная функцией библиотеки sciзy stats.ttest_1samp
.
Мы передаем выборку и популяционное среднее относительно которого выполняется проверка. Так, если мы просто хотим узнать, не отличается ли обновленный веб-сайт значимо от существующего популяционного среднего времени пребывания, равного 90 сек., то подобную проверку можно выполнить следующим образом:
def ex_2_18():
groups = load_data('new-site.tsv').groupby('site')['dwell-time']
b = groups.get_group(1)
return t_test_verbose(b, mean=90, fn=stats.ttest_1samp)
{'p-значение ': 0.13789520958229415,
'степени свободы df ': 15.0,
'интервал уверенности': (78.4815276659039, 165.5184723340961),
'n1 ': 16,
'среднее x ': 122.0,
'дисперсия x ': 6669.866666666667,
't-статистика ': 1.5672973291495713}
Служебная функция t_test_verbose
не только возвращает p-значение для выполненной проверки, но и интервал уверенности для популяционного среднего. Интервал имеет широкий диапазон между 78.5 и 165.5 сек., и, разумеется, перекрывается 90 сек. нашего теста. Как раз он и объясняет, почему мы не смогли отклонить нулевую гипотезу.
Многократные выборки
В целях развития интуитивного понимания относительно того, каким образом t-тест способен подтвердить и вычислить эти статистики из столь малых данных, мы можем применить подход, который связан с многократными выборками, от англ. resampling. Извлечение многократных выборок основывается на допущении о том, что каждая выборка является лишь одной из бесконечного числа возможных выборок из популяции. Мы можем лучше понять природу того, какими могли бы быть эти другие выборки, и, следовательно, добиться лучшего понимания опорной популяции, путем извлечения большого числа новых выборок из нашей существующей выборки.
На самом деле существует несколько методов взятия многократных выборок, и мы обсудим один из самых простых — бутстрапирование. При бустрапировании мы генерируем новую выборку, неоднократно извлекая из исходной выборки случайное значение с возвратом до тех пор, пока не сгенерируем выборку, имеющую тот же размер, что и оригинал. Поскольку выбранные значения возвращаются назад после каждого случайного отбора, то в новой выборке то же самое исходное значение может появляться многократно. Это как если бы мы неоднократно вынимали случайную карту из колоды игральных карт и каждый раз возвращали вынутую карту назад в колоду. В результате время от времени мы будем иметь карту, которую мы уже вынимали.
Бутстраповская выборка, или бутстрап, — синтетический набор данных, полученный в результате генерирования повторных выборок (с возвратом) из исследуемой выборки, используемой в качестве «суррогатной популяции», в целях аппроксимации выборочного распределения статистики (такой как, среднее, медиана и др.).
В библиотеке pandas при помощи функции sample можно легко извлекать бутстраповские выборки и генерировать большое число многократных выборок. Эта функция принимает ряд опциональных аргументов, в т.ч. n
(число элементов, которые нужно вернуть из числового ряда), axis
(ось, из которой извлекать выборку) и replace
(выборка с возвратом или без), по умолчанию равный False
. После этой функции можно задать метод агрегирования, вычисляющий сводную статистику в отношении бутстраповских выборок:
def ex_2_19():
'''Построение графика синтетических времен пребывания
путем извлечения бутстраповских выборок'''
groups = load_data('new-site.tsv').groupby('site')['dwell-time']
b = groups.get_group(1)
xs = [b.sample(len(b), replace=True).mean() for _ in range(1000)]
pd.Series(xs).hist(bins=20)
plt.xlabel('Бутстрапированные средние значения времени пребывания, сек.')
plt.ylabel('Частота')
plt.show()
Приведенный выше пример наглядно показывает результаты на гистограмме:
Гистограмма демонстрирует то, как средние значения изменялись вместе с многократными выборками, взятыми из времени пребывания на обновленном веб-сайте. Хотя на входе имелась лишь одна выборка, состоящая из 16 посетителей, бутстрапированные выборки очень четко просимулировали стандартную ошибку изначальной выборки и позволили визуализировать интервал уверенности (между 78 и 165 сек.), вычисленный ранее в результате одновыборочного t-теста.
Благодаря бутстрапированию мы просимулировали взятие многократных выборок, при том, что у нас на входе имелась всего одна выборка. Этот метод обычно применяется для оценивания параметров, которые мы не способны или не знаем, как вычислить аналитически.
Проверка многочисленных вариантов дизайна
Было разочарованием обнаружить отсутствие статистической значимости на фоне увеличенного времени пребывания пользователей на обновленном веб-сайте. Хотя хорошо, что мы обнаружили это на малой выборке пользователей, прежде чем выкладывать его на всеобщее обозрение.
Не позволяя себя обескуражить, веб-команда AcmeContent берется за сверхурочную работу и создает комплект альтернативных вариантов дизайна веб-сайта. Беря лучшие элементы из других проектов, они разрабатывают 19 вариантов для проверки. Вместе с нашим изначальным веб-сайтом, который будет действовать в качестве контрольного, всего имеется 20 разных вариантов дизайна веб-сайта, куда посетители будут перенаправляться.
Вычисление выборочных средних
Веб-команда разворачивает 19 вариантов дизайна обновленного веб-сайта наряду с изначальным. Как отмечалось ранее, каждый вариант дизайна получает случайные 5% посетителей, и при этом наше испытание проводится в течение 24 часов.
На следующий день мы получаем файл, показывающий значения времени пребывания посетителей на каждом варианте веб-сайта. Все они были промаркированы числами, при этом число 0 соответствовало веб-сайту с исходным дизайном, а числа от 1 до 19 представляли другие варианты дизайна:
def ex_2_20():
df = load_data('multiple-sites.tsv')
return df.groupby('site').aggregate(sp.mean)
Этот пример сгенерирует следующую ниже таблицу:
site |
dwell-time |
0 |
79.851064 |
1 |
106.000000 |
2 |
88.229167 |
3 |
97.479167 |
4 |
94.333333 |
5 |
102.333333 |
6 |
144.192982 |
7 |
123.367347 |
8 |
94.346939 |
9 |
89.820000 |
10 |
129.952381 |
11 |
96.982143 |
12 |
80.950820 |
13 |
90.737705 |
14 |
74.764706 |
15 |
119.347826 |
16 |
86.744186 |
17 |
77.891304 |
18 |
94.814815 |
19 |
89.280702 |
Мы хотели бы проверить каждый вариант дизайна веб-сайта, чтобы увидеть, не генерирует ли какой-либо из них статистически значимый результат. Для этого можно сравнить варианты дизайна веб-сайта друг с другом следующим образом, причем нам потребуется вспомогательный модуль Python itertools, который содержит набор функций, создающих итераторы для эффективной циклической обработки:
import itertools
def ex_2_21():
'''Проверка вариантов дизайна веб-сайта на основе t-теста
по принципу "каждый с каждым"'''
groups = load_data('multiple-sites.tsv').groupby('site')
alpha = 0.05
pairs = [list(x) # найти сочетания из n по k
for x in itertools.combinations(range(len(groups)), 2)]
for pair in pairs:
gr, gr2 = groups.get_group( pair[0] ), groups.get_group( pair[1] )
site_a, site_b = pair[0], pair[1]
a, b = gr['dwell-time'], gr2['dwell-time']
p_val = stats.ttest_ind(a, b, equal_var = False).pvalue
if p_val < alpha:
print('Варианты веб-сайта %i и %i значимо различаются: %f'
% (site_a, site_b, p_val))
Однако это было бы неправильно. Мы скорее всего увидим статистическое расхождение между вариантами дизайна, показавшими себя в особенности хорошо по сравнению с вариантами, показавшими себя в особенности плохо, даже если эти расхождения носили случайный характер. Если вы выполните приведенный выше пример, то увидите, что многие варианты дизайна веб-сайта статистически друг от друга отличаются.
С другой стороны, мы можем сравнить каждый вариант дизайна веб-сайта с нашим текущим изначальным значением — средним значением времени пребывания, равным 90 сек., измеренным на данный момент для существующего веб-сайта:
def ex_2_22():
groups = load_data('multiple-sites.tsv').groupby('site')
alpha = 0.05
baseline = groups.get_group(0)['dwell-time']
for site_a in range(1, len(groups)):
a = groups.get_group( site_a )['dwell-time']
p_val = stats.ttest_ind(a, baseline, equal_var = False).pvalue
if p_val < alpha:
print('Вариант %i веб-сайта значимо отличается: %f'
% (site_a, p_val))
В результате этой проверки будут идентифицированы два варианта дизайна веб-сайта, которые существенно отличаются:
Вариант 6 веб-сайта значимо отличается: 0.005534
Вариант 10 веб-сайта 10 значимо отличается: 0.006881
Малые p-значения (меньше 1%) указывают на то, что существует статистически очень значимые расхождения. Этот результат представляется весьма многообещающим, однако тут есть одна проблема. Мы выполнили t-тест по 20 выборкам данных с уровнем значимости α, равным 0.05. Уровень значимости α определяется, как вероятность неправильного отказа от нулевой гипотезы. На самом деле после 20-кратного выполнения t-теста становится вероятным, что мы неправильно отклоним нулевую гипотезу по крайней мере для одного варианта веб-сайта из 20.
Сравнивая таким одновременным образом многочисленные страницы, мы делаем результаты t-теста невалидными. Существует целый ряд альтернативных технических приемов решения проблемы выполнения многократных сравнений в статистических тестах. Эти методы будут рассмотрены в следующем разделе.
Поправка Бонферрони
Для проведения многократных проверок используется подход, который объясняет увеличенную вероятность обнаружить значимый эффект в силу многократных испытаний. Поправка Бонферрони — это очень простая корректировка, которая обеспечивает, чтобы мы вряд ли совершили ошибки 1-го рода. Она выполняется путем настройки значения уровня значимости для тестов.
Настройка очень простая — поправка Бонферрони попросту делит требуемое значение α на число тестов. Например, если для теста имелось k вариантов дизайна веб-сайта, и α эксперимента равно 0.05, то поправка Бонферрони выражается следующим образом:
Она представляет собой безопасный способ смягчить увеличение вероятности совершения ошибки 1-го рода при многократной проверке. Следующий пример идентичен примеру ex-2-22
, за исключением того, что значение α разделено на число групп:
def ex_2_23():
'''Проверка вариантов дизайна веб-сайта на основе t-теста
против исходного (0) с поправкой Бонферрони'''
groups = load_data('multiple-sites.tsv').groupby('site')
alpha = 0.05 / len(groups)
baseline = groups.get_group(0)['dwell-time']
for site_a in range(1, len(groups)):
a = groups.get_group(site_a)['dwell-time']
p_val = stats.ttest_ind(a, baseline, equal_var = False).pvalue
if p_val < alpha:
print('Вариант %i веб-сайта значимо отличается от исходного: %f'
% (site_a, p_val))
Если вы выполните приведенный выше пример, то увидите, что при использовании поправки Бонферрони ни один из веб-сайтов больше не считается статистически значимым.
Метод проверки статистической значимости связан с поддержанием равновесия — чем меньше шансы совершения ошибки 1-го рода, тем больше риск совершения ошибки 2-го рода. Поправка Бонферрони очень консервативна, и весьма возможно, что из-за излишней осторожности мы пропускаем подлинное расхождение.
Примеры исходного кода для этого поста находятся в моем репо на Github. Все исходные данные взяты в репозитории автора книги.
В заключительном посте, посте №4, этой серии постов мы проведем исследование альтернативного подхода к проверке статистической значимости, который позволяет устанавливать равновесие между совершением ошибок 1-го и 2-го рода, давая нам возможность проверить все 20 вариантов веб-сайта одновременно.
Метрики качества классификаторов
Матрица ошибок (Confusion matrix)
Матрица ошибок — это способ разбить объекты на четыре категории в зависимости от комбинации истинного ответа и ответа алгоритма.
Основные термины:
- TPTP — истино-положительное решение;
- TNTN — истино-отрицательное решение;
- FPFP — ложно-положительное решение (Ошибка первого рода);
- FNFN — ложно-отрицательное решение (Ошибка второго рода).
Интерактивная картинка с большим числом метрик
Accuracy — доля правильных ответов:
Accuracy=TP+TNP+N=TP+TNTP+TN+FP+FN{\displaystyle \mathrm {Accuracy} ={\frac {\mathrm {TP} +\mathrm {TN} }{P+N}}={\frac {\mathrm {TP} +\mathrm {TN} }{\mathrm {TP} +\mathrm {TN} +\mathrm {FP} +\mathrm {FN} }}}
Данная матрика имеет существенный недостаток — её значение необходимо оценивать в контексте баланса классов. Eсли в выборке 950 отрицательных и 50 положительных объектов, то при абсолютно случайной классификации мы получим долю правильных ответов 0.95. Это означает, что доля положительных ответов сама по себе не несет никакой информации о качестве работы алгоритма a(x), и вместе с ней следует анализировать соотношение классов в выборке.
Гораздо более информативными критериями являются точность (precision) и полнота (recall).
Точность показывает, какая доля объектов, выделенных классификатором как положительные, действительно является положительными:
Precision=TPTP+FPPrecision = \frac{TP}{TP+FP}
Полнота показывает, какая часть положительных объектов была выделена классификатором:
Recall=TPTP+FNRecall = \frac{TP}{TP+FN}
Существует несколько способов получить один критерий качества на основе точности и полноты. Один из них — F-мера, гармоническое среднее точности и полноты:
F_\beta = (1 + \beta^2) \cdot \frac{\mathrm{precision} \cdot \mathrm{recall}}{(\beta^2 \cdot \mathrm{precision}) + \mathrm{recall}} = \frac {(1 + \beta^2) \cdot \mathrm{true\ positive} }{(1 + \beta^2) \cdot \mathrm{true\ positive} + \beta^2 \cdot \mathrm{false\ negative} + \mathrm{false\ positive}}
Среднее гармоническое обладает важным свойством — оно близко к нулю, если хотя бы один из аргументов близок к нулю. Именно поэтому оно является более предпочтительным, чем среднее арифметическое (если алгоритм будет относить все объекты к положительному классу, то он будет иметь recall = 1 и precision больше 0, а их среднее арифметическое будет больше 1/2, что недопустимо).
Чаще всего берут β=1\beta=1, хотя иногда встречаются и другие модификации. F2F_2 острее реагирует на recall (т. е. на долю ложноположительных ответов), а F0.5F_{0.5} чувствительнее к точности (ослабляет влияние ложноположительных ответов).
В sklearn есть удобная функция sklearn.metrics.classification_report, возвращающая recall, precision и F-меру для каждого из классов, а также количество экземпляров каждого класса.
import pandas as pd import seaborn as sns from sklearn.model_selection import train_test_split from sklearn.preprocessing import StandardScaler, LabelEncoder from sklearn.linear_model import LogisticRegression from sklearn.svm import SVC from sklearn import datasets import numpy as np import matplotlib.pyplot as plt import seaborn as sns %matplotlib inline %config InlineBackend.figure_format = 'svg'
from sklearn.metrics import classification_report y_true = [0, 1, 2, 2, 2] y_pred = [0, 0, 2, 2, 1] target_names = ['class 0', 'class 1', 'class 2'] print(classification_report(y_true, y_pred, target_names=target_names))
precision recall f1-score support
class 0 0.50 1.00 0.67 1
class 1 0.00 0.00 0.00 1
class 2 1.00 0.67 0.80 3
micro avg 0.60 0.60 0.60 5
macro avg 0.50 0.56 0.49 5
weighted avg 0.70 0.60 0.61 5
Линейная классификация
Основная идея линейного классификатора заключается в том, что признаковое пространство может быть разделено гиперплоскостью на две полуплоскости, в каждой из которых прогнозируется одно из двух значений целевого класса.
Если это можно сделать без ошибок, то обучающая выборка называется линейно разделимой.
Указанная разделяющая плоскость называется линейным дискриминантом.
Логистическая регрессия
Логистическая регрессия является частным случаем линейного классификатора, но она обладает хорошим «умением» – прогнозировать вероятность отнесения наблюдения к классу. Таким образом, результат логистической регрессии всегда находится в интервале [0, 1].
iris = pd.read_csv("https://nagornyy.me/datasets/iris.csv")
sepal_length | sepal_width | petal_length | petal_width | |
---|---|---|---|---|
count | 150.000000 | 150.000000 | 150.000000 | 150.000000 |
mean | 5.843333 | 3.057333 | 3.758000 | 1.199333 |
std | 0.828066 | 0.435866 | 1.765298 | 0.762238 |
min | 4.300000 | 2.000000 | 1.000000 | 0.100000 |
25% | 5.100000 | 2.800000 | 1.600000 | 0.300000 |
50% | 5.800000 | 3.000000 | 4.350000 | 1.300000 |
75% | 6.400000 | 3.300000 | 5.100000 | 1.800000 |
max | 7.900000 | 4.400000 | 6.900000 | 2.500000 |
sns.pairplot(iris, hue="species")
/usr/local/lib/python3.7/site-packages/scipy/stats/stats.py:1713: FutureWarning: Using a non-tuple sequence for multidimensional indexing is deprecated; use `arr[tuple(seq)]` instead of `arr[seq]`. In the future this will be interpreted as an array index, `arr[np.array(seq)]`, which will result either in an error or a different result.
return np.add.reduce(sorted[indexer] * weights, axis=axis) / sumval
sns.lmplot(x="petal_length", y="petal_width", data=iris)
X = iris.iloc[:, 2:4].values y = iris['species'].values
array([‘setosa’, ‘setosa’, ‘setosa’, ‘setosa’, ‘setosa’], dtype=object)
from sklearn.preprocessing import LabelEncoder le = LabelEncoder() le.fit(y) y = le.transform(y) y[:5]
iris_pred_names = le.classes_ iris_pred_names
array([‘setosa’, ‘versicolor’, ‘virginica’], dtype=object)
from sklearn.model_selection import train_test_split X_train, X_test, y_train, y_test = train_test_split( X, y, test_size=0.3, random_state=0)
from sklearn.preprocessing import StandardScaler sc = StandardScaler() sc.fit(X_train) X_train_std = sc.transform(X_train) X_test_std = sc.transform(X_test)
X_train[:5], X_train_std[:5]
(array([[3.5, 1. ],
[5.5, 1.8],
[5.7, 2.5],
[5. , 1.5],
[5.8, 1.8]]), array([[-0.18295039, -0.29318114],
[ 0.93066067, 0.7372463 ],
[ 1.04202177, 1.63887031],
[ 0.6522579 , 0.35083601],
[ 1.09770233, 0.7372463 ]]))
from sklearn.linear_model import LogisticRegression lr = LogisticRegression(C=100.0, random_state=1) lr.fit(X_train_std, y_train)
/usr/local/lib/python3.7/site-packages/sklearn/linear_model/logistic.py:433: FutureWarning: Default solver will be changed to ‘lbfgs’ in 0.22. Specify a solver to silence this warning.
FutureWarning)
/usr/local/lib/python3.7/site-packages/sklearn/linear_model/logistic.py:460: FutureWarning: Default multi_class will be changed to ‘auto’ in 0.22. Specify the multi_class option to silence this warning.
«this warning.», FutureWarning)
LogisticRegression(C=100.0, class_weight=None, dual=False, fit_intercept=True,
intercept_scaling=1, max_iter=100, multi_class=’warn’,
n_jobs=None, penalty=’l2′, random_state=1, solver=’warn’,
tol=0.0001, verbose=0, warm_start=False)
lr.predict_proba(X_test_std[:3, :])
array([[2.77475804e-08, 6.31730607e-02, 9.36826912e-01],
[7.87476628e-03, 9.91707489e-01, 4.17744834e-04],
[8.15542033e-01, 1.84457967e-01, 8.14812482e-12]])
lr.predict_proba(X_test_std[:3, :]).sum(axis=1)
lr.predict_proba(X_test_std[:3, :]).argmax(axis=1)
Предсказываем класс первого наблюдения
lr.predict(X_test_std[0, :].reshape(1, -1))
На основе его коэффициентов:
array([0.70793846, 1.51006688])
X_test_std[0, :].reshape(1, -1)
array([[0.70793846, 1.51006688]])
y_pred = lr.predict(X_test_std)
print(classification_report(y_test, y_pred, target_names=iris_pred_names))
precision recall f1-score support
setosa 1.00 1.00 1.00 16
versicolor 1.00 0.94 0.97 18
virginica 0.92 1.00 0.96 11
micro avg 0.98 0.98 0.98 45
macro avg 0.97 0.98 0.98 45
weighted avg 0.98 0.98 0.98 45
Машина опорных векторов
Основная идея метода — перевод исходных векторов в пространство более высокой размерности и поиск разделяющей гиперплоскости с максимальным зазором в этом пространстве. Две параллельных гиперплоскости строятся по обеим сторонам гиперплоскости, разделяющей классы. Разделяющей гиперплоскостью будет гиперплоскость, максимизирующая расстояние до двух параллельных гиперплоскостей. Алгоритм работает в предположении, что чем больше разница или расстояние между этими параллельными гиперплоскостями, тем меньше будет средняя ошибка классификатора.
На практике случаи, когда данные можно разделить гиперплоскостью, довольно редки. В этом случае поступают так: все элементы обучающей выборки вкладываются в пространство X более высокой размерности, так, чтобы выборка была линейно разделима.
from sklearn.svm import SVC svm = SVC(kernel='linear', C=1.0, random_state=1) svm.fit(X_train_std, y_train)
SVC(C=1.0, cache_size=200, class_weight=None, coef0=0.0,
decision_function_shape=’ovr’, degree=3, gamma=’auto_deprecated’,
kernel=’linear’, max_iter=-1, probability=False, random_state=1,
shrinking=True, tol=0.001, verbose=False)
Kernel (ядро) отвечается за гиперплоскость и может принимать значения linear (для линейной), rbf (для нелинейной) и другие.
С — параметр регуляризации. Он в том числе контролирует соотношение между гладкой границей и корректной классификацией рассматриваемых точек.
gamma — это «ширина» rbf ядра (kernel). Она участвует в подгонке модели и может являться причиной переобучения.
y_pred_svm = svm.predict(X_test_std)
print(classification_report(y_test, y_pred_svm, target_names=iris_pred_names))
precision recall f1-score support
setosa 1.00 1.00 1.00 16
versicolor 1.00 0.94 0.97 18
virginica 0.92 1.00 0.96 11
micro avg 0.98 0.98 0.98 45
macro avg 0.97 0.98 0.98 45
weighted avg 0.98 0.98 0.98 45
Нелинейная классификация
svm_rbf = SVC(kernel='rbf', random_state=1, gamma=0.10, C=10.0) svm_rbf.fit(X_train_std, y_train)
SVC(C=10.0, cache_size=200, class_weight=None, coef0=0.0,
decision_function_shape=’ovr’, degree=3, gamma=0.1, kernel=’rbf’,
max_iter=-1, probability=False, random_state=1, shrinking=True,
tol=0.001, verbose=False)
y_pred_svm_rbf = svm_rbf.predict(X_test_std) print(classification_report(y_test, y_pred_svm_rbf, target_names=iris_pred_names))
precision recall f1-score support
setosa 1.00 1.00 1.00 16
versicolor 1.00 0.94 0.97 18
virginica 0.92 1.00 0.96 11
micro avg 0.98 0.98 0.98 45
macro avg 0.97 0.98 0.98 45
weighted avg 0.98 0.98 0.98 45
Деревья решений
Деревья решений используются в повседневной жизни в самых разных областях человеческой деятельности.
До внедрения масштабируемых алгоритмов машинного обучения в банковской сфере задача кредитного скоринга решалась экспертами. Решение о выдаче кредита заемщику принималось на основе некоторых интуитивно (или по опыту) выведенных правил, которые можно представить в виде дерева решений:
В этом случае можно сказать, что решается задача бинарной классификации (целевой класс имеет два значения: «Выдать кредит» и «Отказать») по признакам «Возраст», «Наличие дома», «Доход» и «Образование».
Дерево решений как алгоритм машинного обучения – по сути то же самое. Огромное преимущество деревьев решений в том, что они легко интерпретируемы, понятны человеку.
В основе популярных алгоритмов построения дерева решений лежит принцип жадной максимизации прироста информации – на каждом шаге выбирается тот признак, при разделении по которому прирост информации оказывается наибольшим. Дальше процедура повторяется рекурсивно, пока энтропия не окажется равной нулю или какой-то малой величине (если дерево не подгоняется идеально под обучающую выборку во избежание переобучения).
Плюсы:
- Порождение четких правил классификации, понятных человеку, например, «если возраст < 25 и интерес к мотоциклам, то отказать в кредите». Это свойство называют интерпретируемостью модели;
- Деревья решений могут легко визуализироваться, как сама модель (дерево), так и прогноз для отдельного взятого тестового объекта (путь в дереве);
- Быстрые процессы обучения и прогнозирования;
- Малое число параметров модели;
- Поддержка и числовых, и категориальных признаков.
Минусы:
- У порождения четких правил классификации есть и другая сторона: деревья очень чувствительны к шумам во входных данных, вся модель может кардинально измениться, если немного изменится обучающая выборка (например, если убрать один из признаков или добавить несколько объектов), поэтому и правила классификации могут сильно изменяться, что ухудшает интерпретируемость модели;
- Разделяющая граница, построенная деревом решений, имеет свои ограничения (состоит из гиперплоскостей, перпендикулярных какой-то из координатной оси), и на практике дерево решений по качеству классификации уступает некоторым другим методам;
from sklearn.tree import DecisionTreeClassifier tree = DecisionTreeClassifier(criterion='gini', max_depth=4, random_state=1) tree.fit(X_train_std, y_train)
DecisionTreeClassifier(class_weight=None, criterion=’gini’, max_depth=4,
max_features=None, max_leaf_nodes=None,
min_impurity_decrease=0.0, min_impurity_split=None,
min_samples_leaf=1, min_samples_split=2,
min_weight_fraction_leaf=0.0, presort=False, random_state=1,
splitter=’best’)
y_pred_tree = tree.predict(X_test_std) print(classification_report(y_test, y_pred_tree, target_names=iris_pred_names))
precision recall f1-score support
setosa 1.00 1.00 1.00 16
versicolor 1.00 0.94 0.97 18
virginica 0.92 1.00 0.96 11
micro avg 0.98 0.98 0.98 45
macro avg 0.97 0.98 0.98 45
weighted avg 0.98 0.98 0.98 45
Collecting pydotplus
[?25l Downloading https://files.pythonhosted.org/packages/60/bf/62567830b700d9f6930e9ab6831d6ba256f7b0b730acb37278b0ccdffacf/pydotplus-2.0.2.tar.gz (278kB)
[K 100% |████████████████████████████████| 286kB 1.5MB/s
[?25hRequirement already satisfied: pyparsing>=2.0.1 in /usr/local/lib/python3.7/site-packages (from pydotplus) (2.3.0)
Building wheels for collected packages: pydotplus
Running setup.py bdist_wheel for pydotplus … [?25ldone
[?25h Stored in directory: /Users/hun/Library/Caches/pip/wheels/35/7b/ab/66fb7b2ac1f6df87475b09dc48e707b6e0de80a6d8444e3628
Successfully built pydotplus
Installing collected packages: pydotplus
Successfully installed pydotplus-2.0.2
from pydotplus import graph_from_dot_data from sklearn.tree import export_graphviz dot_data = export_graphviz(tree, filled=True, rounded=True, class_names=iris_pred_names, feature_names=['petal length', 'petal width'], out_file=None) graph = graph_from_dot_data(dot_data)
Метод ближайших соседей
Метод ближайших соседей (k Nearest Neighbors, или kNN) — тоже очень популярный метод классификации, также иногда используемый в задачах регрессии. Это, наравне с деревом решений, один из самых понятных подходов к классификации. На уровне интуиции суть метода такова: посмотри на соседей, какие преобладают, таков и ты. Формально основой метода является гипотеза компактности: если метрика расстояния между примерами введена достаточно удачно, то схожие примеры гораздо чаще лежат в одном классе, чем в разных.
Для классификации каждого из объектов тестовой выборки необходимо последовательно выполнить следующие операции:
- Вычислить расстояние до каждого из объектов обучающей выборки
- Отобрать k объектов обучающей выборки, расстояние до которых минимально
- Класс классифицируемого объекта — это класс, наиболее часто встречающийся среди k ближайших соседей
Под задачу регрессии метод адаптируется довольно легко – на 3 шаге возвращается не метка, а число – среднее (или медианное) значение целевого признака среди соседей.
Примечательное свойство такого подхода – его ленивость. Это значит, что вычисления начинаются только в момент классификации тестового примера, а заранее, только при наличии обучающих примеров, никакая модель не строится. В этом отличие, например, от ранее рассмотренного дерева решений, где сначала на основе обучающей выборки строится дерево, а потом относительно быстро происходит классификация тестовых примеров.
Качество классификации/регрессии методом ближайших соседей зависит от нескольких параметров:
- число соседей
- метрика расстояния между объектами (часто используются метрика Хэмминга, евклидово расстояние, косинусное расстояние и расстояние Минковского). Отметим, что при использовании большинства метрик значения признаков надо масштабировать. Условно говоря, чтобы признак «Зарплата» с диапазоном значений до 100 тысяч не вносил больший вклад в расстояние, чем «Возраст» со значениями до 100.
- веса соседей (соседи тестового примера могут входить с разными весами, например, чем дальше пример, тем с меньшим коэффициентом учитывается его «голос»)
Плюсы и минусы метода ближайших соседей
Плюсы:
- Простая реализация;
- Неплохо изучен теоретически;
- Как правило, метод хорош для первого решения задачи;
- Можно адаптировать под нужную задачу выбором метрики или ядра (ядро может задавать операцию сходства для сложных объектов типа графов, а сам подход kNN остается тем же);
- Неплохая интерпретация, можно объяснить, почему тестовый пример был классифицирован именно так.
Минусы:
- Метод считается быстрым в сравнении, например, с композициями алгоритмов, но в реальных задачах, как правило, число соседей, используемых для классификации, будет большим (100-150), и в таком случае алгоритм будет работать не так быстро, как дерево решений;
- Если в наборе данных много признаков, то трудно подобрать подходящие веса и определить, какие признаки не важны для классификации/регрессии;
- Зависимость от выбранной метрики расстояния между примерами. Выбор по умолчанию евклидового расстояния чаще всего ничем не обоснован. Можно отыскать хорошее решение перебором параметров, но для большого набора данных это отнимает много времени;
- Нет теоретических оснований выбора определенного числа соседей — только перебор (впрочем, чаще всего это верно для всех гиперпараметров всех моделей). В случае малого числа соседей метод чувствителен к выбросам, то есть склонен переобучаться;
- Как правило, плохо работает, когда признаков много, из-за «прояклятия размерности». Про это хорошо рассказывает известный в ML-сообществе профессор Pedro Domingos – тут в популярной статье «A Few Useful Things to Know about Machine Learning», также «the curse of dimensionality» описывается в книге Deep Learning в главе «Machine Learning basics».
Класс KNeighborsClassifier в Scikit-learn
Основные параметры класса sklearn.neighbors.KNeighborsClassifier:
- weights: «uniform» (все веса равны), «distance» (вес обратно пропорционален расстоянию до тестового примера) или другая определенная пользователем функция
- algorithm (опционально): «brute», «ball_tree», «KD_tree», или «auto». В первом случае ближайшие соседи для каждого тестового примера считаются перебором обучающей выборки. Во втором и третьем — расстояние между примерами хранятся в дереве, что ускоряет нахождение ближайших соседей. В случае указания параметра «auto» подходящий способ нахождения соседей будет выбран автоматически на основе обучающей выборки.
- leaf_size (опционально): порог переключения на полный перебор в случае выбора BallTree или KDTree для нахождения соседей
- metric: «minkowski», «manhattan», «euclidean», «chebyshev» и другие
from sklearn.neighbors import KNeighborsClassifier knn = KNeighborsClassifier(n_neighbors=5, p=2, metric='minkowski') knn.fit(X_train_std, y_train) y_pred_knn = knn.predict(X_test_std) print(classification_report(y_test, y_pred_knn, target_names=iris_pred_names))
precision recall f1-score support
setosa 1.00 1.00 1.00 16
versicolor 1.00 1.00 1.00 18
virginica 1.00 1.00 1.00 11
micro avg 1.00 1.00 1.00 45
macro avg 1.00 1.00 1.00 45
weighted avg 1.00 1.00 1.00 45
Домашнаяя работа
Примените изученные классификаторы для предсказания выживаемости на Титанике и постойте наилучший классификатор. Каковы значения основных его метрик?
Опирайтесь на эти статьи:
- Kaggle и Titanic — еще одно решение задачи с помощью Python
- Основы анализа данных на python с использованием pandas+sklearn
- Титаник на Kaggle: вы не дочитаете этот пост до конца
Рассмотрим основные метрики качества в задачах классификации : доля правильных ответов, точность, полнота, F-мера и матрица ошибок. А также четыре различных комбинации фактических и прогнозируемых значений: истинно отрицательные (TN), ложноотрицательные (FN), истинно положительные (TP) и ложноположительные (FP).
Напомним, что в случае бинарной классификации мы говорим о положительном (positive) классе и отрицательном (negative) классе, подразумевая под положительным классом интересующий нас класс.
True Positive (TP) — это общее количество предсказанных положительных результатов, которые на самом деле являются положительными.
True Negative (TN) — это общее количество предсказанных отрицательных результатов, которые на самом деле являются отрицательными.
Ложное срабатывание (FP), ошибка 1 рода, — это общее количество предсказанных положительных результатов, которые на самом деле являются отрицательными.
Ложноотрицательный (FN), ошибка 2 рода — это общее количество предсказанных отрицательных результатов, которые на самом деле являются положительными.
Матрица ошибок (Сonfusion matrix)
#### Доля правильных ответов (Accuracy Score)
$$Accuracy Score = \frac{Correct Predictions}{All Predictions} = \frac{TP + TN}{TP + TN + FP + FN}$$
#### Точность (Precision)
$$Positive Precision = \frac{True Positive}{Predicted Positive} = \frac{TP}{TP + FP}$$
$$Negative Precision = \frac{True Negative}{Predicted Negative} = \frac{TN}{TN + FN}$$
#### Полнота (Recall)
$$Positive Recall = \frac{True Positive}{Actual Positive} = \frac{TP}{TP + FN}$$
$$Negative Recall = \frac{True Negative}{Actual Negative} = \frac{TN}{TN + FP}$$
### F-мера (f1-score)
$$\frac{2*Precision*Recall}{Precision+Recall}$$
В качестве примера рассмотрим следующую задачу : в период с января по март происходит большая часть увольнений торгового персонала, необходимо создать модель, которая по ряду показателей сотрудника определяла в нем потенциального кандидата на увольнение. В этом примере положительным классом будем считать увольнение.
Загружаем наши данные и смотрим первые и последние пять строк :
df = pd.read_excel(«data\shop_quit.xlsx»)
pd.concat([df.head(), df.tail()])
sex | experience | wd | salary | quit | |
---|---|---|---|---|---|
0 | female | 184.5 | 39.0 | 3738 | 0 |
1 | female | 34.5 | 47.0 | 3079 | 0 |
2 | female | 79.4 | 42.5 | 2887 | 0 |
3 | male | 4.7 | 47.5 | 2742 | 0 |
4 | female | 56.3 | 64.0 | 3483 | 0 |
116 | male | 28.1 | 26.0 | 2651 | 0 |
117 | female | 19.0 | 63.0 | 3177 | 0 |
118 | female | 18.9 | 49.0 | 2239 | 0 |
119 | male | 46.6 | 4.0 | 4731 | 0 |
120 | female | 21.5 | 17.5 | 2163 | 0 |
Расшифруем показатели :
- sex — пол
- experience — стаж работы в месяцах
- wd — количество отработанных смен
- salary — среднесменная заработная плата
- quit — метка 0-работает, 1 — уволился
Перекодируем пол
from sklearn.preprocessing import LabelEncoder
class_le = LabelEncoder ()
df.sex=class_le.fit_transform(df.sex.values)
df.sample(5)
sex | experience | wd | salary | quit | |
---|---|---|---|---|---|
98 | 1 | 93.8 | 40.5 | 6092 | 0 |
38 | 0 | 0.9 | 19.0 | 3734 | 0 |
34 | 0 | 20.6 | 23.0 | 3526 | 0 |
104 | 0 | 13.2 | 44.0 | 3687 | 0 |
115 | 0 | 18.5 | 33.5 | 3494 | 0 |
Присвоим показатели матрице data , а метки увольнения — вектору label
label_col = «quit»
data = df.drop(label_col, axis=1)
label = df[label_col]
Разделим данные на обучающую и тестовую части
d_train, d_test, l_train, l_test = train_test_split(data, label, test_size=0.4, random_state=23, stratify=label)
Через stratify = label используем встроенную поддержку стратификации. Это означает,
что метод train_test_split возвращает обучающий и испытательный
поднаборы, имеющие такие же количественные соотношения меток классов ,
как у входного набора данных. Мы можем проверить , так ли это , с применением функции bincount из NumPy, которая подсчитывает количество вхождений каждого значения в массиве:
print(‘Labels counts in label:’, np.bincount(label))
print(‘Labels counts in l_train:’, np.bincount(l_train))
print(‘Labels counts in l_test:’, np.bincount(l_test))
print(‘Labels % in label:’, round(label.mean(),2))
print(‘Labels % in l_train:’, round(l_train.mean(),2))
print(‘Labels % in l_test:’, round(l_test.mean(),2))
Labels counts in y: [103 18] Labels counts in y_train: [61 11] Labels counts in y_test: [42 7] Labels % in y: 0.15 Labels % in y_train: 0.15 Labels % in y_test: 0.14
Как видим, процент меток класса очень близки в основном наборе и в обучающей и тестовой частях.
Для решения нашей задачи используем классификатор на основе случайного леса из sklearn.ensemble и сразу выведем метрики качества включая матрицу ошибок
from sklearn.ensemble import RandomForestClassifier
from sklearn import metrics
from sklearn.metrics import plot_confusion_matrix
rfr_model = RandomForestClassifier(n_estimators=10, random_state=0)
rfr_model.fit(d_train, l_train)
predicted = rfr_model.predict(d_test)
score = metrics.accuracy_score(l_test, predicted)
print(f’Decision Tree Classification Score = {score:.1%}\n’)
print(f’Classification Report:\n {metrics.classification_report(l_test, predicted)}\n’)
plot_confusion_matrix(rfr_model, d_test, l_test, cmap=’Greys’)
Decision Tree Classification Score = 89.8% Classification Report: precision recall f1-score support 0 0.93 0.95 0.94 42 1 0.67 0.57 0.62 7 accuracy 0.90 49 macro avg 0.80 0.76 0.78 49 weighted avg 0.89 0.90 0.89 49
Как видим, результат неплохой, доля правильных ответов 89.8%
Из реализация случайных лесов в scikit-leam можно вывести показатели важности признаков, через атрибут feature importances после подгонки классификатора RandomForestClassifier.
features_response = d_test.columns.tolist()
feat_imp_df = pd.DataFrame({
‘Importance’:rfr_model.feature_importances_},
index=features_response)
feat_imp_df.sort_values(‘Importance’, ascending=True).plot.barh()
В результате выполнения кода отображается график, который располагает признаки в наборе данных по их относительной важности, при этом показатели важности признаков нормализованы, в сумме давая 1.0.
По степени важности вперед вышли два показателя — стаж и количество отработанных смен, это говорит о том, что для сокращения текучести кадров необходимо оказывать повышенное внимание новым сотрудникам.
Содержание
- Модели классификации в Python
- Метрики качества классификаторов
- Матрица ошибок (Confusion matrix)
- Линейная классификация
- Логистическая регрессия
Модели классификации в Python
Классификация – отнесение объекта к одной из категорий на основании его признаков.
Метрики качества классификаторов
Матрица ошибок (Confusion matrix)
Матрица ошибок — это способ разбить объекты на четыре категории в зависимости от комбинации истинного ответа и ответа алгоритма.
- T P TP TP — истино-положительное решение;
- T N TN TN — истино-отрицательное решение;
- F P FP FP — ложно-положительное решение (Ошибка первого рода);
- F N FN FN — ложно-отрицательное решение (Ошибка второго рода).
Accuracy — доля правильных ответов:
A c c u r a c y = T P + T N P + N = T P + T N T P + T N + F P + F N = +\mathrm >>= +\mathrm > <\mathrm +\mathrm +\mathrm +\mathrm >>> Accuracy = P + N TP + TN = TP + TN + FP + FN TP + TN
Данная матрика имеет существенный недостаток — её значение необходимо оценивать в контексте баланса классов. Eсли в выборке 950 отрицательных и 50 положительных объектов, то при абсолютно случайной классификации мы получим долю правильных ответов 0.95. Это означает, что доля положительных ответов сама по себе не несет никакой информации о качестве работы алгоритма a(x), и вместе с ней следует анализировать соотношение классов в выборке.
Гораздо более информативными критериями являются точность (precision) и полнота (recall).
Точность показывает, какая доля объектов, выделенных классификатором как положительные, действительно является положительными: P r e c i s i o n = T P T P + F P Precision = \frac P rec i s i o n = TP + FP TP
Полнота показывает, какая часть положительных объектов была выделена классификатором: R e c a l l = T P T P + F N Recall = \frac R ec a ll = TP + FN TP
Существует несколько способов получить один критерий качества на основе точности и полноты. Один из них — F-мера, гармоническое среднее точности и полноты: F_\beta = (1 + \beta^2) \cdot \frac <\mathrm \cdot \mathrm><(\beta^2 \cdot \mathrm) + \mathrm> = \frac <(1 + \beta^2) \cdot \mathrm> <(1 + \beta^2) \cdot \mathrm+ \beta^2 \cdot \mathrm + \mathrm>
Среднее гармоническое обладает важным свойством — оно близко к нулю, если хотя бы один из аргументов близок к нулю. Именно поэтому оно является более предпочтительным, чем среднее арифметическое (если алгоритм будет относить все объекты к положительному классу, то он будет иметь recall = 1 и precision больше 0, а их среднее арифметическое будет больше 1/2, что недопустимо).
Чаще всего берут β = 1 \beta=1 β = 1 , хотя иногда встречаются и другие модификации. F 2 F_2 F 2 острее реагирует на recall (т. е. на долю ложноположительных ответов), а F 0.5 F_ F 0.5 чувствительнее к точности (ослабляет влияние ложноположительных ответов).
В sklearn есть удобная функция sklearn.metrics.classification_report, возвращающая recall, precision и F-меру для каждого из классов, а также количество экземпляров каждого класса.
import pandas as pd import seaborn as sns from sklearn.model_selection import train_test_split from sklearn.preprocessing import StandardScaler, LabelEncoder from sklearn.linear_model import LogisticRegression from sklearn.svm import SVC from sklearn import datasets import numpy as np import matplotlib.pyplot as plt import seaborn as sns %matplotlib inline %config InlineBackend.figure_format = 'svg'
from sklearn.metrics import classification_report y_true = [0, 1, 2, 2, 2] y_pred = [0, 0, 2, 2, 1] target_names = ['class 0', 'class 1', 'class 2'] print(classification_report(y_true, y_pred, target_names=target_names))
precision recall f1-score support
class 0 0.50 1.00 0.67 1
class 1 0.00 0.00 0.00 1
class 2 1.00 0.67 0.80 3
micro avg 0.60 0.60 0.60 5
macro avg 0.50 0.56 0.49 5
weighted avg 0.70 0.60 0.61 5
Линейная классификация
Основная идея линейного классификатора заключается в том, что признаковое пространство может быть разделено гиперплоскостью на две полуплоскости, в каждой из которых прогнозируется одно из двух значений целевого класса. Если это можно сделать без ошибок, то обучающая выборка называется линейно разделимой.
Указанная разделяющая плоскость называется линейным дискриминантом.
Логистическая регрессия
Логистическая регрессия является частным случаем линейного классификатора, но она обладает хорошим «умением» – прогнозировать вероятность отнесения наблюдения к классу. Таким образом, результат логистической регрессии всегда находится в интервале [0, 1].
Источник