﻿ 基于正则化秩<i>k</i>矩阵逼近的稀疏主成分分析<sup>*</sup>
 文章快速检索 高级检索

Sparse principal component analysis via regularized rank-k matrix approximation
YANG Qian, LIU Hongying
School of Mathematics and Systems Science, Beijing University of Aeronautics and Astronautics, Beijing 100083, China
Received: 2016-05-31; Accepted: 2016-09-02; Published online: 2016-10-10 09:03
Foundation item: National Natural Science Foundation of China (61172060, 61403011)
Corresponding author. LIU Hongying, E-mail:liuhongying@buaa.edu.cn
Abstract: In calculating the sparse principal components (PCs), attaining k PCs simultaneously can reduce the accumulated error arising from the calculation process. We proposed the sparse principal component model via regularized rank-k matrix approximation and designed a block coordinate descent method (BCD-sPCA-rSVD) to solve this problem. Its main idea is to first divide variables into 2k blocks by coordinates, and then solve sub-problem with respect to each single coordinate block when keeping other 2k-1 variables fixed. By solving these sub-problems with explicit solutions recursively until the stopping criterion is satisfied, the BCD-sPCA-rSVD algorithm can be easily constructed. Its per-iteration complexity is linear in both sample size and variable dimensionality. The algorithm is convergent and easy to implement. Numerical simulation results show that the algorithm is feasible and effective when applied to real and synthetic data sets. The proposed method reduces the accumulated error and has lower computational complexity, which makes it well suited to handling large-scale problems of sparse principal component analysis.
Key words: dimension reduction     sparse principal component     regularization     block coordinate descent method     singular value decomposition     threshold

1 相关的稀疏PCA模型与算法

1.1 SPCA与sPCA-rSVD的模型与算法

Zou和Hastie[1]在2006年提出的SPCA，它的模型是先将PCA表述成回归型的优化问题，并添加关于回归系数的LASSO[13](l1范数)惩罚项，即

 (1)

Shen和Huang[4]在2008年提出了sPCA-rSVD，此模型首先利用了数据矩阵奇异值分解与PCA的联系，由秩1矩阵逼近得到第一主成分，然后在对应的目标函数中添加正则化项(l1l0范数等)来得到稀疏的载荷，即

 (2)

sPCA-rSVD也是利用块坐标下降法将模型(2) 的变量分成uv这2个坐标块，当固定其中一个坐标块的变量时，求解关于另外一个坐标块的子问题，交替地求解关于这2个坐标块的子问题直至终止。当正则化项为l1范数时(sPCA-rSVD-soft)，sPCA-rSVD中子问题的求解较为简单(主要的计算任务包括简单的线性回归和逐分量阈值)，但是sPCA-rSVD依次求解主成分的计算方式会产生累积误差。

1.2 BCD-SPCA的模型与算法

Zhao等[12]在2015年提出了BCD-SPCA，它考虑了含有PC载荷的l1l0范数约束的数据矩阵的奇异值分解模型，即

 (3)

2 秩k矩阵逼近的块坐标下降法

2.1 模型与算法

 (4)

 (5)

 (6)

 (7)

 (8)

 (9)

BCD-sPCA-rSVD就是利用上面所介绍的分块模式，并按照v1, u1, v2, u2, …, vk, uk的次序循环地求解关于vjuj的子问题直至满足终止条件。BCD-sPCA-rSVD的伪码描述见算法1。

1.  给定初始点V0, U0

2.  repeat

3.  for j=1, 2, …, k do

4.   计算

5.  求解关于vj的子问题(5) 来更新vj，即令

6.  计算

7.  求解关于uj的子问题(6) 来更新uj，即令uj=Ejvj

8. end for

9.  until stopping criterion satisfied

1) 与SPCA相比，BCD-sPCA-rSVD的分块模式与模型均与SPCA不同。

2) 与sPCA-rSVD相比，BCD-sPCA-rSVD能同时求多个主成分，并且模型的目标函数是含有l1范数惩罚项的数据矩阵的秩k逼近，约束是要求vj为单位向量而不是uj

3) 与BCD-SPCA相比，BCD-sPCA-rSVD的模型是把控制PC载荷的稀疏性条件作为目标函数的惩罚项，而不是作为约束条件。

