PID参数调节的谱多流形聚类算法研究*

罗养霞，马 迪，常言说

1.西安财经大学 信息学院，西安 710100

2.密西根大学 迪尔伯恩分校 计算机与信息科学系，美国 迪尔伯恩 48124

3.计算机应用与商务智能研究中心，西安 710100

1 引言

随着人类社会的发展，从多个数据源得到的多种形态的数据成指数级爆炸性增长，人们逐渐被数据所“淹没”。这些不断涌现的数据表现出数据量大、维数高、结构混杂以及不为人感知等特点。数据的膨胀给数据分析和处理带来了前所未有的困难，而大数据分析算法和处理能力的滞后，使得复杂数据阅读和分析，挖掘和揭示隐藏在复杂表象下的本质特征和内在规律成为挑战。对数据的分析能力没有得到显著提高，也使得对数据的现实应用成为一个瓶颈问题。

流形学习（manifold learning）可以从流形的角度把握数据的内部结构，对离散数据集合的分析，探求嵌入在高维数据空间中低维流形的不同表现形式，寻求事物产生的内在规律性，得到数据所蕴含的与人类认知一致的结构信息。近年来，流形数据研究与分析越来越受到重视。

首先，Seung和Lee[1]在Science上发表的研究认为感知通常是以流形方式存在的，开启了对流形的探究。目前的流形聚类算法大多数是针对表面无相交的情况设计[2-3]，还有一些聚类算法是区别噪声数据在相交的表面具有不同的维度或密度[4-6]。文献[7]研究解决相交曲线上的数据分析。Souvenir等人[8]研究K-means算法，并改进流形曲面数据分析[9-11]。Guo等人[12]提出了使用tabu搜索方式，最小化包含在局部方向的(组合)信息。较好的算法是基于局部主成分分析（principal component analysis，PCA）方法，这方面较成熟的研究[13]是基于细化的多尺度谱方法研究。文献[14]的研究是基于半监督学习环境，其促进了文献[15-16]的工作研究成果。Di Terlizzi等人[17]研究K-manifolds，通过全局度量将局部分解为a个黎曼乘积。文献[18]基于流行表面光滑的假设，关注非参数设置，研究引入曲率约束路径多曲面交叠区域数据分析。文献[19-20]研究基于嵌入式标签改进局部线性嵌入，研究将采样点分片断的多类流形聚类建模和避免噪声。而Jeub等人[21]研究引入模块化思想，来避免参数增加，通过一致化的聚类过程提高聚类算法的一般性应用价值。

2 PID参数调节

对于五类聚类算法，如K-均值聚类、基于模型的方法、基于密度的聚类、高斯混合聚类、层次聚类方法等，这些算法的聚类结果依赖于初始值、数据维数、聚类数目、滑动窗口大小、不同密度阈值点、领域半径或概率分布等。对于参数的调节如果不当，会严重影响聚类算法的适用性，甚至使算法失效。特别是对于多个参数的调节，会降低算法的分类效果和收敛性，部分算法在研究时，考虑非参数设置进行约束研究[18]，有的对数据特征进行假设[15]，如约束为低平滑多流形。目前在聚类算法参数方面研究可分为两类趋势：一类是研究不用输入任何参考的无指导聚类算法[22-23]；另一类研究聚类参数的自动获取或自适应调节[24-25]。

闭环自动控制技术，调节器控制规律为比例-积分-微分控制（proportional-integral-derivative，PID）方式，调节方法被广泛应用于非线性和随时间变化的系统，用于改进神经网络、免疫算法、蚁群优化等研究[26-28]。参数调节的基本原则是由设定合适的比例、积分和微分三个因素作用，对控制结果进行优化和调节，是一种基于反馈的、根据偏差变化趋势而响应的调节方式，可减少结果的不确定性，具有调节速度快，消除余差等优点。PID控制的一般原则如式（1）所示：

根据系统区域指标设定函数，在系统更新迭代时，根据输出效果的精度或适应度，对参数大小或惯性权值进行适应性调整，以达到全局搜索并达到最优解。本文以混合流形算法为基础，研究参数调节，目的使参数调节随结果反馈和优化，最大程度地避免设置参数的随机性和盲目性，保证聚类算法局部和全局最优。

3 参数约束的谱多流形聚类算法

参数约束谱聚类算法的主要思想是：研究基于参数改变的差异函数谱聚类，通过构造相似性矩阵，多个主成分分析器估计局部切空间，并通过参数传递标签PT_label()来标记参数，基于PID调节和控制模型逼近过程，当逼近局部切空间时，按当前最优参数传递作为初始值，以误差和误差变化率为输出，按三维向量调节模型参数，应用ZN法（Ziegler-Nichols）向两边扩展搜索空间，直到获得合理结果。主要步骤包括：（1）构造相似性矩阵；（2）参数计算与求解；（3）参数传递；（4）PID参数调节。

