СТАТЬЯ

Тематическое исследование: диагностика рака молочной железы

Из Методы ансамбля для машинного обучения Гаутама Кунапули

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

___________________________________________________________

Получите скидку 40% на Методы ансамбля для машинного обучения, введя fcckunapuli в поле кода скидки при оформлении заказа на сайте manning.com.
_____________________________________________________________

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

Эта задача является сложной из-за незначительных различий в данных пациентов, из-за чего для одной модели очень сложно достичь высокой производительности. Бэггинг и случайный лес способны обучать эффективные модели. Это тематическое исследование будет рассмотрено с использованием реализации scikit-learn бэггинга и случайных лесов. Производительность этих изученных ансамблевых моделей будет сравниваться с индивидуальными (не ансамблевыми/традиционными) моделями машинного обучения. Результаты будут использованы, чтобы показать, что ансамблевые методы могут значительно превзойти индивидуальные модели, тем самым предоставив читателю конкретный реальный пример того, как ансамбли могут быть успешными.

Листинг 1. Загрузка и предварительная обработка

In [1]:
 from sklearn.datasets import load_breast_cancer
 dataset = load_breast_cancer()
  
 # Convert to a Pandas DataFrame so we can visualize nicely
 import pandas as pd
 import numpy as np
 df = pd.DataFrame(data=dataset['data'], columns=dataset['feature_names'])
 i = np.random.permutation(len(dataset['target']))
 df = df.iloc[i, :7]
 df['diagnosis'] = dataset['target'][i]
 df = df.reset_index()
 df.head()
  
 Out[1]:

Предварительно обработайте данные, стандартизировав признаки как нулевое среднее значение и единичное стандартное отклонение.

In [2]:
 from sklearn.preprocessing import StandardScaler
 X, y = dataset['data'], dataset['target']
 X = StandardScaler().fit_transform(X)

Бэггинг, случайный лес и дополнительные деревья

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

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

scikit-learn предоставляет ExtraTreesClassifier, который поддерживает нестандартную оценку и распараллеливание, подобно BaggingClassifier и RandomForestClassifier. Обратите внимание, что ExtraTrees обычно не выполняют начальную выборку (bootstrap = False по умолчанию), поскольку мы можем добиться разнообразия базовых оценок за счет экстремальной рандомизации.

После того, как мы предварительно обработаем наш набор данных, мы обучим и оценим пакетирование с помощью деревьев решений, случайных лесов и ExtraTrees, чтобы ответить на следующие вопросы:

  1. Как меняется производительность ансамбля в зависимости от размера ансамбля? То есть, что происходит, когда наши ансамбли становятся все больше и больше?
  2. Как производительность ансамбля меняется с базовой сложностью обучаемого? То есть, что происходит, когда наши индивидуальные базовые оценки становятся все более и более сложными. В этом тематическом исследовании, поскольку все три рассмотренных метода ансамбля используют деревья решений в качестве базовых оценок, одной «мерой» сложности является глубина дерева, причем более глубокие деревья являются более сложными.

Размер ансамбля в сравнении с исполнением ансамбля

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

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

Конкретный набор данных, который мы будем использовать, — это набор данных о раке молочной железы штата Висконсин (WDBC), общий эталонный набор данных в машинном обучении. С 1993 года данные WDBC использовались для оценки производительности десятков алгоритмов машинного обучения.

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

WDBC доступен в scikit-learn. Мы сравниваем поведение трех алгоритмов при увеличении параметра n_estimators.

Листинг 2. Ошибки обучения и тестирования при увеличении размера ансамбля.