2.2 计算复杂度分析

 算法 计算复杂度 n>p n＜p BCD-sPCA-rSVD O(Jk+nk2) O(Jk+ pk2) SPCA O(p3) O(Jk+ pk2) BCD-SPCA O(Jk+nk2) O(Jk+ pk2)

1) 不论对于n>p还是np的数据集，BCD-sPCA-rSVD的计算复杂度关于样本个数n与变量维数p都是线性的。

2) 在高维低样本数据集(n＜p)下，虽然3个算法每次迭代的计算复杂度完全相同，但是BCD-sPCA-rSVD和BCD-SPCA是将变量分成2k块逐次更新，而SPCA是分成两块逐次更新，理论上前面2个算法要比SPCA收敛得更快。后文3.3节和3.4节的数值结果也很好地支持了该事实。进一步，BCD-sPCA-rSVD在每次迭代又比BCD-SPCA少求解一个方程，所以BCD-sPCA-rSVD的计算效率更高。

2.3 收敛性分析

 (10)

3 数值结果与分析

BCD-sPCA-rSVD中，初始点V0的每一列为主成分载荷，U0=XV0，控制最大迭代次数为1 000，UV的更新率均取10-6。如下各实验中，这些参数的取值均保持不变。

3.1 Pitprop数据

 算法 稀疏度 非正交性 相关性 PEV/% RRE/% PCA 0 0 0 86.94 36.06 DSPCA 63 13.63 0.57 72.46 47.71 Gpower 63 17.88 0.51 75.04 45.89 SCoTLASS 27 0.32 0.44 78.24 49.24 ALSPCA 63 0 0.30 73.32 45.37 sPCA-rSVD-soft 53 14.76 0.46 76.59 46.76 SPCA 60 0.86 0.40 75.82 44.48 BCD-SPCA 63 20.05 0.40 75.86 44.19 BCD-sPCA-rSVD 63 1.51 0.28 75.13 44.18

3.2 合成数据

1) 通过观测变量X5X6X7X8来计算对应于V2的第一稀疏主成分PC1。

2) 通过观测变量X1X2X3X4来计算对应于V1的第二稀疏主成分PC2。

3) 由于(X1, X2, X3, X4)与(X5, X6, X7, X8)是相互独立的，则它们所对应的稀疏主成分不相关，并且其载荷是相互正交的。

 变量 PCA BCD-sPCA-rSVD SPCA(BCD-SPCA) PC1 PC2 PC1 PC2 PC1 PC2 X1 -0.116 3 -0.478 4 0 -0.500 0 0 -0.500 0 X2 -0.116 2 -0.478 4 0 -0.500 0 0 -0.500 0 X3 -0.116 2 -0.478 4 0 -0.500 0 0 -0.500 0 X4 -0.116 2 -0.478 3 0 -0.500 0 0 -0.500 0 X5 0.395 1 -0.145 3 0.500 0 0 0.500 0 0 X6 0.395 1 -0.145 3 0.500 0 0 0.500 0 0 X7 0.395 1 -0.145 4 0.500 0 0 0.500 0 0 X8 0.395 1 -0.145 3 0.500 0 0 0.500 0 0 X9 0.400 9 0.009 1 0 0 0 0 X10 0.400 9 0.0091 0 0 0 0 PEV/% 99.72 80.46 80.46

3.3 Colon Cancer数据

 算法 稀疏度 非正交性 相关性 PEV/% RRE/% 计算时间/s 迭代次数 PCA 0 0 0 58.35 64.54 1.98 SPCA 5 370 22.88 0.49 46.12 65.48 4.47 165 BCD-SPCA 5 373 29.15 0.54 48.88 69.73 3.65 97 BCD-sPCA-rSVD 5 376 24.02 0.48 47.24 65.34 2.76 91

3.4 20Newsgroups数据

20Newsgroups是一组低维高样本的数据集，它记录了100个单词在16 242篇报导中出现的频率，即n=16 242，p=100。所有的新闻报导均从全球最大的电子布告栏系统Usenet上取得。该数据集下载自http://cs.nyu.edu/roweis/data.html

 算法 稀疏度 非正交性 相关性 PEV/% RRE/% 时间/s 迭代次数 PCA 0 0 0 10.69 94.50 1.08 SPCA 161 0.07 0.15 8.42 95.76 2.62 192 BCD-SPCA 161 17.64 0.30 8.71 95.34 3.11 58 BCD-sPCA-rSVD 166 0.02 0.09 8.58 95.60 2.73 52