3.1 构造相似性矩阵

相似性矩阵构造基于下述事实：

（1）在局部每个数据点及其邻近点属于一个线性流形，但从全局来看可能是位于或近似属于光滑的非线性流形。

（2）每个数据点的局部切空间提供了非线性流形局部几何结构的优良低维线性近似。

（3）不同流形的相交区域，同一个流形聚类的数据点，有相似的局部切空间；不同流形聚类的数据点，其切空间是不同的。因此可以利用数据点所内含的局部几何结构信息来辅助构造更合适的相似性矩阵（affinity matrix，AM）。

只有当两个数据点满足相互靠近而且具有相似的局部切空间时，才能够断定是来自同一个流形聚类。可举证的两个反例是：（1）数据点和垂直的仿射子空间上的数据点有相似的局部切空间，但它们相互远离；（2）曲面和垂直仿射子空间的相交区域附近的点相互靠近，但它们有不同的局部切空间。

基于以上考虑，在构造相似性矩阵时，既要考虑数据点局部切空间之间的距离相似性LSij（local similarity，LS），又要考虑数据点局部切空间之间的结构相似性SSij̇（structural similarity，SS），这两个相似性融合在一起来决定最后的矩阵相似性：

为了使得构造出的相似性矩阵具有前面分析中所期望的性质，需要构造函数f，使其关于局部切空间相似性LSij单调递减，同时关于两个切空间的结构相似性SSij单调递增。为了细化构造过程，区分流形相叠情况，对两数据点xi和xj的K近邻情况进行区分。关于函数f,LSij和SSij̇详细求解如下所示。

对于数据点xi和xj构造近邻图，Knn（）表示xi或xj的K个近邻数据点，局部相似性LSij定义如下：

两个数据点xi和xj的局部切空间分别为Θi和Θj，则它们之间的结构相似性SSij̇定义如下：

对于式（4）中0≤θ1≤…≤θd≤π/2是两个切空间Θi和Θj之间的主角度，递归定义为：

上述方法不同的是，当含有噪声数据或流形相叠时，调节参数被细化，减小收敛和逼近速度，以得到更好的性能。3.2节将讨论如何计算与求解，递进地逼近局部切空间。

3.2 参数计算与求解

（1）确定x的边际分布函数

其中，uk是数据的均值向量，潜在变量y和噪声εk分别是高斯分布y～N(0,I)和εk～N(0,σ2kI)，在此模型下的联合密度函数为：

其中,X=(x1,x2,…,xk+1)T。

（2）确定对数似然

在逼近训练时，需要确定πk、uk和Σk模型参数，为了防止数据过小而造成浮点数下溢情况，利用EM（expectation maximization）算法最大化观测数据的对数似然来计算。

由于在对数函数里面又有加和，无法直接应用求导解方程来获得最大值。因此利用从GMM（Gaussian mixture model）中随机选点的办法：分成两步逼近，实际上类似于K-means的两步。

（3）迭代计算

估计数据由每个数据样本生成的概率来计算，如式（11）所示。

注意对于每个数据xi来说是生成的概率，而不是被选中的概率。由于要估计uk和Σk的值，采用迭代法，取上一次迭代的值作为初始值。

（4）参数值求解

由于每个数据样本符合标准高斯分布，因此按照下式得到参数值。

通过估计出每个数据点的局部切空间后，计算得到相似性矩阵后，重复迭代前面两步，直到似然函数的值收敛为止。

3.3 基于标签的参数传递

最大对数似然GMM如何能获得全局最优解，若从局部最优达到全局最优，将是最好的建模。但实践中由于GMM和K-means的迭代求解方法都有EM相似，其受初值影响较大，因此也有和K-means同样的问题，即如何保障全局最优。参数值如果调整不好，就有可能得到很差的结果，同时GMM计算工作量要远大于K-means，增加了算法的复杂度。

为了对算法有效控制，增强其在含噪声数据时的鲁棒性，从局部最优收敛达到全局最优。同时，为了减少GMM计算工作量，引入数据点的label（），用参数传递标签PT_label(xt)记录K-means最优参数，然后将其作为初值传递给GMM再进行细致迭代。

传递时将K-means所得的聚类中心（cluster center）传递给GMM函数，GMM所得的结果P(xt)不仅是数据点参数标记PT_label(xt)，同时包含了数据点标记为每个PT_label的概率。再通过估计出每个数据点的局部切空间和计算得到相似性矩阵得到聚类结果，基本过程见算法1。

