Этот учебник является продолжением раздела Методы кластеризации с данными об экспрессии генов, в котором мы обсуждали алгоритм, относящийся к раздельной кластеризации (k-средних) и иерархической кластеризации.

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

Вы можете найти первую часть: здесь.

Структура руководства:

  1. Введение в кластеризацию на основе плотности
  2. Реализация DBSCAN в Python
  3. Введение в модели гауссовских смесей
  4. Реализация гауссовских смешанных моделей в Python
  5. Предыдущие руководства
  6. библиография

Введение в кластеризацию на основе плотности

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

Как мы видели ранее, иерархическая кластеризация и кластеризация k-средних не очень способны формировать кластеры на основе различной плотности или когда они имеют сложную форму. Как мы видим на рисунке 1, интуитивно мы видим три разных концентрических кластера с некоторыми точками данных шума. Иерархическая кластеризация и k-средние в этом случае плохо работают, они плохо обнаруживают ни кластеры, ни шум. В этом случае DBSCAN хорошо идентифицирует кластеры и выделяет шум как выбросы.

Пространственная кластеризация приложений с шумом на основе плотности (DBSCAN) – это алгоритм кластеризации данных, предложенный Мартином Эстером, Хансом-Питером Кригелем, Йоргом Сандером и Сяовей Сюй в 1996 году (1). DBSCAN группирует наблюдения данных в плотном регионе в кластер. Одним из самых больших преимуществ является то, что он очень устойчив к выбросам (точки данных, которые представляют экстремальное значение или сильно отличаются от других точек данных). Общая идея кластеризации на основе плотности заключается в том, что для заданной точки данных вы хотите найти точки данных в пространстве вокруг (ε-окрестность, которая определяется как набор точек данных на расстоянии ε от точки данных p). в двумерном пространстве геометрическая форма, в которой все точки равноудалены от центра, представляет собой круг. В этом случае мы можем описать понятие плотности как количество точек данных в определенной области (радиус ε в данном контексте).

DBSCAN требует, чтобы мы установили только два параметра:

  • Eps или эпсилон (ε) — это радиус круга, который необходимо создать вокруг каждого наблюдения, чтобы проверить плотность. В пространстве более высокого измерения (более 3 объектов или измерений) эпсилон — это радиус гиперсферы. В конкретном случае эпсилон — это расстояние между двумя соседними наблюдениями данных (две точки данных считаются соседями, если расстояние между ними меньше или равно eps).
  • minPts количество наблюдений за соседними данными внутри круга для этого наблюдения, чтобы считаться базовой точкой. Это также можно суммировать как минимальное количество наблюдений данных в кластере (или плотном регионе).

Учитывая два параметра, мы можем определить точки данных как:

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

Давайте посмотрим на рисунок 2, эпсилон представляет радиус кругов, а minPts равно 4. Красные точки — это основные точки, потому что внутри радиуса эпсилон есть как минимум 4 точки. Желтые точки — это пограничные точки, потому что они доступны из центральной точки (внутри круга центральной точки), но внутри их окружающей области нет, по крайней мере, точек. Внутри круга B всего 2 точки, а это означает, что внутри радиуса эпсилона, начиная с B, количество точек данных в окружающей области меньше minPts. Вместо этого N является выбросом, поскольку не является базовой точкой, но также недостижим из любой базовой точки.

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

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

Принимая во внимание эти две концепции, в DBSCAN две точки данных можно назвать:

  • Прямая доступность по плотности
  • Достижимая плотность
  • Подключение по плотности

Точка X считается непосредственно достижимой по плотности из точки Y, если X находится по соседству с Y (расстояние между X и Y меньше эпсилон) и Y является центральной точкой, как на рис. 3А.

Точка X достижима по плотности из Y, если существует цепочка точек (P1, P2, P3 и т. д.), каждая из которых достижима по плотности. Здесь X не является плотностью, достижимой напрямую из Y, но, как мы видим в цепочке, P3 достижима по плотности из Y, P2 достижима по плотности из P3, а X достижима по плотности из P2 (рис. 3 B). В качестве аналогии, две точки данных, достижимые по плотности как «друг друга». Размерность эпсилон определяет, достижимы ли две точки данных по плотности, при выборе меньшего ε плотность достижима меньше точек данных.

