Кластеризация данных на основе гауссовских смесей

Введение

В данном примере представлена методика реализации алгоритма жёсткой кластеризации на синтезированных данных, генерируемых смесью гауссовских распределений.

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

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

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

Исходные данные

Импортируем и присоединим необходимые библиотеки.

import Pkg
Pkg.add(["Random", "Distributions", "LinearAlgebra", "Clustering", "GaussianMixtures"])
using Random, Distributions, LinearAlgebra, Clustering, GaussianMixtures

Смоделируем данные из смеси двух двумерных гауссовских распределений, используя функцию MvNormal.

Random.seed!(123)
μ1 = [1.0, 2.0]
σ1 = [3.0 0.2; 0.2 2.0]
μ2 = [-1.0, -2.0]
σ2 = [2.0 0.0; 0.0 1.0]
X1 = rand(MvNormal(μ1, σ1), 200)'
X2 = rand(MvNormal(μ2, σ2), 100)'
X = vcat(X1, X2)
точек = size(X, 1)

Отобразим точечную диаграмму сгенерированных данных.

график1 = scatter(X[:, 1], X[:, 2], markersize=3, markercolor=:black, marker=:o, legend=false, title="Набор точек")#, xlabel="X₁", ylabel="X₂")
display(график1)

Обучение модели для сгенерированных данных

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

gm = GMM(2, X; kind=:full, nIter=5000, nInit=100, method=:kmeans)

GMM{Float64} with 2 components in 2 dimensions and full covariance
Mix 1: weight 0.343617
 mean: [-0.7377815809293755, -1.942919371363109]
 covariance: 2×2 Matrix{Float64}:
 2.12472   0.294642
 0.294642  0.994616
Mix 2: weight 0.656383
 mean: [0.800774284695651, 2.041802685241795]
 covariance: 2×2 Matrix{Float64}:
 3.29572   0.218722
 0.218722  1.8319

Создадим скрипт для одновременного отображения набора точек и контуров.

function отобразить(X, gm)
    к = size(gm.μ, 1)
    μs = gm.μ
    Σ_ = gm.Σ
    веса = gm.w 
    
    function gmPDF(x, y)
        точка = [x, y]
        сумма = 0.0
        for k in 1:к
            μ = μs[k, :]
            L = Σ_[k]
            Σ = L * L'
            распределение = MvNormal(μ, Σ)
            сумма += веса[k] * pdf(распределение, точка)
        end
        return сумма
    end
    
    x_сетка = range(-5, 6, length=80)
    y_сетка = range(-5, 6, length=80)
    
    plot(size=(800, 600))
    scatter!(X[:, 1], X[:, 2], label="", color=:black, markersize=3, marker=:circle, alpha=0.7)
    contour!(x_сетка, y_сетка, (x,y)->gmPDF(x,y), levels=20, color=:jet, linewidth=3.5, colorbar=false, label="", alpha=0.8)
    title!("Набор точек и контуры")
    #xlabel!("X₁")
    #ylabel!("X₂")
    return current()
end

отобразить (generic function with 1 method)

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

график2 = отобразить(X, gm)
display(график2)

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

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

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

function предсказать(gm::GMM, X::Matrix{Float64})
    точек = size(X, 1)
    к = size(gm.μ, 1)
    метки = zeros(Int, точек)
    for i in 1:точек
        точка = X[i, :]
        макс_вероятность = -Inf
        лучший = 1
        for k in 1:к
            μ = gm.μ[k, :]
            L = gm.Σ[k]
            Σ = L * L'
            распределение = MvNormal(μ, Σ)
            логарифм_плотности = logpdf(распределение, точка)
            логарифм_вероятности = log(gm.w[k]) + логарифм_плотности
            if логарифм_вероятности > макс_вероятность
                макс_вероятность = логарифм_вероятности
                лучший = k
            end
        end
        метки[i] = лучший
    end
    return метки
end

function кластеризация(gm::GMM, X::Matrix{Float64})
    return предсказать(gm, X)
end

кластеризация (generic function with 1 method)

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

idx = кластеризация(gm, X)

300-element Vector{Int64}:
 2
 2
 2
 2
 2
 2
 2
 1
 2
 2
 2
 2
 2
 ⋮
 1
 1
 1
 1
 1
 1
 1
 1
 1
 1
 1
 1

Отобразим результат кластеризации.

