在图像处理和机器学习的应用中，经常须处理高维数据，这些数据存在于一个低维结构中．低维结构可以用高维空间中的线性低维子空间建模．子空间聚类是指对一组来自多个线性子空间的数据进行聚类．现有算法可分为迭代算法[1-2]、代数算法[3]、统计算法[4]和谱聚类算法[5-6]．谱聚类算法首先在数据样本间建立一个相似矩阵，然后对相似矩阵进行谱聚类[7]得到数据聚类[6]．根据构建相似矩阵时使用信息的不同，谱聚类算法可分为局部谱聚类算法和全局谱聚类算法．局部谱聚类算法，如局部子空间相似性(local subspace affinity，LSA)[8]、局部线性流形聚类(locally linear manifold clustering，LLMC)[9]和局部最优谱平面(spectral local best-fit flats，SLBF)[10]，使用每个数据点周围的局部信息构建相似矩阵[6]．全局子空间聚类算法采用全局信息构建更加合理的相似矩阵．谱曲率聚类(spectral curvature clustering，SCC)[11]使用极曲率的概念定义同一子空间中数据点之间的相似性[6]．近年来在稀疏和低秩恢复算法方面出现了一些重要的全局子空间聚类算法，如稀疏子空间聚类(sparse subspace clustering，SSC)[5]和低秩表示(low-rank representation，LRR)[12]．这些算法通过找到一个用数据自身表示数据的表示矩阵，然后构建一个凸优化问题，并求解表示矩阵[6]，最后根据表示矩阵构建相似矩阵，并得到数据的聚类．在低秩表示算法中，须采用一个矩阵范数近似矩阵秩函数．核范数和F范数都被用来近似秩函数，在某些应用中，基于F范数的算法具有更好的性能．为了扩展核范数和F范数的概念，本研究提出采用加权核范数(weighted nuclear norm，WNN)近似矩阵秩函数，以得到更优的聚类结果．1 低秩表示模型受文献[13]的启发，使用加权核范数(weighted nuclear norm，WNN)约束低秩表示模型中的表示矩阵Z．基于加权核范数最小化(weighted nuclear norm minimization，WNNM)的低秩表示模型可以写为minZ,E||Z||w,*+λ||E||2,1;s.t.   [X=XZ+E], (1)式中：||∙||w,*为加权核范数；Z为表示矩阵；w=[w1,w2,…,wn]，wi≥0为施加于矩阵Z奇异值的权重，加入权重可以提升低秩表示模型的灵活性；λ为权重参数；X为数据矩阵；E为噪声，将Z的核范数约束替换为WNN约束．式(1)称为WNNM-LRR模型．类似于LRR模型，使用ADMM求解WNNM-LRR问题，可将式(1)写成如下的等价问题：min||J||w,*+λ||E||2,1;s.t.    [X=XZ+E], [Z=J], (2)式中J为辅助变量．式(2)可以通过最小化增广拉格朗日函数求解    L=||J||w,*+λ||E||2,1+tr(Y1T(X-XZ-E))+tr(Y2T(Z-J))+μ2(||X-XZ-E||F2+||Z-J||F2),式中：Y1和Y2为拉格朗日乘子；μ0为惩罚因子．上述问题可以通过固定其他变量，分别对J，Z和E最小化，然后更新拉格朗日乘子Y1和Y2实现．求解过程见算法1．算法1 WNNM-LRR输入 数据矩阵X，参数λ，权重向量w初始化 Z(0)=J(0)=0，E(0)=0，Y1,(0)=0，Y2,(0)=0，μ(0)=1×10-6，μmax=1×106，ρ=1.1，ε=1×10-8，k=0while not converged do固定其他变量，更新J，J(k+1)=argmin            J1μ(k)||J||w,*+12||J-(Z(k)+Y2,(k)/μ(k))||F2.固定其他变量，更新Z，    Z(k+1)=(I+XTX)-1(XT(X-E(k))+J(k+1)+(XTY1,(k)-Y2,(k))/μ(k)).固定其他变量，更新E，    E(k+1)=argmin            Eλμ(k)||E||2,1+12||E-(X-XZ(k+1)+Y(k)/μ(k))||F2.更新拉格朗日乘子：Y1,(k+1)=Y1,(k)+μ(k)(X-XZ(k+1)-E(k+1))；Y2,(k+1)=Y2,(k)+μ(k)(Z(k+1)-J(k+1))．更新惩罚因子μ，μ(k+1)=min(ρμ(k),μmax)，检查是否满足收敛条件：||X-XZ(k+1)-E(k+1)||∞ε；||Z(k+1)-J(k+1)||ε；k=k+1．end while输出 矩阵Z=Z(k+1)，E=E(k+1)．算法1中，ρ为惩罚因子相邻迭代时的放大倍数；k为迭代次数；ε为小常数，用于检验收敛条件是否满足．在算法1的第1步中，须求解WNNM问题．使用权重按照降序排列，权重wi=σiγ(X)，式中γ为参数，决定了权重的分布．采用奇异值的某个指数作为权重，希望扩展核范数与F范数的概念．新的范数将是核范数与F范数的自然扩展．如果γ=0，即所有权重取1，那么加权核范数变为核范数；如果γ=1，即权重等于奇异值本身，那么加权核范数变为F范数；如果0γ1，那么可以得到核范数与F范数之间的矩阵范数．根据基于F范数低秩表示的实验结果可以发现，该算法的聚类性能比LRR高．使用一个更加灵活的矩阵范数来约束表示矩阵，通过取合适的γ值，可以期望得到更高的聚类性能．因为奇异值按照降序排列，所以权重也按照降序排列．此时第1步中的优化问题是一个凸优化问题，该问题有闭式解，可以写为J(k+1)=US1/μ(k)wΣVT，式中UΣVT是矩阵Z(k)+Y2,(k)/μ(k)的奇异值分解．2 LADMM求解模型当使用ADMM求解低秩表示模型时，须将表示矩阵Z替换为矩阵J，并引入一个等式约束．若直接通过ADMM求解，则增广拉格朗日函数可以写为    ℒ'(Z,E,μ)=||Z||*+λ||E||2,1+tr(Y'(X-XZ-E))+μ||X-XZ-E||F2/2.然后将ℒ'的最小化分解为对(Z，E)最小化的2个子问题．ADMM的迭代过程为    Z(k+1)=argmin            Zℒ'(Z,E(k),μ'(k))=argmin           Z||Z||*+μ'(k)||X-XZ-E(k)+Y(k)/μ'(k)||F2/2; (3)               E(k+1)=argmin             Eℒ'(Z(k+1),E,μ'(k))=           argmin             Eλ||E||2,1+μ'(k)||X-XZ(k+1)-            E+Y(k)/μ'(k)||F2/2; (4)Y'(k+1)=Y'(k)+μ'(k)(X-XZ(k+1)-E(k+1))，式中：Y'为拉格朗日乘子；μ'为惩罚因子．式(4)有闭式解，但是式(3)没有闭式解，因为E前面的线性变换是单位矩阵，而Z前面的线性变换是X，不是单位矩阵，此时必须迭代求解式(3)．为了克服这一困难，通常的做法是引入辅助变量J，然后变为式(4)中的等价形式．文献[14]引入了LADMM来求解待优化变量前不是单位矩阵的子问题．对式(3)中的二次项在Z(k)处线性化，并加入了邻近项．因为无须增加辅助变量，所以LADMM可以改善ADMM的收敛性能．LADMM[14]是一种通用方法，可以用于求解其他凸问题．所提的WNNM-LRR模型也可以通过LADMM求解．相应的优化过程见算法2．算法2 通过LADMM求解WNNM-LRR输入 数据矩阵X，参数λ，权重向量w初始化 Z(0)=0，E(0)=0，μ'maxμ'(0)0，ρ01，ε10，ε20，ηXσmax2(X)，k=0while not converged do固定其他变量，更新E    E(k+1)=argmin            Eλμ'(k)||E||2,1+12||E-(X-XZ(k)+Y'(k)/μ'(k))||F2.固定其他变量，更新Z，    Z(k+1)=argmin              Z||Z||w,*+μ'(k)ηX2||Z-Z(k)-XT(Y(k)+μ'(k)(X-XZ(k)-E(k+1)))/(μ'(k)ηX)||F2.更新拉格朗日乘子，Y'(k+1)=Y'(k)+μ'(k)(X-XZ(k+1)-E(k+1))．更新惩罚因子μ'，μ'(k+1)=min(ρ'μ'(k),μ'max)，ρ'=    ρ0   (μ'(k)max(||Z(k+1)-Z(k)||F||E(k+1)-E(k)||F)/||X||Fε2);    1   (其他),检查是否满足收敛条件：||X-XZ(k+1)-E(k+1)||∞ε1；max(||Z(k+1)-Z(k)||F;||E(k+1)-E(k)||F)/||X||Fε2；k=k+1．end while输出 矩阵Z=Z(k+1)，E=E(k+1)．通过LADMM求解的WNNM-LRR称为WNNM-LRR(L)．基于WNNM的低秩表示模型可以有效进行优化．3 实验结果与分析将所提算法和LSA[8]，SLBF[10]，LLMC[9]，SCC[11]，MSL[4]，LRR[12]和SSC[5]算法进行比较．在WNNM-LRR和WNNM-LRR(L)中，权重参数λ取为1/logn[15]．参数γ取为1/3．对于LSA，SLBF，LLMC，SCC和MSL，须事先知道子空间的数量，而且通过图拉普拉斯的特征谱估计子空间的数量通常并不可靠．对于所有算法，子空间的数量都作为已知值．对于人脸聚类实验，使用Extended Yale B数据集[16]，包含38张人脸，每个人脸存在低维子空间中．数据集测试样本数为2 432，维度为192×168，类别数为38．3.1　聚类实验比较WNNM-LRR，WNNM-LRR(L)和其他算法在Extended Yale B数据集[16]上的性能．数据集包含38人(n=38)的192×168像素的图像，每个人都有Ni=64幅在不同光照条件下拍摄的正面人脸图像．和文献[5]中的做法类似，将图像降采样至48×42像素以减少所有算法的计算复杂度和存储要求．每个数据点是2 016维的向量．与文献[5]类似，将38人分为4组，其中1~3组人员对应于编号1~10，11~20和21~30，4组人员对应于编号31~38．所有算法都直接对原数据进行聚类．不同算法的聚类精度示于表1中．10.13245/j.hust..T001表1Extended Yale B数据集上不同算法的聚类精度算法组别1234LSA25.0025.6825.3130.44SLBF37.2638.3634.0640.12LLMC40.0943.1832.3438.54SCC29.2520.3930.3130.63MSL64.9459.5571.4142.69SSC83.1889.8994.2278.06LRR95.4485.2396.7290.32WNNM-LRR96.5484.9197.3493.08WNNM-LRR(L)94.5097.5998.1395.46%局部谱聚类算法使用数据点的K近邻估计该点的局部线性子空间，得到了较低的聚类精度．对于Extended Yale B数据集，有很多点处于不同子空间的交集中，数据点的局部邻域包含了很多其他子空间中的数据点，从而估计的局部线性子空间不准确，导致了基于该子空间建立的相似性不准确．SCC得到了较低的聚类精度，尽管全局谱聚类算法不依赖近邻点，数据中存在的噪声和异常值将会影响极曲率的估计，从而基于极曲率建立的相似性不能反映数据点间的相似性．LRR和SSC得到了比较高的聚类精度，其中LRR的结果比SSC好．LRR和SSC都估计了数据集中的全局结构，并且显式地处理噪声和异常值．得到的表示矩阵可以反映数据点间的关系，相似矩阵也是准确的．MSL得到了较低的聚类精度，该方法须对数据集进行初始聚类，并使用一个迭代过程更新初始聚类结果，很容易在初始聚类的附近陷入局部最小，因此聚类精度高度依赖于初始聚类．这里使用形状作用矩阵作为初始聚类方法．使用组1数据集，LRR和WNNM表示矩阵的奇异值(σ(z))分布情况如图1所示，i为奇异值按照从大到小的顺序排列的编号．从图1中可以看出：WNNM-LRR有更多小的奇异值，从大奇异值向小奇异值的过渡更加平缓．WNNM-LRR使用权重按照降序排列，大奇异值被更多地抑制，小奇异值被更少地抑制，从而在奇异值阈值处理过程中，保留了更多的小奇异值．而LRR在奇异值阈值处理过程中，所有的奇异值都按照相同数值抑制，造成了小奇异值的丢失．10.13245/j.hust..F001图1LRR和WNNM-LRR表示矩阵的奇异值分布对于WNNM-LRR(L)，在2~4上聚类精度进一步得到提高．当使用线性化技术时，对应于Z的表示矩阵被直接闭式求解．在ADMM中必须采用另一个变量替换Z，并引入额外的等式约束．直接求解Z可以得到表示矩阵更加准确的信息，并且使用加权核范数可以保留较小的奇异值，从而得到更加高的聚类精度．3.2　WNNM-LRR和WNNM-LRR(L)的参数分析参数γ的取值范围选择为0.1~1.0，步长是0.1．此时λ=0.18．使用WNNM-LRR和WNNM-LRR(L)对4组数据进行聚类得到的聚类精度如表2和表3所示．10.13245/j.hust..T002表2当γ=0.1~1.0时，WNNM-LRR在Extended Yale B数据集上的聚类精度组别γ0.10.20.30.40.50.60.70.80.91.0196.2395.9196.3896.7096.3897.0197.0196.8696.8695.76285.3985.5585.3985.0785.0785.3987.8087.9687.9687.64397.0396.2597.3496.5697.5096.7298.2898.5978.7579.22490.9191.7091.1193.8794.8695.8582.6183.4072.3371.74%10.13245/j.hust..T003表3当γ=0.1~1.0时，WNNM-LRR(L)在Extended Yale B数据集上的聚类精度组别γ0.10.20.30.40.50.60.70.80.91.0191.9892.7792.4594.1894.5094.3494.8195.1395.2895.28297.5997.5998.0797.7597.7597.5988.1288.1288.1288.28397.3497.9798.1398.2898.1398.1398.5998.5998.5998.44494.4794.8695.4695.2695.6595.8595.4696.0595.4683.60%从表2和3可以看出：所提算法在不同组中的表现不同；1组中WNNM-LRR的聚类精度高于LRR，而WNNM-LRR(L)的聚类精度低于LRR；WNNM-LRR的聚类精度随着参数γ的变化，幅度较小，而WNNM-LRR(L)的聚类精度随着参数γ的增加而增加；2组中WNNM-LRR和WNNM-LRR(L)的聚类精度高于LRR，WNNM-LRR(L)的精度高于WNNM-LRR；WNNM-LRR(L)的精度在γ从0.6到0.7时，下降较多；3组和4组的结果比较类似．WNNM-LRR(L)的聚类精度随着参数γ的变化，幅度较小．当参数γ靠近1.0时，WNNM-LRR的精度下降较多．在大多数情况下，WNNM-LRR和WNNM-LRR(L)的精度都高于LRR，说明所提的加权核范数是有效的．比较WNNM-LRR和WNNM-LRR(L)可以发现：在大多数情况下，WNNM-LRR(L)的聚类精度较高，并且对于不同的γ值更加鲁棒．4 结语提出基于低秩表示和加权核范数最小化的子空间聚类方法．通过ADMM和LADMM最优化算法求解，在计算矩阵范数的过程中通过加入一个参数γ，扩展了现有的矩阵范数．当0γ1时，可以得到一个介于矩阵核范数和F范数之间的范数．与基于核范数最小化的子空间聚类算法相比，所提算法在优化过程中可以保留更多小的奇异值，性能更好．