In [3]:
 import pickle
 import os
  
 from sklearn.model_selection import train_test_split
 from sklearn.tree import DecisionTreeClassifier
 from sklearn.ensemble import BaggingClassifier, RandomForestClassifier, ExtraTreesClassifier
 from sklearn.metrics import accuracy_score
  
 max_leaf_nodes = 8
 n_runs = 20
 n_estimator_range = range(2, 20, 1)
  
 bag_trn_error = np.zeros((n_runs, len(n_estimator_range)))  # Train error for the bagging classifier
 rf_trn_error = np.zeros((n_runs, len(n_estimator_range)))   # Train error for the random forest classifier
 xt_trn_error = np.zeros((n_runs, len(n_estimator_range)))   # Train error for the extra trees classifier
  
 bag_tst_error = np.zeros((n_runs, len(n_estimator_range)))  # Test error for the bagging classifier
 rf_tst_error = np.zeros((n_runs, len(n_estimator_range)))   # Test error for the random forest classifier
 xt_tst_error = np.zeros((n_runs, len(n_estimator_range)))   # Test error for the extra trees classifier
  
 for run in range(0, n_runs):
     print('Run {0}'.format(run))
  
     # Split into train and test sets
     X_trn, X_tst, y_trn, y_tst = train_test_split(X, y, test_size=0.25)
  
     for j, n_estimators in enumerate(n_estimator_range):
         # Train using bagging
         bag_clf = BaggingClassifier(base_estimator=DecisionTreeClassifier
                                                    (max_leaf_nodes=max_leaf_nodes),
                                     n_estimators=n_estimators, max_samples=0.5, n_jobs=-1)
         bag_clf.fit(X_trn, y_trn)
         bag_trn_error[run, j] = 1 - accuracy_score(y_trn, bag_clf.predict(X_trn))
         bag_tst_error[run, j] = 1 - accuracy_score(y_tst, bag_clf.predict(X_tst))
  
         # Train using random forests
         rf_clf = RandomForestClassifier(max_leaf_nodes=max_leaf_nodes,
                                         n_estimators=n_estimators, n_jobs=-1)
         rf_clf.fit(X_trn, y_trn)
         rf_trn_error[run, j] = 1 - accuracy_score(y_trn, rf_clf.predict(X_trn))
         rf_tst_error[run, j] = 1 - accuracy_score(y_tst, rf_clf.predict(X_tst))
  
         # Train using extra trees
         xt_clf = ExtraTreesClassifier(max_leaf_nodes=max_leaf_nodes, bootstrap=True,
                                       n_estimators=n_estimators, n_jobs=-1)
         xt_clf.fit(X_trn, y_trn)
         xt_trn_error[run, j] = 1 - accuracy_score(y_trn, xt_clf.predict(X_trn))
         xt_tst_error[run, j] = 1 - accuracy_score(y_tst, xt_clf.predict(X_tst))
  
     results = (bag_trn_error, bag_tst_error, \
                rf_trn_error, rf_tst_error, \
                xt_trn_error, xt_tst_error)
  
  
 In [4]:
 %matplotlib inline
  
 import matplotlib.pyplot as plt
  
 fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(10, 4))
 n_estimator_range = range(2, 20, 1)
  
 # Plot the training error
 m = np.mean(bag_trn_error*100, axis=0)
 ax[0].plot(n_estimator_range, m, linewidth=3)
  
 m = np.mean(rf_trn_error*100, axis=0)
 ax[0].plot(n_estimator_range, m, linestyle='--', linewidth=3)
  
 m = np.mean(xt_trn_error*100, axis=0)
 ax[0].plot(n_estimator_range, m, linestyle='-.', linewidth=3)
  
 ax[0].legend(['Bagging', 'Random Forest', 'Extra Trees'], fontsize=16)
 ax[0].set_xlabel('Number of Estimators', fontsize=16)
 ax[0].set_ylabel('Training Error (%)', fontsize=16)
 ax[0].axis([0, 21, 1, 9])
 ax[0].grid()
  
 # Plot the test error
 m = np.mean(bag_tst_error*100, axis=0)
 ax[1].plot(n_estimator_range, m, linewidth=3)
  
 m = np.mean(rf_tst_error*100, axis=0)
 ax[1].plot(n_estimator_range, m, linestyle='--', linewidth=3)
  
 m = np.mean(xt_tst_error*100, axis=0)
 ax[1].plot(n_estimator_range, m, linestyle='-.', linewidth=3)
  
 ax[1].legend(['Bagging', 'Random Forest', 'Extra Trees'], fontsize=16)
 ax[1].set_xlabel('Number of Estimators', fontsize=16)
 ax[1].set_ylabel('Hold-out Test Error (%)', fontsize=16)
 ax[1].grid()
 ax[1].axis([0, 21, 1, 9]);
 plt.tight_layout()
 Out[4]:

Базовая сложность учащегося и совокупная производительность

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

Хорошей практикой всегда является стандартизация обучающих данных, тем более что масштабы функций сильно различаются. Стандартизация — это этап предварительной обработки, при котором функции масштабируются так, чтобы среднее значение равнялось 0, а стандартное отклонение — 1.

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

