谱聚类

这篇文章主要介绍的是谱聚类。这是一种目前常见的聚类算法,试图通过无监督的方式将数据点划分为若干类或者簇。谱聚类作为一种聚类算法,具有诸多优点:(1)效果较好,诸多数据集上的对比试验表明其性能优于一般算法(2)使用范围广泛,可以在任意形状的样本空间上聚类且收敛于全局最优解(而传统聚类算法如Kmeans和EM算法都是建立在凸球形的样本空间上的,否则容易陷入局部最优)(3)计算简单:其最终求解是划归成了矩阵特征值的求解,实现简单(4)聚类仅和数据点数目有关,而与维数无关。这个方法在学校的数据挖掘课上老师也专门花了一大节三小节课来介绍,在UW做暑研的老板的一个主要方向就是谱聚类(然而我并不是这个方向的)。在这里也想系统地学习和介绍一下这一算法。本文主要参考的是上课的笔记和讲义。

Kmeans聚类回顾

聚类分析的目的是将诸多数据按照一定的技术非常较为集中的几个类别。我们想要将相似的数据划分为一类。这意味着对于两个样本点来说,我们要用一定的方法对它们是否相似进行定义。一般而言这其实需要特定领域的知识。在不同领域的数据可能具有不同的特征,其算法可能也会需要有针对性地进行调整。

在介绍具体的聚类算法之前我们需要先看一看相似性矩阵(proximity matrices)的概念。这是所有聚类算法的中心概念。一般相似矩阵可以表示成 \(N\times N\)的方阵\(D\). 相似矩阵\(D\)一般是一个对称矩阵,对角线元素为0,其余元素非负. 每个元素\(d_{ij}\)衡量元素\(i\)和元素\(j\)之间的相似程度。