Точка плотно связана из точки Y, если существует точка O такая, что и X, и Y достижимы по плотности из O

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

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

Из этой информации мы можем сделать вывод, что DBSCAN очень чувствителен к значению epsilon и MinPts, вариации этих значений значительно изменят вывод.

минуты. Он должен быть установлен выше 3 (если вы установите равным 1, каждая точка будет считаться отдельным кластером). Интересно, что выбор minPts равным 2 подобен выполнению иерархического кластера с метрикой одной ссылки (и поскольку дендрограмма была сокращена на расстоянии эпсилон). Как правило, minPts выбирается как удвоенное количество функций. Большие значения minPts лучше подходят для набора данных с шумом, и, как правило, чем больше ваш набор данных, тем больше должно быть minPts.

Эпсилон обычно выбирается в качестве колена на графике k-расстояний (подробнее о нем будет рассказано позже). Если эпсилон мал, создается большее количество кластеров, и многие другие точки данных будут рассматриваться как шум. Слишком маленький эпсилон будет означать, что большая часть ваших данных не будет сгруппирована и будет считаться выбросами. Наоборот, если эпсилон слишком велик, многие маленькие кластеры будут объединены в один больший кластер. Если вы выберете слишком большой эпсилон, вы можете иметь большинство точек данных только в одном кластере).

Реализация DBSCAN в python

Если мы используем вариант по умолчанию, результаты не обнадеживают:

from sklearn.cluster import DBSCAN
dbscan=DBSCAN()
dbscan.fit(df)
pca_df["dbscan_labels"] = dbscan.labels_
sns.set(style="whitegrid", palette="muted")
ax = sns.scatterplot(x="PCA1", y="PCA2", hue="dbscan_labels",  data=pca_df)
plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)

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

from sklearn.neighbors import NearestNeighbors
neigh = NearestNeighbors(n_neighbors=2)
nbrs = neigh.fit(df)
distances, indices = nbrs.kneighbors(df)
# Plotting K-distance Graph
distances = np.sort(distances, axis=0)
distances = distances[:,1]
#plt.figure(figsize=(20,10))
plt.plot(distances)
plt.title(‘K-distance Graph’,fontsize=20)
plt.xlabel(‘Data Points sorted by distance’,fontsize=14)
plt.ylabel(‘Epsilon’,fontsize=14)
plt.show()

Мы выбираем значение локтя, значение эпсилон, где есть максимальная кривизна. Здесь около 60. minPts также зависит от знания домена (здесь мы выбираем 50)

dbscan=DBSCAN(eps=60,min_samples=50)
dbscan.fit(df)
pca_df[“dbscan_labels”] = dbscan.labels_
sns.set(style=”whitegrid”, palette=”muted”)
ax = sns.scatterplot(x=”PCA1", y=”PCA2", hue=”dbscan_labels”, data=pca_df)
plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)

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

from sklearn.neighbors import NearestNeighbors
neigh = NearestNeighbors(n_neighbors=2)
nbrs = neigh.fit(X)
distances, indices = nbrs.kneighbors(X)
# Plotting K-distance Graph
distances = np.sort(distances, axis=0)
distances = distances[:,1]
#plt.figure(figsize=(20,10))
plt.plot(distances)
plt.title(‘K-distance Graph’,fontsize=20)
plt.xlabel(‘Data Points sorted by distance’,fontsize=14)
plt.ylabel(‘Epsilon’,fontsize=14)
plt.show()

dbscan=DBSCAN(eps=40,min_samples=30)
dbscan.fit(X)
pca_df[“dbscan_labels”] = dbscan.labels_
pca_df[‘dbscan_labels’] = pca_df.dbscan_labels.astype(‘category’)
sns.set(style=”whitegrid”, palette=”muted”)
ax = sns.scatterplot(x=”PCA1", y=”PCA2", hue=”dbscan_labels”, data=pca_df)
plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)

