python谱聚类算法_从零开始的谱聚类（Spectral Clustering），使用Python实现

机器学习的主要领域之一是无监督学习领域。主要思想是在我们的数据中找到一种模式，而不需要像监督学习那样的标签的先验知识。

它通常通过将我们的数据聚类成组并尝试从聚类中推断出意义来实现。一种比较流行的算法是K均值算法(以及熟悉的EM算法)。在这个算法中，我们在迭代过程中调整K个质心来找到我们的clusters。听起来不错吧？

但主要问题是：

1)它假设数据的形状(圆球，径向基)。

2)有时需要多次重启才能找到局部最小值(即最佳聚类)。

谱聚类算法有助于解决这两个问题。该算法依赖于图的能力和数据点之间的接近度以便对它们进行聚类，从而可以避免K均值算法迫使我们假设的球形聚类。

连接到图表的数据点

该算法的各个阶段是：

1)构造最近邻图(KNN图)或基于半径的图。

2)利用图拉普拉斯算子的特征向量，将数据点嵌入低维空间(谱嵌入)中，其中clusters更明显。

3)使用最低特征值来选择clusters的特征向量。

Python代码如下：

%load_ext watermark

import warnings

warnings.filterwarnings("ignore")

from IPython.core.display import display, HTML

import time

import pandas as pd

import pandas_datareader.data as web

import numpy as np

import matplotlib as mpl

import matplotlib.pyplot as plt

%matplotlib inline

%watermark

from sklearn.datasets import make_moons

random_state = 21

X_mn, y_mn = make_moons(150, noise=.07, random_state=random_state)

cmap = 'viridis'

dot_size=50

fig, ax = plt.subplots(figsize=(9,7))

ax.set_title('Data with ground truth labels - linear separation not possible', fontsize=18, fontweight='demi')

ax.scatter(X_mn[:, 0], X_mn[:, 1],c=y_mn,s=dot_size, cmap=cmap)

好的，我们导入了必要的包并创建了我们的数据点。它看起来像这样：

Data with the ground labels, not a very friendly c

很明显，K均值算法在这里没有帮助，因为群集与算法的假设不一致。

以下是尝试使用K-means解决此聚类问题的示例。

使用kmeans算法的误分类数据样本

让我们尝试遵循谱聚类的各个阶段。

让我们从创建我们需要的最近邻图开始(Python代码实现)：

from sklearn.neighbors import radius_neighbors_graph

from sklearn.neighbors import kneighbors_graph

A = radius_neighbors_graph(X_mn,0.4,mode='distance', metric='minkowski', p=2, metric_params=None, include_self=False)

# A = kneighbors_graph(X_mn, 2, mode='connectivity', metric='minkowski', p=2, metric_params=None, include_self=False)

A = A.toarray()

我们创建了一个邻接矩阵A，它描述了数据点之间的图形。A的定义为：

视觉表现是：

点1和2之间的权重是W12, 1和4之间的权重为0

在我们的例子中(Python代码实现)：

A.shape

(150, 150)

A[:5,:5]

array([[0. , 0. , 0. , 0. , 0.23955829],

[0. , 0. , 0. , 0. , 0. ],

[0. , 0. , 0. , 0. , 0. ],

[0. , 0. , 0. , 0. , 0. ],

[0.23955829, 0. , 0. , 0. , 0. ]])

fig, ax = plt.subplots(figsize=(9,7))

ax.set_title('5 first datapoints', fontsize=18, fontweight='demi')

ax.set_xlim(-1, 2)

ax.set_ylim(-1,1)

ax.scatter(X_mn[:5, 0], X_mn[:5, 1],s=dot_size, cmap=cmap)

for i in range(5):

ax.annotate(i, (X_mn[i,0],X_mn[i,1]))

所以我们在这里有我们的邻接矩阵'A'(150x150)，它描绘了每个点到另一个点之间的距离，我们设置了一个限制，只连接它们之间距离至少为0.4的点，所以所有其他权重之间点将为零。

我们来看看前5个数据点的关系。

矩阵中的对角线为零，因为每个点都没有连接到自身。请注意，它与我们定义的性质也是对称的(即第1点和第2点之间的权重与2和1相同)。

每个数字代表两点之间的权重(距离)，因为你可以看到点0和4非常接近，权重是0.23955829。

但是点0和点3相距很远，因此权重为0(简单来说，点0从点3开始不低于0.4)。

好的，所以我们完成了构建最近邻图(由我们的邻接矩阵描绘)，

让我们继续我们的拉普拉斯算子。

拉普拉斯算子被定义为L = DA，其中A是我们刚看到的我们的成对矩阵，D是称为“度矩阵”的对角矩阵。

在度矩阵中，对角线中的每个单元格是该点的权重之和。

例如，D(1,1)是点1和所有其他点之间所有权重的总和(因此它是邻接矩阵中行号1的总和)。

D(2,2)等同样如此。当我与J不同时，D(I，J)为0。

矩阵D对角线中每个单元的公式

from scipy.sparse import csgraph

L = csgraph.laplacian(A, normed=False)

L[:5,:5]

由A矩阵描绘的图的拉普拉斯矩阵的前五行和列

所以我使用scipy的内置函数进行拉普拉斯算法，因为你可以看到对角线上的拉普拉斯算子值是度矩阵，其余的是来自邻接矩阵的权重，由于我们的公式而带有减号。

接下来，我们将数据点嵌入低维空间(谱嵌入)，其中使用图拉普拉斯的特征向量使聚类更明显。

由于没有深入了解线性代数的细节，通过获得特征值(n个数据点的n个值)和特征向量(n个特征向量大小nx 1)的值，我们可以得到图的最佳分区。

特征值和向量方程

拉普拉斯矩阵的特征值

如果图完全连接，则λ2(特征值编号2)大于0，表示图的代数连通性。特征值(λ2)越大，图形越连接。

eigval, eigvec = np.linalg.eig(L)

np.where(eigval == np.partition(eigval, 1)[1])# the second smallest eigenvalue

output:

(array([1], dtype=int64),)) # the second eigenvalue

因此，对应于第二个最小特征值，为了进行图的Bipartition，我们需要采用相应的特征向量和阈值。特征向量中大于0的每个单元将标记为“1”，而每个低于0的单元将标记为“0”。

y_spec =eigvec[:,1].copy()

y_spec[y_spec < 0] = 0

y_spec[y_spec > 0] = 1

type(y_spec),y_mn.shape,y_spec.shape

output:

(numpy.ndarray, (150,), (150,))

fig, ax = plt.subplots(figsize=(9,7))

ax.set_title('Data with ground truth labels - linear separation not possible', fontsize=18, fontweight='demi')

ax.scatter(X_mn[:, 0], X_mn[:, 1],c=y_spec ,s=dot_size, cmap=cmap)

为了比较，我们可以使用SKlearn中的内置算法。

from sklearn.cluster import SpectralClustering

model = SpectralClustering(n_clusters=2, affinity='nearest_neighbors',

assign_labels='kmeans')

labelsS = model.fit_predict(X_mn)

fig, ax = plt.subplots(figsize=(9,7))

ax.set_title('kernal transform to higher dimension\nlinear separation is possible', fontsize=18, fontweight='demi')

plt.scatter(X_mn[:, 0], X_mn[:, 1], c=labelsS, s=dot_size, cmap=cmap)

与我们刚刚构建的算法结果相同

另一种选择是在我们得到的所有特征向量上使用Kmeans算法。