CN102184349A - 基于流形学习的基因表达数据的聚类系统及方法 - Google Patents

Abstract

本发明公开一种基于流形学习的基因表达数据的聚类方法，步骤为：通过获取系统获取基因表达数据矩阵A，采用局部线性平滑算法对其进行预处理；引入处理过的数据矩阵A，在三维空间中进行构建加权邻域图G；用两点间最短路径近似两点间测地线距离；用MDS计算二维嵌入坐标，将三维空间数据矩阵A映射到二维可视化空间；采用K均值聚类算法在映射后的二维可视化空间上进行聚类，得出聚类结果。此种聚类方法具有计算代价小、可以消除高阶冗余、适合于模式分类任务等特点，可以根据聚类结果有效区分细胞的当前状态、药物对恶性细胞的有效性等，获取生物体中基因的具体功能和表达调控信息。本发明还提供一种基于流形学习的基因表达数据的聚类系统。

Description

基于流形学习的基因表达数据的聚类系统及方法

技术领域

本发明属于数据挖掘处理领域，特别涉及一种基于流形学习的基因表达数据的聚类系统及方法。

背景技术

科学的进步，特别是生物信息学的迅速发展，把我们带入了一个崭新的时代，作为生命科学以及自然科学的核心、前沿领域之一，人类对于基因序列的奥秘还知之甚少，破解其生物信息以及相关的内容还有很长的路要走。高通量的微阵列检测技术是近几年生物学实验技术的一个显著突破，利用该技术可以同时对成千上万个基因的转录水平进行平行分析，根据转录产物的含量信息，从而产生了大规模的基因表达实验数据。由于细胞种类繁多，同时基因表达具有时空特异性，基因表达数据更为复杂，数据量更大，增长速度更快，如何将海量的数据中包含的丰富生物信息和知识进行迅速有效的组织、浓缩、提取、分类处理和存储就成了一个瓶颈。

面对海量的数据库以及大量繁杂的生物信息，机器学习、统计数据分析等方法在这个背景下迅速发展，巨大的计算量、复杂的噪声模式、海量的时变数据给传统的分析带来了巨大的困难。数据挖掘中的聚类分析是一种有效的数据分析工具，已广泛应用于数据处理、信号处理、图像分析等领域，而且可以用来提取、分类基因表达数据中有价值的生物信息，从而进一步提高生物信息的利用率。为了对基因表达数据进行有效的聚类，近些年来，大量的聚类算法(包括传统的聚类算法、双聚类算法等)相继被提出并应用到生物信息学的研究当中，但是目前的聚类方法主要是在二维空间上从相似性度量准则等方面进行聚类，而且假设基因表达数据存在的结构为全局线性结构，即各变量之间是独立无关的，虽然这些聚类方法的使用有助于对基因功能、细胞过程以及基因间相互调控关系的研究，但是，面对现有基因表达数据的高数据量、高维数、高增长率以及非结构化等特点，目前的这些聚类方法很难对高维的基因表达数据进行有效的聚类，发现其真实内在结构。

经对现有文献技术的检索发现，L.Zhao等人在Proceedings of the 2005 ACM SIGMOD International Conference on Management of Data，ACM，New York，NY，USA(ISBN：1-59593-060-4)中首次阐述了三维聚类方法Tricluster，并将其应用在基因表达数据上。Tricluster算法通过构建基因与样本的多重图得到双聚类结果从而再次构建双聚类与时间的多重图，进而得到三聚类结果。由此可以看出，Tricluster并不是对基因表达数据从三维上同时进行聚类操作。近年来的研究发现，真实世界中的数据往往是高维的，难以表示、理解和处理，需要采用数据降维以获得较容易分析理解的低维数据，尤其是非线性降维成为数据挖掘的一个重要手段，而且高维的基因表达数据属于本质上的低维流形，为此，研究人员将目光投向了非线性流形学习。