Некоторые соображения:

  • DBSCAN не требует, чтобы мы указывали количество кластеров перед использованием.
  • Он хорошо работает с кластерами произвольной формы. DBSCAN может находить нелинейно разделяемые кластеры. он также может идентифицировать кластеры, которые полностью окружены другими кластерами (но не соединены). Варьируя minPts, можно также добиться хороших результатов в разделении в случае однозвенного эффекта (когда два разных кластера соединены тонкой линией точек и могут рассматриваться как одиночные кластеры).
  • Он устойчив к выбросам и может использоваться для обнаружения выбросов.
  • Требуется только установка двух параметров. Однако иногда бывает непросто определить подходящее значение параметров, и для этого требуются знания предметной области.
  • Он работает хорошо, если кластеры сильно различаются по плотности. Это связано с тем, что трудно найти комбинацию epsilon и minPts, подходящую для всех кластеров.
  • Это не совсем детерминировано, граничная точка, достижимая из более чем одного кластера, может быть вставлена ​​в тот или иной кластер, в зависимости от порядка обработки данных. В то время как это полностью детерминировано для основной точки и выбросов.
  • Качество DBSCAN зависит от выбранной меры расстояния (обычно евклидово расстояние). Евклидово расстояние не является масштабно-инвариантным, что означает, что оно чувствительно к масштабированию функций (нормализация ваших данных - лучшая практика). С данными высокой размерности (где количество признаков намного превышает количество наблюдений) евклидово расстояние не дает хороших результатов. С увеличением количества размерностей снижается качество результата (проклятие размерности).

Введение в модели гауссовских смесей

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

Начнем с распределения Гаусса (также известного как нормальное распределение). Это кривая в форме колокола, где точки данных симметрично распределены вокруг среднего значения. Распределение Гаусса определяется средним (µ) и дисперсией (σ2), чем выше дисперсия, тем более разбросана кривая (влияющая на ширину колокола). ). Другое среднее значение означает, что гауссиана смещена на рисунке ниже, есть 4 разных гауссовских распределения с разными средними значениями и дисперсией.

Изображение у нас есть три распределения Гаусса (которые мы можем определить как GD1, GD2 и GD3), и каждое из них имеет свое собственное среднее значение (µ1, µ2, µ3) и свою собственную дисперсию (σ21, σ22, σ23, σ — стандартное отклонение и квадрат стандартного отклонения – это дисперсия, σ2). При рассмотрении набора данных (набора точек данных) гауссовая смешанная модель будет определять вероятность того, что точка данных принадлежит каждому из этих распределений. Это связано с тем, что GMM являются вероятностными моделями.

В случае, если у нас есть три кластера (обозначенных разными цветами) и мы выбираем точку данных, мы можем рассчитать вероятность принадлежности к определенному кластеру. Как показано на рис. 5А, вероятность точки данных, обведенной красным, для синего кластера равна 1. Если мы выберем другую точку данных, мы увидим, что она с большей вероятностью принадлежит голубому кластеру. Обе точки данных не имеют вероятности принадлежать зеленому кластеру.

В одномерном пространстве для распределения Гаусса функция плотности вероятности (функция плотности вероятности является мерой того, насколько вероятно, что переменная будет иметь определенное значение):

Если мы посмотрим на этот пример (рис. 6), у нас есть гауссово распределение со средним значением с центром в 0 и стандартным отклонением 1, и это описано функцией выше. При 0 вероятность наблюдать вход x с этим конкретным гауссовым распределением составляет 40%. Чем дальше мы удаляемся от центра, тем ниже вероятность, поскольку ключевое свойство гауссовского распределения состоит в том, что наибольшая вероятность находится в центре и уменьшается по мере удаления от центра, приближаясь к нулю. Сумма значений под кривой равна 1 (вероятность 100 %).

Это верно, когда у нас есть только одна переменная, как на рис. 4, а если у нас есть две переменные вместо 2-мерной кривой, мы имеем 3-мерную кривую. В этом случае вход больше не скаляр, а вектор, и у нас также есть вектор для среднего (среднее значение представляет собой центр наших данных, поэтому оно имеет ту же размерность, что и вход). Вместо этого ковариация представляет собой матрицу (в данном случае 2 x 2). Ковариационная матрица также хранит информацию о связи между x и y, о том, как изменяется y при изменении x. На следующем рисунке мы видим, что чем трехмерная кривая, интенсивность цвета представляет более высокую вероятность.