In [5]:
 n_estimators = 10
 max_leaf_nodes = 8
 n_runs = 20
 n_leaf_range = [2, 4, 8, 16, 24, 32]
  
 bag_trn_error = np.zeros((n_runs, len(n_leaf_range)))  # Train error for the bagging classifier
 rf_trn_error = np.zeros((n_runs, len(n_leaf_range)))   # Train error for the random forest classifier
 xt_trn_error = np.zeros((n_runs, len(n_leaf_range)))   # Train error for the extra trees classifier
  
 bag_tst_error = np.zeros((n_runs, len(n_leaf_range)))  # Test error for the bagging classifier
 rf_tst_error = np.zeros((n_runs, len(n_leaf_range)))   # Test error for the random forest classifier
 xt_tst_error = np.zeros((n_runs, len(n_leaf_range)))   # Test error for the extra trees classifier
  
 for run in range(0, n_runs):
     print('Run {0}'.format(run))
  
     # Split into train and test sets
     X_trn, X_tst, y_trn, y_tst = train_test_split(X, y, test_size=0.25)
  
     for j, max_leaf_nodes in enumerate(n_leaf_range):
         # Train using bagging
         bag_clf = BaggingClassifier(base_estimator=DecisionTreeClassifier
                                                    (max_leaf_nodes=max_leaf_nodes),
                                     n_estimators=n_estimators, max_samples=0.5, n_jobs=-1)
         bag_clf.fit(X_trn, y_trn)
         bag_trn_error[run, j] = 1 - accuracy_score(y_trn, bag_clf.predict(X_trn))
         bag_tst_error[run, j] = 1 - accuracy_score(y_tst, bag_clf.predict(X_tst))
  
         # Train using random forests
         rf_clf = RandomForestClassifier(max_leaf_nodes=max_leaf_nodes,
                                         n_estimators=n_estimators, n_jobs=-1)
         rf_clf.fit(X_trn, y_trn)
         rf_trn_error[run, j] = 1 - accuracy_score(y_trn, rf_clf.predict(X_trn))
         rf_tst_error[run, j] = 1 - accuracy_score(y_tst, rf_clf.predict(X_tst))
  
         # Train using extra trees
         xt_clf = ExtraTreesClassifier(max_leaf_nodes=max_leaf_nodes, bootstrap=True,
                                       n_estimators=n_estimators, n_jobs=-1)
         xt_clf.fit(X_trn, y_trn)
         xt_trn_error[run, j] = 1 - accuracy_score(y_trn, xt_clf.predict(X_trn))
         xt_tst_error[run, j] = 1 - accuracy_score(y_tst, xt_clf.predict(X_tst))
            
     results = (bag_trn_error, bag_tst_error, \
                rf_trn_error, rf_tst_error, \
                xt_trn_error, xt_tst_error)
            
    
  
  
  
 In [6]:
 %matplotlib inline
  
 fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(10, 4))
 n_leaf_range = [2, 4, 8, 16, 24, 32]
  
 # Plot the training error
 m = np.mean(bag_trn_error*100, axis=0)
 ax[0].plot(n_leaf_range, m, linewidth=3, marker='o')
  
 m = np.mean(rf_trn_error*100, axis=0)
 ax[0].plot(n_leaf_range, m, linestyle='--', linewidth=3, marker='o')
  
 m = np.mean(xt_trn_error*100, axis=0)
 ax[0].plot(n_leaf_range, m, linestyle='-.', linewidth=3, marker='o')
  
 ax[0].legend(['Bagging', 'Random Forest', 'Extra Trees'], fontsize=16)
 ax[0].set_xlabel('Depth of Trees in Ensemble', fontsize=16)
 ax[0].set_ylabel('Training Error (%)', fontsize=16)
 ax[0].set_xticks(n_leaf_range)
 ax[0].grid()
  
 # Plot the test error
 m = np.mean(bag_tst_error*100, axis=0)
 ax[1].plot(n_leaf_range, m, linewidth=3, marker='o')
  
 m = np.mean(rf_tst_error*100, axis=0)
 ax[1].plot(n_leaf_range, m, linestyle='--', linewidth=3, marker='o')
  
 m = np.mean(xt_tst_error*100, axis=0)
 ax[1].plot(n_leaf_range, m, linestyle='-.', linewidth=3, marker='o')
  
 ax[1].legend(['Bagging', 'Random Forest', 'Extra Trees'], fontsize=16)
 ax[1].set_xlabel('Depth of Trees in Ensemble', fontsize=16)
 ax[1].set_ylabel('Hold-out Test Error (%)', fontsize=16)
 ax[1].set_xticks(n_leaf_range)
 ax[1].grid();
 plt.tight_layout()
 Out[6]:

Важность функций в случайных лесах

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

Как меняется производительность ансамбля в зависимости от размера ансамбля? То есть, что происходит, когда наши ансамбли становятся все больше и больше?

Важность функций с корреляциями меток

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

Во-первых, мы рассмотрим, как производительность обучения и тестирования меняется в зависимости от размера ансамбля. То есть мы сравниваем поведение трех алгоритмов при увеличении параметра n_estimators.

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