流形学习算法的研究是涉及拓扑学、机器学习、图论、信号处理等多个领域交叉的结果，作为一种新的非监督学习方法，它的目的是寻找嵌入在高维数据中的低维流形结构，发现高维数据分布的内在几何结构，挖掘隐藏在高维数据中的本征信息和内在规律，甚至找出产生数据集的内在变量，从而减少数据空间变化的自由度，达到维数约简或数据可视化的目的，然后通过传统的聚类算法聚类出结果，自2000年以来，学者们提出了很多流形学习算法，比如局部线性嵌入(Locally Linear Embedding，LLE)，拉普拉斯映射(Laplacian Eigenmap，LE)，局部保存投影(Locality Preserving Projection，LPP)等，可以直接获得数据点的低维映射结构，可视化效果非常出色。

上述的基因表达数据聚类算法中，传统聚类算法从单个方向进行聚类，无法发现数据的局部信息；双聚类算法可以有效地挖掘出局部信息，但无法适用于更高维的基因表达数据；而流形学习算法是一个比较新且应用前景非常广泛的技术，相比较而言，其能够克服基因表达数据的高维特征造成的影响。随着研究的不断深入及相关应用的推广，非线性流形学习也已经从最初的非监督学习发展到了监督和半监督学习，已成为机器学习相关领域的一个研究热点，针对流形学习的研究和应用也会越来越多，对生物信息学领域的研究也将更为有效及深入。

基于以上分析，本发明人将现有的聚类算法与流形算法进行有效结合，以期提高对数据的聚类效果，本案由此产生。

发明内容

本发明所要解决的技术问题，是针对前述背景技术中的缺陷和不足，提供一种基于流形学习的基因表达数据的聚类系统及方法，具有计算代价小、可以消除高阶冗余、适合于模式分类任务等特点，可以根据聚类结果有效区分细胞的当前状态(正常还是恶化)、药物对恶性细胞的有效性等，获取生物体中基因的具体功能和表达调控信息。

本发明为解决以上技术问题，所采用的技术方案是：

一种基于流形学习的基因表达数据的聚类系统，包括通过数据线相互连接的获取系统和计算机，其中，获取系统包括外壳、支架、透光载玻片、微阵列芯片、激光扫描仪、滴管和清洗装置，支架设于外壳内的底部，透光载玻片设于支架上，而微阵列芯片置于透光载玻片上；外壳内的顶部固定有激光扫描仪和滴管，而清洗装置嵌设于外壳的侧面。

上述激光扫描仪的头部嵌入外壳内。

一种基于前述基于流形学习的基因表达数据的聚类系统的方法，包括如下步骤：

(1)通过基因表达数据的获取系统获取基因表达数据矩阵A；

(2)采用局部线性平滑算法对数据矩阵A进行预处理；

(3)引入步骤(2)中处理过的数据矩阵A，在三维空间中进行构建加权邻域图G；

(4)计算最短路径，用两点间最短路径近似两点间测地线距离，距离较近的点直接用邻接矩阵D_G中的权值距离代替，而距离比较远的点则通过迭代的Dijkstra算法方法计算得到；

(5)用多维尺度变换方法计算二维嵌入坐标，将三维空间数据矩阵A映射到二维可视化空间；

(6)采用k均值聚类算法在映射后的二维可视化空间上进行聚类，得出聚类结果。

上述步骤(2)的过程为：

21)确定样本点x_i的k个(取值范围为20-30)近邻点x_ij，j＝1，2......k，记作

22)采用迭代权值选取方法计算权值向量w_i，并将其正交化∑_jw_i，j＝1；

23)采用加权PCA方法计算样本点x_i处切线空间的估计值：

${\mathrm{\Σ}}_j{w}_{i,j}{\left|\right|x_{i_j}-({\stackrel{\‾}{x_i}}^w+U_i{\mathrm{\θ}}_j^{\left(i\right)})\left|\right|}_2^2=\underset{c,U,{\mathrm{\θ}}_j}{\min }{\mathrm{\Σ}}_j{w}_{i,j}{\left|\right|x_{i_j}-(c+U{\mathrm{\θ}}_j)\left|\right|}_2^2$