function диаграмма(X, idx; colors=[:red, :blue], markers=[:+, :o])#, titles=["Кластер 1", "Кластер 2"])
    p = plot(size=(800, 600), legend=:best)
    уникальные = unique(idx)
    for (i, кластер) in enumerate(уникальные)
        маска = idx .== кластер
        x_класт = X[маска, 1]
        y_класт = X[маска, 2]
        форма_маркера = markers[i] == :+ ? :circle : :dtriangle
        scatter!(x_класт, y_класт, markershape=форма_маркера, markercolor=colors[i], markersize=5, label="")#, label=titles[i])
    end
    return p
end

график3 = диаграмма(X, idx)
title!(график3, "Кластеризация")
#xlabel!(график3, "X₁")
#ylabel!(график3, "X₂")
display(график3)

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

Оценка вероятностей принадлежности к кластерам

Создадим специальную функцию для непосредственного расчёта матрицы оценочных вероятностей

function оценка(gm, X)
    точек = size(X, 1)
    к = size(gm.μ, 1)
    P = zeros(точек, к)
    for i in 1:точек
        логарифмы = zeros(к)
        for k in 1:к
            μ = gm.μ[k, :]
            L = gm.Σ[k]
            Σ = L * L'
            распределение = MvNormal(μ, Σ)
            логарифм_плотности = logpdf(распределение, X[i, :])
            логарифмы[k] = log(gm.w[k]) + логарифм_плотности
        end
        макс_лог = maximum(логарифмы)
        смещенные = логарифмы .- макс_лог
        эксп_вероятности = exp.(смещенные)
        P[i, :] = эксп_вероятности ./ sum(эксп_вероятности)
    end
    return P
end

оценка (generic function with 1 method)

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

P = оценка(gm, X)
кластер1 = idx .== 1
кластер2 = idx .== 2

график4 = plot(size=(800, 600))
scatter!(X[кластер1, 1], X[кластер1, 2], marker_z = P[кластер1, 2], markersize = 5, marker = :dtriangle, linewidth = 1.5, clims = (0, 1), color = :jet, colorbar = true, colorbar_title = "Оценочная вероятность", label="")#, label = "Кластер 1")
scatter!(X[кластер2, 1], X[кластер2, 2], marker_z = P[кластер2, 2], markersize = 5, marker = :circle, linewidth = 1.5, clims = (0, 1), color = :jet, label="")#, label = "Кластер 2", colorbar = true, colorbar_title = "Оценочная вероятность")
title!("Оценочные вероятности")
#xlabel!("X₁")
#ylabel!("X₂")
display(график4)

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

Назначения кластерам новых данных

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

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

mix = MixtureModel([MvNormal(μ1, σ1), MvNormal(μ2, σ2)], [0.75, 0.25])
X0 = collect(rand(mix, 75)')

Назначим кластерам новые данные и визуализируем результат.

idx0 = кластеризация(gm, X0)

график5 = plot(size=(800, 600))
к = size(gm.μ, 1)
μs = gm.μ
Σ_ = gm.Σ
веса = gm.w

function gmPDF(x, y)
    точка = [x, y]
    сумма = 0.0
    for k in 1:к
        μ = μs[k, :]
        L = Σ_[k]
        Σ = L * L'
        распределение = MvNormal(μ, Σ)
        сумма += веса[k] * pdf(распределение, точка)
    end
    return сумма
end

x_min, x_max = minimum(X0[:, 1])-1, maximum(X0[:, 1])+1
y_min, y_max = minimum(X0[:, 2])-1, maximum(X0[:, 2])+1
x_сетка = range(x_min, x_max, length=80)
y_сетка = range(y_min, y_max, length=80)

contour!(x_сетка, y_сетка, (x,y)->gmPDF(x,y), levels=20, color=:jet, linewidth=1.5, alpha=0.7, label="", colorbar=false)

for (i, кластер) in enumerate(unique(idx0))
    маска = idx0 .== кластер
    x_класт = X0[маска, 1]
    y_класт = X0[маска, 2]
    форма_маркера = i == 1 ? :dtriangle : :circle
    цвет_маркера = i == 1 ? :blue : :red
    #номер_кластера = "Кластер $i"
    scatter!(x_класт, y_класт, markershape=форма_маркера, markercolor=цвет_маркера, markersize=5, alpha=0.7, label="")#, label=номер_кластера
end

title!("Кластеризация новых данных")
#xlabel!("X₁")
#ylabel!("X₂")
#legend!(position=:best)
display(график5)

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

Заключение

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

Основные выводы:

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

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

Кластеризация данных