TNNLS-2024《Manifold Peaks Nonnegative Matrix Factorization (MPNMF)》

2. 核心思想

这篇论文的核心思想是解决传统非负矩阵分解（NMF）在聚类任务中一个关键问题：分解得到的聚类中心（centroids）往往会偏离数据的真实流形（manifold）结构，尤其是在处理具有复杂几何结构的数据集时，这会严重影响聚类性能。

为了解决这个问题，MPNMF提出了以下核心策略：

• 引入流形峰（Manifold Peaks, MPs）: 首先识别出数据流形上的关键点——流形峰（MPs）。这些点被认为是数据流形的“骨干”（backbone），能够较好地表征流形的主要结构，同时排除了噪声和离群点的干扰。
• 基于MPs构建聚类中心: 不是让NMF直接学习聚类中心，而是强制要求每个聚类中心必须是附近流形峰（nearby MPs）的锥组合（conic combination）。锥组合意味着组合系数是非负的，但不要求像凸组合那样和为1。这使得聚类中心能够被约束在由MPs张成的空间内，从而更有可能落在原始数据流形上。
• 利用广义图平滑正则化（Generalized Graph Smoothness Regularization）: 引入一种新的正则化方法来指导图的构建。具体来说，它同时考虑了MPs之间的流形距离和它们在构建聚类中心时的权重分布相似性，确保邻近的MPs在构建不同聚类中心时具有相似的权重分布，进一步促使聚类中心位于流形上。
• 稀疏性约束: 对每个聚类中心仅由少量（附近的）MPs组合而成施加稀疏性约束，增强了聚类中心构建的局部性。

3. 目标函数

MPNMF的目标函数（公式15）旨在最小化数据重构误差，同时引入上述约束：

$$\underset{Z,V,Y}{\min }\frac{1}{2}\Vert X-PZ{V}^T{\Vert }_F^2+\frac{\alpha }{2}⟨Q_P,Z{Z}^T⟩+\frac{\beta }{2}⟨Q_Y,V{V}^T⟩$$

$$\text{s.t. }Z\ge 0,V\ge 0,\forall j\Vert z_j{\Vert }_0\le {\tau }_z,\forall i\Vert v_i{\Vert }_0\le {\tau }_v,Y^{T}Y=I$$

其中：

• $X\in {\mathbb{R}}_{+}^{m\times n}$ 是输入的非负数据矩阵（$m$个特征，$n$个样本）。
• $P\in {\mathbb{R}}^{m\times p}$ 是选出的$p$个流形峰（MPs）组成的矩阵。
• $Z\in {\mathbb{R}}_{+}^{p\times r}$ 是MPs的锥组合权重矩阵，$z_j$是第$j$个聚类中心（列）的组合权重向量，$r$是聚类中心的数量。
• $V\in {\mathbb{R}}_{+}^{n\times r}$ 是系数矩阵，$v_i$是第$i$个数据样本（行）与聚类中心的关联向量。
• $Y\in {\mathbb{R}}^{n\times c}$ 是最终的聚类指示矩阵（$c$个簇），满足正交性约束 $Y^{T}Y=I$。
• $\Vert \cdot {\Vert }_F$ 是Frobenius范数。
• $⟨A,B⟩=\text{tr}(A^{T}B)$ 是矩阵内积。
• $Q_P=[\frac{1}{2}(d_g^{(p)}(i,j){)}^2{]}_{p\times p}$，其中 $d_g^{(p)}(i,j)$ 是第$i$个和第$j$个MPs之间的测地距离。这一项 $⟨Q_P,Z{Z}^T⟩$ 是作用于MPs的广义图平滑正则化，鼓励邻近MPs具有相似的组合权重。
• $Q_Y=[\frac{1}{2}\Vert y_i-y_j{\Vert }^2{]}_{n\times n}$，其中 $y_i$ 是第$i$个样本的聚类指示向量。这一项 $⟨Q_Y,V{V}^T⟩$ 是作用于数据样本的广义图平滑正则化，鼓励关联相似聚类中心的样本属于相似的簇。
• $\alpha ,\beta$ 是正则化参数，控制正则化项的强度。
• ${\tau }_z,{\tau }_v$ 是稀疏性参数，分别限制每个聚类中心使用的MPs数量和每个数据样本关联的聚类中心数量。
• $Z\ge 0,V\ge 0$ 是非负约束。

4. 目标函数的优化过程

该问题是一个非凸优化问题。论文通过交替迭代的方式优化变量 $Z$, $V$, 和 $Y$。