其中R^m×d是嵌入在m维欧几里得空间中的d维流形，d＜m，U∈R^m×d是构成映射子空间的标准正交基，c∈R^m为映射空间的置换，||·||₂为向量的欧几里得空间标准化，θ_j代表映射子空间中点的局部坐标；

24)将样本点x_i投影到估计的切线空间内，得到样本点x_i的估计值

其中

噪音数据所对应的权值较小。

上述步骤(3)中确定邻域图G的边的方法为：计算矩阵A中每个样本点x_i同其余样本点之间的欧氏距离，当x_j是x_i的最近的k个邻近点中一个时，即认为它们彼此是相邻的，图G中存在边x_ix_j。

上述步骤(3)中确定邻域图G的边的方法为：当x_i和x_j的欧式距离小于固定值ε时，认为图G中有边x_ix_j，并以样本点间的欧式距离作为权重赋予边x_ix_j，由此可以得到邻接矩阵D_G。

上述步骤(4)的具体步骤为：

41)当邻域图G中存在边x_ix_j时，设最短路径d_G(x_i，x_j)＝d(x_i，x_j)，否则设d_G(x_i，x_j)＝∞；

42)对于空间中的所有样本点l＝1，2......N，N为样本的个数，最短路径通过下式计算：

d_G(x_i，x_j)＝min{d_G(x_i，x_j)，d_G(x_i，x_l)+d_G(x_l，x_j)}

43)通过上式可得到描述基因表达数据全局流形结构的最短路径的距离矩阵

它由图G中所有样本点之间的最短路径的平方组成。

上述步骤(5)的具体步骤为：

51)假设f为三维空间到二维空间的映射函数，y_i为样本点x_i的二维投影：y_i＝f(x_i)，则该映射的目标函数即为保持映射前后样本点之间的最短路径，即：

$f=\arg \min \underset{i,j}{\mathrm{\Σ}}{(d_G(x_i,x_j)-d(f\left(x_i\right),f\left(x_{\text{j}}\right)))}^2$

其中d(f(x_i)，f(x_j)为基因表达数据映射后二维空间中对应的欧式距离；

52)多维尺度变换方法把距离矩阵D_G作为差异性矩阵，不失一般性，假设样本点被中心化，即

则有

其中e为单位列向量，使用D_Y表示降维后二维空间中的欧式距离矩阵，τ(D_Y)为对应的内积矩阵，则上一步中目标函数等价于：

其中

53)假设映射函数f(x)＝u^Tx，则τ(D_Y)＝Y^TY＝X^Tuu^TX，前述目标函数可转化为求解特征值与特征向量：X||τ(D_G)||X^Tu＝λXX^Tu，内积矩阵τ(D_G)的两个最大特征值λ₁和λ₂以及对应的特征向量u₁和u₂所构成的映射转换矩阵为U＝[u₁，u₂]，那么T＝diag(λ₁ ^1/2，λ₂ ^1/2)U^T是二维嵌入的映射结果。

采用上述方案后，本发明采用非线性降维方法等距映射(ISOMAP)将通过微阵列技术获取的三维基因表达数据进行降维，将流形上邻近的点映射到二维可视化空间中的邻近点，不但可以保持样本的内在结构，直接得到三维基因表达数据在二维可视化空间的投影矩阵，大大降低计算成本，还可以更加准确的发现基因之间的共调控关系，提高聚类精确度。其优点主要体现在：

(1)去除了基因表达数据的奇异值以及噪声数据，从而使等距映射方法能够更有效的提取数据中的鉴别特征；

(2)保持了样本数据的流形结构信息，相比于传统的线性算法如PCA等基于全局欧式距离的假设，本发明克服了传统线性算法在处理非线性时问题时的缺点，能更好的保持样本数据的内在几何结构；

(3)可以有效处理三维甚至更高维的基因表达数据，能较好的获得二维空间的映射，计算代价小。

附图说明

图1是本发明聚类系统的结构示意图；

图2是本发明聚类方法的流程图；

图3是本发明聚类方法中基因表达数据预处理的流程图。

具体实施方式

以下将结合附图，对本发明的结构及工作流程进行详细说明。