算法1基于参数传递的谱多流形聚类（spectral multi-manifold clustering_parameter transfer，SMMC_PT）

输入：数据集X，子空间个数N，局部化模型数M,K-邻近个数K。

1.将数据集X划分为i类；

2.训练N个i维的局部线性模型来逼近潜在的流形数据；

3.初始化i=0，当i≤n时，重复如下步骤：

①对于每一点xi，若xi∈PT_label(xt),t∈[1,i],选择k个最邻近的点组成集合Knn(xi)；

②根据式（3）计算任意两点之间的相似度LSij；

③利用式（4）计算两个局部切空间之间的结构相似性SSij；

④根据式（16）确定xi的局部切空间Θi；

⑤PT_label(xt)中将调整后的最优值进行传递；

⑥由以上数据得到调整之后的相似性矩阵AMij；

⑦更新调整过程i=i+1，更新PT_label(xi)为PT_label(xt)；

5.从更新的对角矩阵和相似矩阵计算广义特征矩阵(DM-AM)x=λDMx，并计算前k个特征值对应的特征向量x1,x2,…,xk；

6.利用K-means将 {u1,u2,…,uk}∈ℜN×K的行向量分组为k个聚类。

输出：数据集X对应的聚类结果。

但这仅仅是参考区间范围，针对一般的数据集而言，需要进一步改进参数值随精度、误差的优化调整。

3.4 基于PID参数调节

由于多流形上存在噪声时，数据点并不是精确地位于某一个具体流形上，会影响聚类划分和流形学习，为了减少算法对噪声的敏感性，达到更快收敛，减少参数调节的盲目性，应用基于反馈的PID参数调节。该方法的调节参数基本思想是，按当前参数传递作为初始值，以误差和误差变化率为输出，按三维向量调节模型参数，应用ZN法向两边扩展搜索空间，直到获得合理结果，变换后的调节函数如下计算。

其中，u(t)=u(t-1)+Δu(t)，而 Δu(t)=K(e(t)-e(t-1))+Me(t)+γ[e(t)-2e(t-1)+e(t-2)]。三个参数组成自适应优化函数X(K,M,γ)，构成一个三维行向量，应用ZN法获得参数，并向两边扩展搜索空间。

式中，K∗、M∗、γ∗为整定值，α、β为扩展系数，控制原则如下：

（1）当误差较大时，为了使聚类过程有更好的全局跟踪性能，取较大的K与较小的M，同时为了避免系统响应超调，通常取γ=8。

（2）当误差和误差变化率中等大小时，K应取得小些，M的取值受系统的影响较大，就取得小一些，γ的取值应适当。

（3）当误差较小时，为了使系统具有较好的稳定性，K和M均应取得大些，为了防止系统受干扰振荡，当误差变化率较大时，M可取得小一些；当误差变化率较小时，M可取大一些。

具体过程如算法2所示。

算法2基于PID的参数调节方法（spectral multimanifold clustering_proportional-integral-derivative，SMMC_PID）

输入：原始数据集X

3.输出聚类错误、错误率和聚类精度；

4.If在PID参数调节控制中出现合理的聚类精度，聚类误差以及变化率在阈值范围内，输出参数和聚类结果。

5.Else ZN法扩展搜索空间，重新设定K，进行下一次M的选择，直到找到合理的结果，输出设定的参数。

6.End if

7.While没有找到合理的分类结果调节γ，重复上面的过程。

输出:原始数据的聚类结果。

4 分析与实验比较

4.1 算法复杂度分析

SMMC_PID算法复杂度分为计算相似性矩阵AMij，参数计算和求解，以及参数调节三个主要过程。假设数据维数为D，近邻点数为K，有N个切空间，通过M个局部分析器进行求解。

（1）计算相似性矩阵

对于任意两个数据点的局部切空间分别为Θi和Θj，在递归求解局部切空间相似性时，复杂度为O(N2Dd2)；搜索每个数据点的近邻点的复杂度为O((D+K)N2)。对AMij利用谱方法将数据投影到k维嵌入空间，并在该空间中利用K-means将数据分组为k个聚类，其中广义特征分析的复杂度为O((N+k)N2),则复杂度之和为O(N3+(Dd2+D+K+k)N2)。

（2）参数计算和求解过程

由于计算过程包括确定x的边际分布函数，确定对数似然，迭代计算和参数传递，因此N个局部切空间，应用M个分析器进行学习时，应用PT_label（t）将模型参数传递来优化初始化值，这个过程的复杂度为O(N(n+dm)MD)，其中n和m分别为K-means和EM过程收敛所需的迭代步数，这里由于PT_label参数传递使t≪n。而K-means在k维投影数据上的复杂度为O(Nk2q)（q为该过程中K-means收敛所需的迭代步数），则此过程的复杂度合计为O(N((n+dm)MD+k2q))。