• 更新 $Z$ (固定 $V,Y$):

  • 优化问题（公式25）：
    $$\underset{Z}{\min }J(Z)=\frac{1}{2}\Vert X-PZ{V}^T{\Vert }_F^2+\frac{\alpha }{2}⟨Q_P,Z{Z}^T⟩\text{s.t. }Z\ge 0,\forall j\Vert z_j{\Vert }_0\le {\tau }_z$$
  • 梯度计算：$\nabla J(Z)=-P^{T}XV+P^{T}PZ{V}^{T}V+\alpha Q_{P}Z$
  • 优化方法：将矩阵问题向量化，转化为一个带有组$l_0$范数约束的二次正则化非负最小二乘（Quadratic Regularized Group $l_0$-NNLS）问题（公式27）。
  • 解决方案：论文推导了一个通用的求解此类问题的迭代近端梯度定理（Theorem 2），并将其应用于更新 $Z$（见Algorithm 2）。具体更新步骤为：
    $$Z\leftarrow P_{+,{\tau }_z}[Z-\eta \nabla J(Z)]$$
    其中 $\eta$ 是步长，$P_{+,{\tau }_z}[\cdot ]$ 是一个投影算子，它对矩阵的每一列独立操作，将每个列向量投影到非负且 $l_0$ 范数不超过 ${\tau }_z$ 的空间（即保留最大的 ${\tau }_z$ 个正值，其余置零，并保持非负性）。

• 更新 $V$ (固定 $Z,Y$):

  • 优化问题（公式29）：
    $$\underset{V}{\min }O(V)=\frac{1}{2}\Vert X-PZ{V}^T{\Vert }_F^2+\frac{\beta }{2}⟨Q_Y,V{V}^T⟩\text{s.t. }V\ge 0,\forall i\Vert v_i{\Vert }_0\le {\tau }_v$$
  • 梯度计算：$\nabla O(V)=-X^{T}PZ+V{Z}^T{P}^{T}PZ+\beta Q_{Y}V$
  • 优化方法：同样向量化并转化为Group $l_0$-NNLS问题（公式30）。
  • 解决方案：应用Theorem 2进行更新（见Algorithm 3）。具体更新步骤为：
    $$V\leftarrow P_{+,{\tau }_v}[(V-\eta \nabla O(V){)}^T{]}^T$$
    （注意这里因为向量化的是 $V^T$，所以更新后需要转置回来）。

• 更新 $Y$ (固定 $Z,V$):

  • 优化问题：
    $$\underset{Y^{T}Y=I}{\min }⟨Q_Y,V{V}^T⟩=\frac{1}{2}\sum _{i,j}⟨v_i,v_j⟩\Vert y_i-y_j{\Vert }^2$$
  • 优化方法：这个问题等价于最小化 $\text{tr}(Y^T{L}_{V}Y)$，其中 $L_V=\text{diag}(W_V{1}_n)-W_V$ 是基于当前 $V$ 构建的图 $G_V$ 的拉普拉斯矩阵 ($W_V=V{V}^T$)。
  • 解决方案：标准的谱聚类方法，取 $L_V$ 的前 $c$ 个最小特征值对应的特征向量作为新的 $Y$（见Algorithm 4第8步）。
  • 更新 $Q_Y$：使用更新后的 $Y$ 按照公式(13) $Q_Y=[\frac{1}{2}\Vert y_i-y_j{\Vert }^2{]}_{n\times n}$ 计算新的 $Q_Y$。

整个优化过程通过交替执行上述三个步骤，直到收敛。

5. 主要贡献点

1. 提出流形峰（MPs）概念： 创新性地引入了MPs来表征数据流形的骨干，作为构建聚类中心的基础，有效避免了噪声和离群点的干扰。
2. 基于MPs的锥组合约束： 强制聚类中心是附近MPs的锥组合，从根本上解决了传统NMF聚类中心偏离数据流形的问题。
3. 广义图平滑正则化： 提出了广义图平滑正则化框架，能够指导图的构建，并将MPs的几何距离与其权重分布联系起来，增强了聚类中心在流形上的定位。
4. 解决Group $l_0$-NNLS问题： 理论上分析并提供了一种高效的算法来求解目标函数中出现的带组$l_0$范数约束的二次正则化非负最小二乘问题。
5. 显著提升聚类性能： 通过在合成和真实世界数据集上的大量实验证明，MPNMF相比多种先进的NMF聚类算法，在NMI和F-measure等指标上取得了显著的性能提升。

6. 算法实现过程 (Algorithm 4)