如图1所示，本发明提供一种基于流形学习的基因表达数据的聚类系统，包括基因表达数据的获取系统S和计算机C，基因表达数据的获取系统S通过USB数据线8与计算机C连接，系统S包括一个中空的外壳7，所述外壳7的底部设置有支架6，支架6上放置有透光载玻片5，该透光载玻片5用于放置微阵列芯片2；而外壳7的顶部设置有激光扫描仪4，且激光扫描仪4的头部嵌入外壳7内；外壳7的顶部还设有滴管1，而外壳7的侧面则嵌入一个清洗装置3。

前述聚类系统的工作原理为：基因表达数据获取时，通过将滴管1中带有荧光标记的待测样品加入到微阵列芯片2上，与探针在65℃下杂交16小时，杂交反应完成后，使用装有清洗液的清洗装置3对杂交反应后的微阵列芯片2在室温下清洗5分钟，然后用激光扫描仪4扫描微阵列芯片2并获取荧光图像，将获得的图像通过USB数据线8输入到计算机C中，由计算机C对获得的图像进行分析，得到微阵列芯片2上每一点的荧光强度值，即基因表达数据，最后使用计算机C对基因表达数据进行预处理、聚类，得到聚类结果。

再请参考图2所示，基于以上聚类系统的一种基于流形学习的基因表达数据的聚类方法，其具体步骤如下：

(1)通过系统S获取基因表达数据矩阵A；

(2)采用局部线性平滑算法对基因表达数据矩阵A进行预处理，主要是采取局部加权主成分分析(Principal Component Analysis，PCA)方法以及迭代的权值选取规则达到实现鲁棒的局部线性平滑处理，去除奇异值与噪声数据，详细过程按下述步骤进行：

21)确定样本点x_i的k(取值范围为20-30)个近邻点x_ij，其中j＝1，2......k，则k个近邻点分别记作

22)采用迭代权值选取方法计算权值向量w_i，并将其正交化∑_jw_i，j＝1；

23)采用加权PCA方法计算样本点x_i处切线空间的估计值：

${\mathrm{\Σ}}_j{w}_{i,j}{\left|\right|x_{i_j}-({\stackrel{\‾}{x_i}}^w+U_i{\mathrm{\θ}}_j^{\left(i\right)})\left|\right|}_2^2=\underset{c,U,{\mathrm{\θ}}_j}{\min }{\mathrm{\Σ}}_j{w}_{i,j}{\left|\right|x_{i_j}-(c+U{\mathrm{\θ}}_j)\left|\right|}_2^2$

其中R^m×d是嵌入在m维欧几里得空间中的d维流形，d＜m，U∈R^m×d是构成映射子空间的标准正交基，c∈R^m为映射空间的置换，||·||₂为向量的欧几里得空间标准化，θ_j代表映射子空间中点的局部坐标；

24)将样本点x_i投影到估计的切线空间内，得到样本点x_i的估计值

其中

噪音数据所对应的权值较小；

(3)引入步骤(2)中处理过的数据矩阵A，在三维空间中进行构建加权邻域图G；

主要有两种方法确定邻域图G的边：一是采用k近邻方法，计算矩阵A中每个样本点x_i同其余样本点之间的欧氏距离，当x_j是x_i的最近的k个邻近点中一个时，即认为它们彼此是相邻的，图G中存在边x_ix_j；二是采用ε阈值法，当x_i和x_j的欧式距离小于固定值ε时，认为图G中有边x_ix_j，并以样本点间的欧式距离作为权重赋予边x_ix_j，由此可以得到邻接矩阵D_G，本方案采用了k近邻方法；

(4)计算最短路径，用两点间最短路径近似两点间测地线距离，具体可以分为两种情况：距离较近的点直接用邻接矩阵D_G中的权值距离代替，而距离比较远的点则通过迭代的Dijkstra算法方法计算得到。具体步骤如下：

41)当图G中存在边x_ix_j时，设最短路径d_G(x_i，x_j)＝d(x_i，x_j)，否则设d_G(x_i，x_j)＝∞；