в этом случае функция плотности вероятности будет:

Как уже говорилось, здесь µ — двумерный вектор, а Σ — ковариационная матрица (2x2), а для двумерной кривой ковариация определяет форму кривой. Мы можем обобщить эту концепцию на многомерный набор данных (скажем, число измерений N), где гауссовская модель будет иметь µ в качестве вектора (длина этого вектора будет N) и ковариационную матрицу (N x N как матрица). С учетом набора данных это будет описываться смесью распределений Гаусса (количество GD будет равно количеству кластеров) и каждого GD с вектором µ и матрицей σ.

Под капотом Gaussian Mixture Model используется алгоритм максимизации ожидания (с которым мы уже сталкивались при обсуждении k-средних). Алгоритм максимизации ожидания (EM) представляет собой статистическую модель, которая определяет скрытые переменные (в данном случае средний вектор и ковариационную матрицу гауссовых распределений, представляющих кластеры). Идея заключается в том, что EM использует существующие данные для поиска оптимальных значений скрытых переменных.

Мы хотим назначить k кластеров: так, у нас есть K гауссовского распределения, каждый из которых имеет вектор µ и матрицу Σ (или ковариационную матрицу). Кроме того, у нас есть еще один параметр, который представляет количество точек, присутствующих в распределении, или плотность распределения (Π). Нам нужно найти значение этих параметров для определения каждого распределения Гаусса. На первом шаге мы присваиваем случайное значение каждому из этих параметров.

На первом этапе (ожидание) мы вычисляем вероятность для каждой точки того, что она принадлежит каждому кластеру (или распределению) c1, c2… ck. GMM использует теорему Байеса для расчета этой вероятности.

Значение высокое, когда точка относится к правому кластеру.

После этого шага обновляются значения µ, ∏ и Σ. Матрицы µ и Σ обновляются пропорционально значениям вероятности для точки данных (в частности, если точка данных имеет более высокую вероятность быть частью этого гауссовского распределения, она будет вносить больший вклад).

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

Прежде чем приступить к коду, предусмотрительно, GMM был описан как алгоритм кластеризации, но по сути это алгоритм оценки плотности. Когда вы подгоняете данные к GMM, технически генерируется не кластеризация, а генеративная вероятностная модель, описывающая распределение данных (2). Мы вернемся к этому последнему пункту.

Реализация гауссовских смешанных моделей в Python

Основные параметры:

  • n_components: по умолчанию 1, количество компонентов гауссовской смеси, которые мы выбираем
  • covariance_type: по умолчанию «полный», описывает тип параметра ковариации.

Подгоняем данные по компонентам PCA.

from sklearn.decomposition import PCA
pca = PCA(n_components=50)
X = pca.fit(df).transform(df)
#we fit the data
from sklearn.mixture import GaussianMixture
gmm = GaussianMixture(n_components=7)
gmm.fit(X)
#we plot the results
pca_df["GMM_labels"] = gmm.predict(X)
pca_df['GMM_labels'] = pca_df.GMM_labels.astype('category')
sns.set(style="whitegrid", palette="muted")
ax = sns.scatterplot(x="PCA1", y="PCA2", hue="GMM_labels",  data=pca_df)
plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)

Мы выбрали 7 компонентов на основе наших предварительных знаний, но это не всегда возможно, поэтому есть и другие методы. Мы будем использовать 2, информационный критерий Акаике (AIC) или байесовский информационный критерий (BIC).

  • Информационный критерий Акаике (AIC) — это оценка ошибки прогнозирования и, следовательно, показатель качества статистической модели для набора данных. Применительно к коллекции моделей позволяет оценить качество каждой модели, предоставляя инструмент для выбора модели. AIC оценивает относительный объем информации, которая теряется каждой моделью, из чего можно сделать вывод, что чем меньше информации теряет модель, тем выше ее качество. На практике каждая модель теряет информацию по отношению к так называемой «истинной модели», и мы хотим выбрать модель, которая теряет меньше информации.
  • Байесовский информационный критерий (BIC) — еще один критерий выбора модели. Увеличение параметров приводит к переоснащению, AIC и BIC описывают этот компромисс и добавляют штраф к увеличению параметра, чтобы избежать этой проблемы (штраф в BIC больше, чем в AIC).