1. 输入: 数据矩阵 $X$，流形峰数量 $p$，聚类中心数量 $r$，簇数量 $c$，正则化参数 $\alpha ,\beta$，稀疏性参数 ${\tau }_z,{\tau }_v$。
2. 选择流形峰 (MPs):
   • 调用 MP(X, p) 算法（Algorithm 1）：
     • 计算所有样本点对之间的测地距离矩阵 $D_g$ (使用Dijkstra算法)。
     • 根据公式(5)计算每个样本的局部测地密度 ${\rho }_i^g$。
     • 根据公式(6)计算每个样本到更高密度样本的最小测地距离 ${\delta }_i^g$ (对于最高密度点，取其到所有点的最大距离)。
     • 根据公式(7)计算综合指标 ${\gamma }_i^g=\frac{{\rho }_i^g}{\max ({\rho }^g)}\times \frac{{\delta }_i^g}{\max ({\delta }^g)}$ 进行归一化。
     • 选择 ${\gamma }_i^g$ 值最大的前 $p$ 个样本作为流形峰，得到矩阵 $P$ 和MPs之间的测地距离矩阵 $D_g^p$。
3. 初始化:
   • 初始化迭代次数 $t=0$。
   • 初始化 $Z^0,V^0,Y^0,Q_Y^0$：
     • $Z^0$：将前 $r$ 个MPs作为初始中心，每个中心的权重设为最近 ${\tau }_z$ 个MPs的组合（权重为 $1/{\tau }_z$ 加上小扰动）。
     • $V^0$：根据每个样本与最近 ${\tau }_v$ 个初始中心的相似性初始化。
     • $Y^0$：使用k-means对 $X$ 进行聚类，结果用one-hot编码表示。
     • $Q_Y^0$：根据 $Y^0$ 和公式(13)计算。
4. 迭代优化:
   • 循环开始 (直到收敛):
     • 更新 $Z^{t+1}$: 使用Algorithm 2（基于近端梯度法求解Group $l_0$-NNLS）更新 $Z$。
     • 更新 $V^{t+1}$: 使用Algorithm 3（基于近端梯度法求解Group $l_0$-NNLS）更新 $V$。
     • 计算中间变量 $W_V^{t+1}$: $W_V^{t+1}=V^{t+1}(V^{t+1}{)}^T$。
     • 计算拉普拉斯矩阵 $L_V^{t+1}$: $L_V^{t+1}=\text{diag}(W_V^{t+1}1)-W_V^{t+1}$。
     • 更新 $Y^{t+1}$: 对 $L_V^{t+1}$ 进行特征分解，取前 $c$ 个最小特征值对应的特征向量作为新的 $Y^{t+1}$。
     • 更新 $Q_Y^{t+1}$: 根据更新后的 $Y^{t+1}$ 和公式(13)计算新的 $Q_Y^{t+1}$。
     • 迭代计数: $t=t+1$。
   • 循环结束条件: 检查目标函数值或变量更新幅度是否满足收敛条件。
5. 输出: 最终的 $Z^t,V^t,Y^t$。$Y^t$ 即为最终的聚类结果。

通过以上步骤，MPNMF实现了在复杂流形结构数据上的高质量聚类。

好的，我们来详细介绍一下流形峰（Manifold Peaks, MPs）的实现过程和背后的数学原理。

正如论文所述，流形峰的核心目的是识别出数据流形上的关键点，这些点能够表征流形的“骨干”（backbone），排除噪声和离群点的干扰，从而为后续的聚类中心构建提供一个高质量、有代表性的基础。

实现过程 (Algorithm 1)

流形峰的识别过程主要在论文的Algorithm 1中描述，具体步骤如下：

1. 输入: 原始数据矩阵 $X$（维度为 $m\times n$，$m$是特征数，$n$是样本数）和期望选出的流形峰数量 $p$。
2. 计算全局测地距离矩阵 ($D^g$):
   • 使用Dijkstra算法（或类似算法）计算数据集中任意两点 $x_i$ 和 $x_j$ 之间的测地距离（Geodesic Distance） $d_{ij}^g$。
   • 测地距离是指数据点在流形上的最短路径距离，而不是简单的欧氏距离。这能更好地反映数据在低维流形结构上的真实距离。
   • 输出一个 $n\times n$ 的测地距离矩阵 $D^g=[d_{ij}^g]$。
3. 计算局部测地密度 (${\rho }_i^g$):
   • 对于每个数据点 $x_i$，根据其与其他所有点的测地距离 $d_{ij}^g$ 计算其局部测地密度 ${\rho }_i^g$。
   • 计算公式（论文公式5）：
     $${\rho }_i^g=\sum _j{e}^{-(d_{ij}^g/d_c^g{)}^2}$$
   • 其中，$d_c^g$ 是一个截断距离（cutoff distance）。只有当 $d_{ij}^g$ 小于 $d_c^g$ 时，点 $x_j$ 对点 $x_i$ 的密度贡献才比较显著。$d_c^g$ 通常根据经验设定，例如取所有测地距离中较小的百分位数（如1%-2%）。
   • 直观理解：密度高的点周围（在测地距离意义上）聚集了更多的邻居点。