42)对于空间中的所有样本点l＝1，2，…，N，N为样本的个数，采用迭代的Dijkstra算法计算最短路径，最短路径通过下式计算：

d_G(x_i，x_j)＝min{d_G(x_i，x_j)，d_G(x_i，x_l)+d_G(x_l，x_j)}

43)通过上式可得到描述基因表达数据全局流形结构的最短路径的距离矩阵它由图G中所有样本点之间的最短路径的平方组成；

(5)用多维尺度变换方法计算二维嵌入坐标，将三维空间数据矩阵A映射到二维可视化空间。通过保持映射前后基因表达数据样本点之间的最短路径而给出目标函数，通过求解此目标函数求取基因表达数据矩阵二维映射矩阵，使目标函数的求解转化为特征值及其特征向量的求解问题；

51)假设f为三维空间到二维空间的映射函数，y_i为样本点x_i的二维投影：y_i＝f(x_i)，则该映射的目标函数即为保持映射前后样本点之间的最短路径，即：

$f=\arg \min \underset{i,j}{\mathrm{\Σ}}{(d_G(x_i,x_j)-d(f\left(x_i\right),f\left(x_{\text{j}}\right)))}^2$

其中d(f(x_i)，f(x_j)为基因表达数据映射后二维空间中对应的欧式距离；

52)多维尺度变换方法把样本之间的距离矩阵D_G作为差异性矩阵，不失一般性，假设样本点被中心化，即

则有

其中e为单位列向量，使用D_Y表示降维后二维空间中的欧式距离矩阵，τ(D_Y)为对应的内积矩阵，则上一步目标函数等价于：

其中

53)假设映射函数f(x)＝u^Tx，则τ(D_Y)＝Y^TY＝X^Tuu^TX，上述目标函数可转化为：X||τ(D_G)||X^Tu＝λXX^Tu，τ(D_G)的两个最大特征值λ₁和λ₂以及对应的特征向量u₁和u₂所构成的映射转换矩阵为U＝[u₁，u₂]，那么T＝diag(λ₁ ^1/2，λ₂ ^1/2)U^T是二维嵌入的映射结果；

(6)采用k均值聚类算法在映射后的二维可视化空间上进行聚类，得出聚类结果。

以上实施例仅为说明本发明的技术思想，不能以此限定本发明的保护范围，凡是按照本发明提出的技术思想，在技术方案基础上所做的任何改动，均落入本发明保护范围之内。

Claims (8)

1.一种基于流形学习的基因表达数据的聚类系统，其特征在于：包括通过数据线相互连接的获取系统和计算机，其中，获取系统包括外壳、支架、透光载玻片、微阵列芯片、激光扫描仪、滴管和清洗装置，支架设于外壳内的底部，透光载玻片设于支架上，而微阵列芯片置于透光载玻片上；外壳内的顶部固定有激光扫描仪和滴管，而清洗装置嵌设于外壳的侧面。

2.如权利要求1所述的基于流形学习的基因表达数据的聚类系统，其特征在于：所述激光扫描仪的头部嵌入外壳内。

3.一种基于如权利要求1所述的基于流形学习的基因表达数据的聚类系统的方法，其特征在于包括如下步骤：

(1)通过基因表达数据的获取系统获取基因表达数据矩阵A；

(2)采用局部线性平滑算法对数据矩阵A进行预处理；

(3)引入步骤(2)中处理过的数据矩阵A，在三维空间中进行构建加权邻域图G；

(4)计算最短路径，用两点间最短路径近似两点间测地线距离，距离较近的点直接用邻接矩阵D_G中的权值距离代替，而距离比较远的点则通过迭代的Dijkstra算法方法计算得到；

(5)用多维尺度变换方法计算二维嵌入坐标，将三维空间数据矩阵A映射到二维可视化空间；

(6)采用k均值聚类算法在映射后的二维可视化空间上进行聚类，得出聚类结果。

4.如权利要求3所述的基于流形学习的基因表达数据的聚类方法，其特征在于所述步骤(2)的过程为：