（3）PID参数调节复杂度

在多次蒙特卡洛重复过程和标记，涉及对K、M的调节，这个过程的复杂度为O(N2+N(M+K)+t)，其中t为迭代次数。

因此，算法SMMC_PID的复杂度为三者之和，又由于K-means和EM过程收敛所需的迭代步数通常较低，并且d＜D，K≪N，k≪N，因此这三者之和可以化简为O(N3+N2D+ND)，因此SMMC_PID复杂度主要由数据点数N和数据维数D来决定。

4.2 实验结果对比

根据上面叙述的规则和算法过程，对不同数据类型进行实验和对比，如图1所示。显示了应用PID参数调节前（左图）、后（右图）谱多流形聚类分类结果性能，从多个实验结果可以看出，基于PID反馈式聚类参数约束机制能显著提高聚类性能。

4.3 聚类精确性及复杂度对比

从精确性和复杂度方面对算法进行比较，这里聚类精确性（clustering accuracy）指所有可能分类中的最大分类精度，度量标准依据以下公式：

Fig.1 Comparison of optimal clustering performance of different datasets图1 不同数据集最优聚类性能对比

其中，ti指正确标签，ci是xi的聚类标签，δ为脉冲函数。在考虑错误时包括N个局部线性分析器逼近潜在流形的重构误差。复杂度（complexity）按平均运行时间（以秒计）来进行算法的复杂度估计。

对比算法的选择取决于算法特点、代码的可用性，以及与此算法的相关性。由于算法GPCA（generalized principal component analysis）、K-planes、LSA（latent semantic analysis）、SCC（spectral curvature clustering）是几种典型的适用于简单、线性流形聚类方法，与本文提出的算法相关性差异较大。本文选择对比算法时，选择以下几种非线性流形聚类算法进行比较：K-means[9]、SC算法[2]、K-manifolds算法[17]、SMMC[15]。分别在图1（a）～（e）合成数据上进行实验，计算均值和标准差、最高精度，以及运行的平均时间，结果如表1所示。

实验中发现，SMMC_PID算法与同类算法相比，在聚类精度上从可视化图形和数值都具有明显优势。结果反馈和参数传递约束聚类算法，减少参数设置的盲目性，缩短了系统迭代过程，但在运行平均时间消耗方面，优势不突出，分析是由于增加了参数传递和部分计算时间的原因。

4.4 真实数据中的实验分析

为了验证算法在真实数据中的可用性，选择四类数据集分别为：（a）Caltech101数据集的一个子集；（b）Yale Face数据集 B 的一个子集；（c）2-维图像COIL-20数据集；（d）3-维运动分割数据库Hopking155数据集（http://www.vision.caltech.edu/Image_datasets/Caltech101,http://vision.ucsd.edu/～leekc/ExtYale-Database/ExtYaleB.html,http://www.cs.columbia.edu/labs/cave/software/softlib/coil-20,http://www.vision.jhu.edu/data/hopkins155)。

在以上数据集进行重复实验，获得平均聚类精度、平均时间如图2、图3所示。

在此实验统计中可以看出，SMMC_PID分类的平均聚类性能优于其他同类算法。聚类的复杂度，除K-manifolds算法总体较高，其他SC、SMMC和SMMC_PID总体接近且偏低。分析认为SMMC_PID算法在参数控制方面增加了计算量和迭代时间，但由于增加了总体聚类收敛速度，减少系统计算复杂度。

Table 1 Clustering performance comparison of different algorithms on different datasets(mean±std.,highest accuracy,computation time)表1 不同聚类算法在不同数据集聚类性能比较（精度均值±标准差，最佳性能，运行时间）

Fig.2 Comparison of average accuracy of different algorithms图2 不同算法的平均聚类精度实验对比

Fig.3 Comparison of average time of different algorithms图3 不同算法的平均时间实验对比

5 总结

参数调节在聚类中是一个难点问题，直接影响着聚类性能。本文研究混合流形聚类参数调节问题，通过构造相似性矩阵，多个主成分分析器估计局部切空间，并通过参数传递和PID调节，控制模型逼近过程，按三维向量调节模型参数，应用ZN法扩展搜索空间，使聚类过程与结果进行反馈，改进聚类算法的精确性和复杂度。经在合成数据和真实数据不同类型数据特征集进行检测，可获得较好的聚类性能。

下一步将整合项目研究工作，结合动态嵌入标签方式[20]，使谱聚类能有效应用于较复杂的实践环境分析。