4. 计算最小距离 (${\delta }_i^g$):
   • 对于每个数据点 $x_i$，计算它到任何密度比它高的点的最小测地距离 ${\delta }_i^g$。
   • 计算公式（论文公式6）：
     $${\delta }_i^g=\{\begin{array}{ll}{\min }_{j:{\rho }_j^g>{\rho }_i^g}{d}_{ij}^g,& \text{如果存在密度比 }{x}_i\text{ 高的点}\\ {\max }_j{d}_{ij}^g,& \text{如果 }{x}_i\text{ 是密度最高的点}\end{array}$$
   • 直观理解：对于普通点，${\delta }_i^g$ 衡量它到更密集区域的远近；对于密度最高的点，用它到最远点的距离 ${\delta }_i^g$ 来衡量其“孤立性”。
5. 计算综合指标 (${\gamma }_i^g$):
   • 结合密度 ${\rho }_i^g$ 和距离 ${\delta }_i^g$，为每个点计算一个无量纲的综合指标 ${\gamma }_i^g$。
   • 计算公式（论文公式7）：
     $${\gamma }_i^g=\frac{{\rho }_i^g}{\max ({\rho }^g)}\times \frac{{\delta }_i^g}{\max ({\delta }^g)}$$
   • 通过将 ${\rho }_i^g$ 和 ${\delta }_i^g$ 分别除以各自的最大值进行归一化，消除了它们因量纲或数值范围不同带来的影响，使得指标更具鲁棒性。
   • 直观理解：${\gamma }_i^g$ 值高的点，要么本身密度高（可能是簇的中心），要么距离高密度区域很远（可能是两个簇之间的桥接点或边界点）。这两种点都可能是流形上的关键点（峰值）。
6. 选择流形峰:
   • 根据计算出的 ${\gamma }_i^g$ 值对所有数据点进行排序。
   • 选择 ${\gamma }_i^g$ 值最高的前 $p$ 个点作为流形峰（MPs）。
   • 输出：选出的 $p$ 个流形峰组成的矩阵 $P$（维度为 $m\times p$）以及这些流形峰之间的测地距离矩阵 $D_g^p$（维度为 $p\times p$）。

数学原理

流形峰算法的数学原理主要基于以下几点：

1. 测地距离 ($d_{ij}^g$): 这是核心。传统的欧氏距离在高维空间中，尤其是在数据实际分布在一个低维流形上时，往往不能准确反映点与点之间的真实相似性。测地距离通过在数据点构成的近邻图上寻找最短路径，能够近似地反映数据点在潜在流形上的距离，更好地捕捉流形的几何结构。
2. 密度 (${\rho }_i^g$) 的定义: 使用高斯核函数 $e^{-(d_{ij}^g/d_c^g{)}^2}$ 来衡量点 $x_j$ 对点 $x_i$ 密度的贡献。这使得距离 $x_i$ 较近的点贡献更大，距离较远的点贡献迅速衰减至0（当 $d_{ij}^g>d_c^g$ 时）。这种定义方式使得密度能够反映点在流形局部区域的集中程度。
3. 距离 (${\delta }_i^g$) 的定义: 通过计算到更高密度点的最小距离，${\delta }_i^g$ 区分了不同类型的点：
   • 高密度且大 ${\delta }_i^g$: 这些点通常是簇的中心（Cluster Centers），因为它们局部密度高且远离其他高密度区域。
   • 低密度且大 ${\delta }_i^g$: 这些点可能是簇间的桥接点或噪声点。
   • 高密度且小 ${\delta }_i^g$: 这些点是簇中心附近的普通点。
   • 低密度且小 ${\delta }_i^g$: 这些点是簇内部或边缘的普通点或噪声点。
4. 综合指标 (${\gamma }_i^g$) 的设计: ${\rho }_i^g\times {\delta }_i^g$（或其归一化形式）的乘积形式是一个关键设计。它有效地结合了密度和距离信息，使得同时具有高密度和大距离（或相对较大的距离）的点能够脱颖而出。这正是簇中心点（流形峰）的典型特征：它们在局部区域内密度高，并且与其他高密度区域（可能是其他簇中心）保持一定的距离。这种乘积形式增强了算法识别真正“峰值”点的能力。
5. 与传统密度峰值聚类（DP）的对比: 论文指出，传统的DP算法主要在欧氏空间中计算密度和距离，容易受到复杂流形结构和噪声的影响。MP算法通过将这些计算转移到测地距离空间，并使用归一化的综合指标，显著提高了在复杂流形数据上识别关键点（骨干点）的鲁棒性和准确性（如论文中的Fig. 3所示）。

总结来说，流形峰算法通过在测地距离空间中计算密度和距离，并利用一个乘积形式的无量纲指标来识别数据流形上的关键点，这些点构成了数据流形的骨干，为后续的、基于这些点构建聚类中心的MPNMF算法提供了坚实的基础。