21)确定样本点x_i的k个(取值范围为20-30)近邻点x_ij，j＝1，2......k，记作

22)采用迭代权值选取方法计算权值向量w_i，并将其正交化∑_jw_i，j＝1；

23)采用加权PCA方法计算样本点x_i处切线空间的估计值：

${\mathrm{\Σ}}_j{w}_{i,j}{\left|\right|x_{i_j}-({\stackrel{\‾}{x_i}}^w+U_i{\mathrm{\θ}}_j^{\left(i\right)})\left|\right|}_2^2=\underset{c,U,{\mathrm{\θ}}_j}{\min }{\mathrm{\Σ}}_j{w}_{i,j}{\left|\right|x_{i_j}-(c+U{\mathrm{\θ}}_j)\left|\right|}_2^2$

其中R^m×d是嵌入在m维欧几里得空间中的d维流形，d＜m，U∈R^m×d是构成映射子空间的标准正交基，c∈R^m为映射空间的置换，||·||₂为向量的欧几里得空间标准化，θ_j代表映射子空间中点的局部坐标；

24)将样本点x_i投影到估计的切线空间内，得到样本点x_i的估计值

其中

噪音数据所对应的权值较小。

5.如权利要求3所述的基于流形学习的基因表达数据的聚类方法，其特征在于所述步骤(3)中确定邻域图G的边的方法为：计算矩阵A中每个样本点x_i同其余样本点之间的欧氏距离，当x_j是x_i的最近的k个邻近点中一个时，即认为它们彼此是相邻的，图G中存在边x_ix_j。

6.如权利要求3所述的基于流形学习的基因表达数据的聚类方法，其特征在于所述步骤(3)中确定邻域图G的边的方法为：当x_i和x_j的欧式距离小于固定值ε时，认为图G中有边x_ix_j，并以样本点间的欧式距离作为权重赋予边x_ix_j，由此可以得到邻接矩阵D_G。

7.如权利要求3所述的基于流形学习的基因表达数据的聚类方法，其特征在于所述步骤(4)的具体步骤为：

41)当邻域图G中存在边x_ix_j时，设最短路径d_G(x_i，x_j)＝d(x_i，x_j)，否则设d_G(x_i，x_j)＝∞；

42)对于空间中的所有样本点l＝1，2......N，N为样本的个数，最短路径通过下式计算：

d_G(x_i，x_j)＝min{d_G(x_i，x_j)，d_G(x_i，x_l)+d_G(x_l，x_j)}

43)通过上式可得到描述基因表达数据全局流形结构的最短路径的距离矩阵它由图G中所有样本点之间的最短路径的平方组成。

8.如权利要求3所述的基于流形学习的基因表达数据的聚类方法，其特征在于所述步骤(5)的具体步骤为：

51)假设f为三维空间到二维空间的映射函数，y_i为样本点x_i的二维投影：y_i＝f(x_i)，则该映射的目标函数即为保持映射前后样本点之间的最短路径，即：

$f=\arg \min \underset{i,j}{\mathrm{\Σ}}{(d_G(x_i,x_j)-d(f\left(x_i\right),f\left(x_{\text{j}}\right)))}^2$

其中d(f(x_i)，f(x_j)为基因表达数据映射后二维空间中对应的欧式距离；

52)多维尺度变换方法把距离矩阵D_G作为差异性矩阵，不失一般性，假设样本点被中心化，即则有

其中e为单位列向量，使用D_Y表示降维后二维空间中的欧式距离矩阵，τ(D_Y)为对应的内积矩阵，则上一步中目标函数等价于：

其中

53)假设映射函数f(x)＝u^Tx，则τ(D_Y)＝Y^TY＝X^Tuu^TX，前述目标函数可转化为求解特征值与特征向量：X||τ(D_G)||X^Tu＝λXX^Tu，内积矩阵τ(D_G)的两个最大特征值λ₁和λ₂以及对应的特征向量u₁和u₂所构成的映射转换矩阵为U＝[u₁，u₂]，那么T＝diag(λ₁ ^1/2，λ₂ ^1/2)U^T是二维嵌入的映射结果。

Images