这类数据在有些社会科学的问卷调查中是直接给出的(问卷中可能会要求被试对两者东西的相似度给出度量). 值得指出的是,这类相似矩阵未必符合距离定义(不满足三角不等式:\(d_{ii'} \leq d_{ik} +d_{i'k},\forall i,i',k\in\{1,\cdots,N\}\))有些基于距离假设定义的聚类算法就不能在这样的数据集上使用。

而对于大部分数据集,我们看到的都是原始数据。对于连续变量,一般情况下我们都把不相似性定义成了距离。欧氏距离,曼哈顿距离等等都有使用。对于有序分类变量(例如成绩A,B,C,D,E)其不相似性一般被归一化到[0,1]区间上的间隔为\(\frac{1}{M}\)的等间隔点(\(M\)是该等级变量总共的等级数),然后被当作是一般的连续变量。而分类变量的不相似性一般直接按照是否在同类进行0-1编码。

对于所有的特征的不相似性计算完毕后,如何对样本点的不相似性进行计算呢?可以按照

\[D(x_i,x_{i'})=\sum_{j=1}^p w_j d_j(x_{ij},x_{i'j})\]

进行计算,这里的权重可以主观的选取(在这里权重取相同值并不意味着每个特征对相似度计算的结果有相同程度的影响!)

我们假设数据集中有\(N\)个样本点,我们希望聚成\(K\)类,即希望给每个样本点一个标签\(k \in \{1,\cdots,K\}\), 使得每个样本点都被分配特定的一类 \(k=C(i)\)中. 在这里我们从一个优化的角度来看待Kmeans方法。为了求解结果,我们可以提出一个自然的损失函数,即求解问题

\[\min W(C)=\frac{1}{2}\sum_{k=1}^K\sum_{C(i)=k}\sum_{C(i')=k}d_{ii'}\]

这个优化函数虽然长得不好看,但是其思想是简单的。对于每个聚出来的类,我们计算每个该类中的点到该类中其它所有点的不相似度的和,再计算每个类里的所有不相似度之和。基于类似的思路我们也可以求解问题来最大化类间的不相似性度量

\[\max B(C)=\frac{1}{2}\sum_{k=1}^K\sum_{C(i)=k}\sum_{C(i')\neq k}d_{ii'}\]

当然我们容易发现

\[W(C)+B(C)=\frac{1}{2}\sum_{i=1}^N\sum_{i'=1}^Nd_{ii'}\]

是矩阵\(D\)所有元素的结果的和,是一个定值。因此两个优化问题是等价的。非常明显的问题是这个优化问题是很难解的。因为对几何做分划的情况数随着样本量增长的速度是非常快的。对于\(N\)样本聚\(K\)类的情形,结果有

\[S(N,K)=\frac{1}{K!}\sum_{k=1}^K(-1)^{K-k}\binom{K}{k}k^N\]

种不同的分配方式。对于一个简单的情形,\(S(19,4)\approx 10^{10}\),可见其增长速度之快. Kmeans算法主要适用于所有变量都是连续型,欧氏距离作为不相似性度量的情况。这是它的损失函数可以较为简单地写成

\[\begin{align}W(C)& =\frac{1}{2}\sum_{k=1}^K\sum_{C(i)=k}\sum_{C(i')=k}\parallel x_i-x_{i'} \parallel ^2 \\\\ &=\sum_{k=1}^K N_k \sum_{C(i)=k} \parallel x_i - \overline{x}_k \parallel ^2 \end{align} \]

而这一优化问题就有着一个迭代的解法:

  1. 对每个样本点随机分配一个1到\(K\)的随机数,作为聚类的初始化.
  2. 迭代直到收敛:
    1. 对于每一个类,计算类中心,即每个特征的平均数所在的点;
    2. 把每个样本点分配到和它最近的类中心所在的类别中.

这是解Kmeans问题的Lloyd算法,当然Kmeans聚类也有别的算法,例如Hartigan-Wong,MacQueen等,R语言中kmeans的默认方法是Hartigan-Wong,每次只分配一个点.

以上算法限于欧氏距离,这会导致离群值的权重被放大。事实上我们也可以把它推广到更一般的不相似性度量情形。只要在计算类中心时,选择最小化到该类别所有点的不相似性的点作为类中心,并把这种不相似性度量作为样本点分配的机制。

谱聚类

现在正式开始介绍谱聚类。谱聚类算法是一种基于图论的聚类算法. 首先我们用一个无向带权图\(G=(V,E)\)来表示数据集,每个顶点\(v_i\)表示一个数据点\(x_i\),每条边带一个权。当权为0就表示顶点之间并没有连接。在这样的设定下,我们把聚类问题转化成了 图的划分问题:我们希望找到数据图的一个分划,使得不同类的边之间具有较低的权重,而类内部有较大的权重。一种常用的建图方法是所谓的\(\epsilon-\)近邻图。这类图把所有距离小于等于\(\epsilon\)的点看作是连通的,赋给权重1;而其他点对之间都认为是不连通的。另一个常见的图是\(K\)最近邻图. 当且仅当点\(x_i\)\(x_i'\)\(k\)个最近邻点之一(且/或)反方向也成立时认为它们是连通的。这保证了图的对称性。在\(K\)最近邻图中,我们常用热核来赋给权重,即\(w_{ii'}=\exp(-\parallel x_i-x_{i'} \parallel ^2/\sigma^2)\). 最终我们可以用邻接矩阵\(W\)来表示这个图.

完成建图后,我们定义顶点\(v_i\)的度是 \[ d_i=\sum_{i'=1}^N w_{ii'} \] 并构造一个度矩阵 \[ D=diag(d_1,\cdots,d_N) \] 由此即可得到未归一化的图Laplace矩阵 \[ L=D-W\] 图Laplace算子有很多良好的性质. 我们考虑一个定义在\(\mathbb{R}^N\)上的函数\(f\),它满足 \[ f^TLf = \frac{1}{2}\sum_{i,i'=1}^N w_{ii'}(f_i-f_{i'})^2\] 易见:

  • \(L\)是对称半正定矩阵
  • \(L\)\(N\)个非负实特征值
  • \(L\)的最小特征值是0,对应的特征向量是所有元素都是1的向量.

对于一个顶点集的一个子集\(A\subset V\),我们定义指示向量 \[ \boldsymbol{1}_A = (f_1,\cdots,f_N)^T\] 其中 \(f_i\)用0-1表示顶点\(v_i\)是否在子集\(A\)内. 可以证明:\(L\)的零特征值的重数等于图的连通分支的个数,而零特征值的特征空间是由这些连通分支的指示向量 \(\boldsymbol{1}_{A_1},\cdots,\boldsymbol{1}_{A_k}\)张成的.

证明的思路非常简单:当\(K=1\)时,设\(f\)是零特征值对应的特征向量 \[0=f^TLf=\frac{1}{2}\sum_{i,i'=1}^N w_{ii'}(f_i-f_{i'})^2\] 这说明\(f\)的元素全部相同. 这就说明了所有元素都是1的向量是特征向量.此时由于图连通,指示函数也是常1向量.而当\(K>1\)时,不妨设点是按连通分枝的顺序排列的,此时的\(L\)矩阵事实上可以写成分块对角阵 \[L=\begin{pmatrix}L_1 & & & \\ & L_2 & & \\ & & \ddots & \\ & & & L_K \end{pmatrix} \] 容易验证性质成立.

在利用图进行聚类时,实际情况可能并不是特别理想:连通分枝和聚成的类并不一定完全相同. 但是同类内部的权重之和应当大于外部权重之和. 可以认为它和理想情况(连通分枝刚好和聚的类相同)相差的不是很大。利用摄动理论的想法,我们给\(L\)矩阵一些小的扰动,并可以给出谱聚类算法:

  1. 计算\(L\)矩阵
  2. \(L\)的最小的\(K\)个特征值:\(0= \lambda_1\leq \lambda_2 \cdots \leq \lambda_K\)和相应的特征向量 \(u_1,\cdots,u_K\)
  3. 定义\(z_i=(u_{i1},\cdots,u_{iK})^T,i=1,\cdots N\),事实上我们把原数据集转换到了Laplace算子的特征空间.
  4. \(\{z_i\}\)进行Kmeans聚类,把\(x_i\)聚到对应类别中.

最后我们用一个图像分割的例子来说明谱聚类算法的效果:图(a)是原始图片,(b)-(h)是将它聚成7类以后的效果.可以看到这一算法在图像分割应用中有良好的效果.