n_components = np.arange(1, 21)
models = [GaussianMixture(n, covariance_type=’full’, random_state=0).fit(X) for n in n_components]
plt.plot(n_components, [m.bic(X) for m in models], label=’BIC’)
plt.plot(n_components, [m.aic(X) for m in models], label=’AIC’)
plt.legend(loc=’best’)
plt.xlabel(‘n_components’)

Оптимальное значение этого параметра находится там, где AIC и BIC минимальны. AIC говорит, что около 7 — это хороший выбор, но чем больше компонентов, тем лучше. BIC предлагает диапазон от 4 до 7, потому что BIC имеет более высокий штраф и предлагает более простую модель. Итак, кажется, что семь - хороший выбор. Предостережение: этот график измеряет, насколько хорошо GMM работает как средство оценки плотности, а не как хорошо он работает как алгоритм кластеризации (2).

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

new_data = gmm.sample(1000)
#since our data are returned in a tuple
gmm_components = new_data[:][0]
gmm_new_labels = new_data[:][1]
new_df = pd.DataFrame(columns = ["x", "y", "gmm_new_labels"])
new_df["PCA1"] = gmm_components[:, 0]
new_df["PCA2"] = gmm_components[:, 1]
new_df["gmm_new_labels"] = gmm_new_labels
new_df['gmm_new_labels'] = new_df.gmm_new_labels.astype('category')
sns.set(style="whitegrid", palette="muted")
ax = sns.scatterplot(x="PCA1", y="PCA2", hue="gmm_new_labels",  data=new_df)
plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)

График похож на наш исходный график (даже если таких пациентов не существует).

Вы можете проверить здесь лучшее использование генеративной модели здесь, где они использовали цифры ручной работы.

некоторые соображения:

  • GMM также является генеративной моделью, это может быть удобно для некоторых задач.
  • GMM могут изучать кластеры любой эллиптической формы. GMM позволяет проявлять гибкость в отношении ковариации кластеров. K-средние — это особый случай GMM, когда ковариация каждого кластера по всем измерениям приближается к нулю.
  • GMM предоставляют вероятности того, как каждая точка данных связана с каждым кластером. На практике GMMS допускает, что каждое наблюдение каким-то образом принадлежит к большему, чем кластер, в пределах определенной степени неопределенности. GMM в основном представляют кластеры как распределение вероятностей. Смешанное членство в кластере (точки данных, принадлежащие к разным кластерам в разной степени) полезно для некоторых задач (вы можете изобразить новостную статью, принадлежащую нескольким тематическим кластерам). В k-средних вероятность того, что одна точка принадлежит одному кластеру, равна 1 (0 для всех остальных кластеров, что еще раз объясняет, почему k-средние практически являются частным случаем GMM), и это причина, по которой k-средние генерируют только сферические кластеры. . Это также называется мягкой классификацией по сравнению с жесткой классификацией (k-means).
  • Нам все еще нужно заранее определить количество кластеров, которые мы хотим.
  • Алгоритм оптимизации функции потерь GMM нетривиален и сложен (наиболее популярен метод максимизации ожидания). Алгоритм может столкнуться с локальными максимумами, которые не являются лучшими в глобальном масштабе.
  • GMM трудно инициализировать, когда у вас есть набор данных с большой размерностью (гораздо больше функций, чем наблюдений). По этой причине мы начали с компонентов PCA.
  • Центральная предельная теорема утверждает, что по мере того, как мы собираем все больше и больше выборок (или наблюдений) для набора данных, он будет напоминать распределение Гаусса.

Доступность кода

Код урока и все изображения вы найдете здесь. Вы можете использовать предоставленный блокнот Colab без необходимости установки каких-либо пакетов на свой компьютер.

если вам было интересно:

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

Предыдущие руководства

О техниках снижения сложности и визуализации: здесь

Для получения дополнительной информации:

Об остром миелоидном лейкозе: здесь

Об анализе медицинских изображений ИИ: здесь

Об искусственном интеллекте на биологических данных: здесь