4 结论

1) 该算法减少了计算的累积误差。

2) 该算法是收敛的。

3) 该算法每次迭代的计算复杂度关于样本个数n与变量维数p都是线性的，因此它可以有效地求解大规模稀疏PCA问题。

 [1] ZOU H, HASTIE T. Sparse principal component analysis[J]. Journal of Computational and Graphical Statistics, 2006, 15(2): 265–286. DOI:10.1198/106186006X113430 [2] TRENDAFILOV N T, JOLLIFFE I T. Projected gradient approach to the numerical solution of the SCoTLASS[J]. Computational Statistics and Data Analysis, 2006, 50(1): 242–253. DOI:10.1016/j.csda.2004.07.017 [3] D'ASPREMONT A, GHAOUI L E, JORDAN M I, et al. A direct formulation for sparse PCA using semidefinite programming[J]. SIAM Review, 2007, 48(3): 434–448. [4] SHEN H, HUANG J Z. Sparse principal component analysis via regularized low rank matrix approximation[J]. Journal of Multivariate Analysis, 2008, 99(6): 1015–1034. DOI:10.1016/j.jmva.2007.06.007 [5] MOGHADDAM B, WEISS Y, AVIDAN S.Spectral bounds for sparse PCA:Exact and greedy algorithms[C]//Advances in Neural Information Processing Systems.Montreal:Neural Information Processing System Foundation, 2006:915-922. [6] D'ASPREMONT A, BACH F R, GHAOUI L E. Optimal solutions for sparse principal component analysis[J]. Machine Learning, 2008, 9(7): 1269–1294. [7] LUSS R, TEBOULLE M. Conditional gradient algorithms for rank-one matrix approximations with a sparsity constraint[J]. SIAM Review, 2013, 55(1): 65–98. DOI:10.1137/110839072 [8] LUSS R, TEBOULLE M. Convex approximations to sparse PCA via Lagrangian duality[J]. Operations Research Letters, 2011, 39(1): 57–61. DOI:10.1016/j.orl.2010.11.005 [9] SIGG C, BUHMANN J.Expectation-maximization for sparse and nonnegative PCA[C]//Proceedings of the 25th International Conference on Machine Learning.NewYork:ACM, 2008:960-967. [10] JOURNEE M, NESTEROV Y, RICHTARIK P, et al. Generalized power method for sparse principal component analysis[J]. Journal of Machine Learning Research, 2010, 11(2): 517–553. [11] LU Z S, ZHANG Y. An augmented Lagrangian approach for sparse principal component analysis[J]. Math Program Series A, 2012, 135(1-2): 149–193. DOI:10.1007/s10107-011-0452-4 [12] ZHAO Q, MENG D Y, XU Z B, et al. A block coordinates descent approach for sparse principal component analysis[J]. Neurocomputing, 2015, 153(4): 180–190. [13] TIBSHIRANI R. Regression shrinkage and selection via the LASSO[J]. Journal of the Royal Statistical Society Series B, 1996, 58(3): 267–268. [14] TSENG P. Convergence of a block coordinates descent method for nondifferentiable minimization[J]. Journal of Optimization Theory Apply, 2001, 109(3): 475–494. DOI:10.1023/A:1017501703105 [15] JEFFERS J N R. Two case studies in the application of principal component analysis[J]. Applied Statistics, 1967, 16(3): 225–236. DOI:10.2307/2985919 [16] ALON U, BARKAI N, NOTTERMAN D A, et al. Broad patterns of gene expression revealed by clustering of tumor and normal colon tissues probed by oligonucleotide arrays[J]. Proceedings of the National Academy of Sciences of the United States of America, 1999, 96(12): 6745–6750. DOI:10.1073/pnas.96.12.6745

文章信息

YANG Qian, LIU Hongying

Sparse principal component analysis via regularized rank-k matrix approximation

Journal of Beijing University of Aeronautics and Astronsutics, 2017, 43(6): 1239-1246
http://dx.doi.org/10.13700/j.bh.1001-5965.2016.0462