非线性主成分分析中数据非线性特征的检验方法*
2011-03-11山东大学公共卫生学院250012高青松薛付忠
山东大学公共卫生学院(250012) 高青松 薛付忠
自从Karl Pearson〔1〕提出主成分分析(PCA)以来,这种多元统计分析方法已经被广泛应用于各个领域,如图像处理、故障诊断、综合分析等。主成分分析的基本思想是通过降维(线性变换)的思想,在损失很少信息的前提下把原始数据中的多个变量化为少数几个综合指标(即主成分),从而达到简化系统结构,抓住问题本质的目的。
传统的线性主成分分析可以有效地处理变量间的线性关系,但是现实中有很多数据是非线性的,而线性主成分分析并不能检测出数据之间的非线性结构〔2〕。另外,线性主成分分析也不能准确地代表不满足高斯分布的数据〔3〕。因此,很多文献提出用非线性的主成分分析方法来提取数据之间的非线性特征,主要包括核主成分分析(KPCA)〔4〕和基于神经网络(neural network)〔2,5〕的主成分分析等。尤其是核主成分分析,不仅可以提取数据的非线性特征,而且其非线性主成分计算简单(求解特征值问题),意义明确(特征空间中的线性主成分),因此在很多领域得到广泛应用。
在进行核主成分分析时,如果选择线性核函数,则核主成分分析等价于线性主成分分析,即在选择正确的核函数的前提下,核主成分分析可用于分析线性数据。然而,到目前为止尚未提出一种可以有效地选择核函数及其参数的方法,如果盲目地对线性数据进行核主成分分析,不仅核函数的选择及其参数的设置费时费力,而且会大大增加计算复杂度;如果选择非线性的核函数,无法代表原始数据的非线性结构,从而降低其效能。所以,在进行主成分分析(线性或非线性)之前,必须对数据是线性结构还是非线性结构进行检验。
Uwe Kruger等在过程控制(process control)领域提出了一个非线性检验的度量〔6〕,但此方法需要对数据有一定的先验知识,或者可以根据其前几个主成分的散点图来进行分区,然而我们在进行数据处理时,绝大部分情况下对数据的了解还处于探索阶段,不大可能有明确的先验知识,而且大部分数据很难根据其主成分的散点图来进行分区,除非所选择的样本有很强的层次结构。基于此,本文在上述非线性检验度量的基础上,提出用主成分聚类分析法来进行分区,减少了主观因素对分区的影响,并可有效地处理主成分的散点图无法分区的情形,从而得到一个实用性更强的非线性检验度量。
原理与方法
1.原理 首先将原始数据中的所有样本分成若干区域,然后任意选择一个区域,计算该区域内变量之间的相关系数矩阵,对于相关系数矩阵中的每一个元素,根据变量的均值和方差计算出其置信限;并根据此置信限估计出主成分分析的残差(也就是丢弃的特征值之和)的最大值和最小值,或简称为精确边界;最后计算每一个区域按照选定区域的均值和方差进行标准化之后,进行主成分分析时的残差。如果所有的残差都落在估计出的精确边界以内,则说明该数据变量之间是线性关系;当至少有一个残差落在精确边界以外时,就判断变量之间是非线性关系。
由于最初是任意选择一个区域来计算精确边界,为了避免可能的偏倚,本文采用交叉验证的原理,对每个区域都估计其精确边界并计算残差,然后按照上述方法来判断变量之间的关系。
2.方法 具体的检验过程分为以下四步:
第一步:定义不连续区域 首先对原始数据进行主成分分析,选出前几个累积贡献率达到80%的主成分,然后对选出的主成分进行K均值聚类(K-means clustering)。
第二步:计算相关系数矩阵中每个元素的置信限给定原始数据集有N个变量,则其相关系数矩阵为
假设将原始数据集分成m个区域,则第h(h=1,…,m)个区域的相关系数矩阵同样可根据该区域中的数据得到,记为R(h)。根据每个变量的均值和方差可计算出R(h)中每个元素的置信限,用矩阵表示如下
第三步:确定精确边界 根据式(1)中相关系数矩阵的置信限可直接估计精确边界。由主成分分析中残差矩阵的弗罗贝尼乌斯范数(frobenius norm)可以证明,残差σ等价于主成分模型中丢弃的特征值之和:
其中σj代表第j个变量的残差,λk代表相关系数矩阵R(h)的第k个特征值,n(h)是该区域中的样本数,n是保留的主成分的个数。
由于特征值 λn+1,…,λN依赖于相关系数矩阵R(h)的元素,所以精确边界的估计可转化为求解下面的最优化问题:
其中,ΔRmax和ΔRmin中的元素分别是R中非对角元素的波动,用于决定最大值λkmax和最小值λkmin,然后根据公式(2)可得第h个区域的精确边界为
本文对上述最优化问题采用粒子群优化算法来求解(particle swarm optimization,PSO)〔7〕。
第四步:定义非线性度量对每一个h(h=1,2,…,m),按照前面三步得到其精确边界σhmax,σhmin和残差以后,对其进行比较,如果所有的残差都落在精确边界之内,则该数据为线性数据;当至少有一个落在精确边界之外,则判断为非线性数据。
实例分析
为验证该非线性检验度量的有效性,我们首先用该度量分析了文献〔6〕中的两份模拟数据,包括一个线性的例子和一个非线性的例子。接下来,又用该度量对北美风湿性关节炎全基因组关联分析的(SNP)数据〔8〕的基因PTPN22基因区域内9个单核苷酸多态性(SNPs)之间的线性或非线性特征进行判断,为我们进一步选择核主成分分析构建基因组区域化关联分析检验统计量提供依据。
1.线性例子 首先产生1000个服从均匀分布的样本,范围在-4到4之间,作为随机变量的值,然后产生两个相互独立且均服从均值为0,方差为0.005的正态分布的随机变量和,则线性情况下可根据如下公式产生两个变量z1和z2:
z1=e1+t,z2=t+e2
图1 线性例子的散点图及其分区结果
图1给出了线性情况下变量z1和z2的散点图,并给出了主成分聚类分析法的分区结果。表1给出了每个区域中相关系数的值及其置信限,并计算了每个区域的精确边界。表2给出了每个区域中主成分分析的残差,也就是相关系数矩阵的第二个特征值。其中对角元素精确边界是根据该区域的数据估计得到的。通过对表1和表2的结果进行比较可得,没有落在精确边界以外的数值,这跟预期结果一样,即我们假设的变量z1和z2之间的关系是线性的。
表1 线性例子中每个区域内的相关系数及其置信限和精确边界
表2 每个区域中主成分分析的残差
2.非线性例子 同样,首先产生1000个服从均匀分布的样本,范围在-4到4之间,作为随机变量t的值,然后产生两个相互独立且均服从均值为0,方差为0.005的正态分布的随机变量e1和e2,非线性情况下变量z1和z2可根据下面的公式产生:
图2 非线性例子的散点图及其分区结果
图2给出了非线性情况下变量z1和z2的散点图,同样给出了主成分聚类分析法的分区结果。表3给出了每个区域中相关系数的值及其置信限,并计算了每个区域的精确边界。表4给出了每个区域中主成分分析的残差。与线性情况不同的是,每个区域都有落在精确边界以外的数值,在表4中用粗体表示,这是因为我们假定的变量z1和z2之间的关系是非线性的。
3.实际数据 在全基因组关联分析中,为了减少多重检验造成的假阳性、数据共线性特征等所致的困扰,需要构建基因组区域化关联分析检验统计量,本课题组已经构建了一种基于线性主成分bootstrap可信区间法的检验统计量〔8〕。然而,基因组SNP之间的线性关系的假设往往不成立,需要进一步采用非线性主成分分析(核主成分分析)构建统计效能更高的基因组区域化关联分析统计量。为了验证本文所提出的非线性检验度量的有效性,我们对北美风湿性关节炎全基因组关联分析的(SNP)数据〔9〕的基因PTPN22基因区域内9个单核苷酸多态性(SNPs)之间的线性或非线性特征进行判断,从而为进一步选择核主成分分析构建基因组区域化关联分析检验统计量提供依据。该数据包括868个病例和1 194个对照。PTPN22基因区域中的9个SNPs(rs114157960-rs114215857)的连锁不平衡区域结构(LD block)见图3。
表3 非线性例子中的相关系数及其置信限和精确边界
表4 每个区域中主成分分析的残差
图3 PTPN22基因中9个SNP的LD图
对该数据进行主成分分析,其前三个主成分的累积贡献率为80.5%,故将其他6个主成分舍去。图4(A)给出了前三个主成分的散点图,从该图中,我们无法得到任何分区结构,而根据主成分聚类分析法可以很清晰地将原始数据分成四个区域(图4(B))。
表5给出了这9个SNP的精确边界及特征值。精确边界与预期的一样,这些值落在区域之内。然而,非对角元素至少一个值落在精确边界之外(用粗体表示),这就意味着这些SNPs之间是非线性关系,从而表明在构建基因组全区域化关联分析检验统计量时采用非线性主成分分析(核主成分分析)的必要性。
讨 论
图4 PTPN22基因中9个SNP前三个主成分的散点图(A)及其分区结果(B)
表5 PTPN22基因中9个SNP的非线性检验结果
非线性的核主成分分析由于其强大的非线性特征提取能力以及其原理简单而被广泛应用于图像识别〔10〕、特征提取〔11〕和综合评价〔12〕等多个领域。然而,核主成分分析的核函数选择及其参数调试问题仍未能得到很好的解决。所以,在进行核主成分分析之前,应该判断是否有必要对数据进行这种复杂的分析,也就是说,应判断数据是否为非线性数据。目前,对数据进行非线性主成分分析,绝大部分研究都只是简单地将数据的非线性作为一个“前提条件(a prior)”〔6〕。为了研究数据的结构,Uwe Kruger等人提出了一个非线性检验的度量,但该度量并不能很好地应用于其他领域,比如说全基因组关联研究(GWAS)得到的单核苷酸多态性(SNP)数据,目前对SNP的研究还处于探索阶段,所以不可能根据其先验知识来进行分区,而且大部分的SNP数据并不能根据其前几个主成分的散点图来进行分区,除非原始数据有很好的群体结构。
本文提出的基于主成分聚类分析法的非线性检验度量,有以下优点:首先,该度量避免了先验知识匮乏或前几个主成分的散点图不能给出任何分区结构的不足。主成分聚类分析法是对主成分分析与聚类分析方法的综合利用,利用主成分分析的结果作为聚类分析的样本矩阵,不仅减少了数据的冗余信息,使得聚类计算较为简单,而且原理清晰,所得结论客观实际,可靠性强。其次,在根据相关系数矩阵每一个元素的置信限估计其精确边界时,本文采用粒子群优化算法来求解最优化问题〔7〕。粒子群优化算法属于进化算法的一种,和遗传算法相似,它也是从随机解出发,通过迭代寻找最优解,并通过适应度来评价解的优良。但粒子群优化算法比遗传算法更容易实现,并且没有许多参数需要调整,目前这种算法已广泛应用于函数优化,神经网络训练,模糊系统控制以及其他遗传算法的应用领域。
两份模拟数据分别验证了线性模型和非线性模型中该方法的有效性。在线性例子中,所有的残差都落在估计出的精确边界之内,说明该模型可以很好地检测出数据之间的线性结构;而在非线性的例子中,至少有一个残差落在精确边界之外,从而说明该非线性检验度量对非线性模型的有效性。对实际基因组关联分析数据PTPN22基因中的9个SNP的分析结果表明,该数据是非线性数据,所以应该采用非线性主成分分析构建基因组区域化关联检验统计量;如果采用线性的主成分分析构建,则不能代表原来的数据结构,从而不能充分地挖掘SNP数据中所蕴含的大量信息。
通过本文提出的非线性检验度量,我们可以更多地了解数据的结构,从而选择最优的模型来挖掘数据中所蕴含的信息。
1.Pearson K.Principal components analysis.The London,Edinburgh,and Dublin Philosophical Magazine and Journal of Science,1901,6(2):559.
2.Botelho SSC,Bem RA,Almeida IL,et al.C-NLPCA:extracting non-linear principal components of image dataset.Simpoio Brasileiro de Sensoriamento Remoto,2005:3495-3502.
3.Heo G,Gader P,Frigui H.2009 Special Issue:RKF-PCA:robust kernel fuzzy PCA.Neural Networks,2009,22(5-6):642-650.
4.Schlkopf B,Smola A,Muller K.Kernel principal component analysis.Artificial Neural Networks-ICANN'97,1997:583-588.
5.Kramer MA.Nonlinear principal component analysis using autoassociative neural networks.AIChE Journal,1991,37(2):233-243.
6.Kruger U,Antory D,Hahn J,et al.Introduction of a nonlinearity measure for principal component models.Computers & Chemical Engineering,2005,29(11-12):2355-2362.
7.Kennedy J,Eberhart R.Particle swarm optimization.IEEE International Conference on Neural Networks,1995:1942-1948.
8.Peng Q,Zhao J,Xue F.PCA-based bootstrap confidence interval tests for gene-disease association involving multiple SNPs.BMC Genet,2010,11:6.
9.Plenge RM,Seielstad M,Padyukov L,et al.TRAF1-C5 as a risk locus for rheumatoid arthritis-A genomewide study.New England Journal of Medicine,2007,357(12):1199-1209.
10.李翊,张静,吴凌华,等.一种基于改进主成分分析的SAR图像识别方法研究.海军航空工程学院学报,2009,24(3):307-312.
11.李岳,温熙森,吕克洪.基于核主成分分析的铁谱磨粒特征提取方法研究.国防科技大学学报,2007,29(2):113-116.
12.祝荣欣,权龙哲,乔金友,等.核主成分分析在农业机械性能综合评价中的应用.中国造船,2009(9):187-189.