In [7]:
 %matplotlib inline
  
 import seaborn as sea
  
 df = pd.DataFrame(data=dataset['data'], columns=dataset['feature_names'])
 df['diagnosis'] = dataset['target']
  
 fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(12, 8))
 cor = np.abs(df.corr())
 sea.heatmap(cor, annot=False, cbar=False, cmap=plt.cm.Reds, ax=ax[0])
  
  
 f = ['mean radius', 'mean perimeter', 'mean area', \
      'worst radius', 'worst perimeter', 'worst area', \
      'radius error', 'perimeter error', 'area error', 'diagnosis']
 cor_zoom = np.abs(df[f].corr())
 sea.heatmap(cor_zoom, annot=True, cbar=False, cmap=plt.cm.Reds, ax=ax[1])
 fig.tight_layout()
  
 Out[7]:

Есть несколько признаков, которые также сильно коррелируют с меткой, то есть диагнозом доброкачественности или злокачественности. Давайте определим 10 признаков, наиболее коррелирующих с диагнозом.

In [8]:
 label_corr = cor.iloc[:, -1]
 label_corr.sort_values(ascending=False)[1:11]
 Out[8]:
 worst concave points    0.793566
 worst perimeter         0.782914
 mean concave points     0.776614
 worst radius            0.776454
 mean perimeter          0.742636
 worst area              0.733825
 mean radius             0.730029
 mean area               0.708984
 mean concavity          0.696360
 worst concavity         0.659610
 Name: diagnosis, dtype: float64

Важность функций с использованием случайных лесов

Случайные леса также могут обеспечивать важность функций. Листинг ниже иллюстрирует это.

Листинг 3. Важность функций в наборе данных WDBC с использованием случайных лесов

In [9]:
 X_trn, X_tst, y_trn, y_tst = train_test_split(X, y, test_size=0.15)
 n_features = X_trn.shape[1]
  
 rf = RandomForestClassifier(max_leaf_nodes=24, n_estimators=50, n_jobs=-1)
 rf.fit(X_trn, y_trn)
 err = 1 - accuracy_score(y_tst, rf.predict(X_tst))
 print('Prediction Error = {0:4.2f}%'.format(err*100))
  
 importance_threshold = 0.02
 for i, (feature, importance) in enumerate(zip(dataset['feature_names'],
                                               rf.feature_importances_)):
  
     if importance > importance_threshold:
         print('[{0}] {1} (score={2:4.3f})'.format(i, feature, importance))
 Prediction Error = 4.65%
 [0] mean radius (score=0.077)
 [1] mean texture (score=0.020)
 [2] mean perimeter (score=0.045)
 [3] mean area (score=0.041)
 [5] mean compactness (score=0.026)
 [6] mean concavity (score=0.036)
 [7] mean concave points (score=0.120)
 [13] area error (score=0.020)
 [20] worst radius (score=0.193)
 [22] worst perimeter (score=0.138)
 [23] worst area (score=0.054)
 [26] worst concavity (score=0.027)
 [27] worst concave points (score=0.096)

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

In [10]:
 fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(8, 4))
  
 # Identify the important features
 importance_threshold = 0.02
 idx = np.array(range(n_features))
 imp = np.where(rf.feature_importances_ >= importance_threshold)  # important features
 rest = np.setdiff1d(idx, imp)  # remaining features
  
 # Plot the important features and the rest on a bar chart
 plt.bar(idx[imp], rf.feature_importances_[imp], alpha=0.65)
 plt.bar(idx[rest], rf.feature_importances_[rest], alpha=0.45)
  
 # Print feature names on the bar chart
 for i, (feature, importance) in enumerate(zip(dataset['feature_names'], rf.feature_importances_)):
     if importance > importance_threshold:
         plt.text(i, 0.015, feature, ha='center', va='bottom', rotation='vertical', fontsize=16, fontweight='bold')
         print('[{0}] {1} (score={2:4.3f})'.format(i, feature, importance))
     else:
         plt.text(i, 0.01, feature, ha='center', va='bottom', rotation='vertical', fontsize=16, color='gray')
    
 # Finish the plot   
 fig.axes[0].get_xaxis().set_visible(False)
 plt.ylabel('Feature Importance Score', fontsize=16)
 plt.xlabel('Features for Breast Cancer Diagnosis', fontsize=16);
 plt.tight_layout()
  
 Out[10]:
  
 [0] mean radius (score=0.077)
 [1] mean texture (score=0.020)
 [2] mean perimeter (score=0.045)
 [3] mean area (score=0.041)
 [5] mean compactness (score=0.026)
 [6] mean concavity (score=0.036)
 [7] mean concave points (score=0.120)
 [13] area error (score=0.020)
 [20] worst radius (score=0.193)
 [22] worst perimeter (score=0.138)
 [23] worst area (score=0.054)
 [26] worst concavity (score=0.027)
 [27] worst concave points (score=0.096)

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

Это все для этой статьи.

Если вы хотите узнать больше о книге, вы можете ознакомиться с ней на нашей браузерной платформе liveBook здесь.