CN1710558A - 基于主集分割的基因芯片表达谱数据聚类方法 - Google Patents

Abstract

本发明属于数据挖掘技术领域，具体涉及一种基于主集分割的基因芯片表达谱数据聚类方法。其步骤是，先利用主集方法对原始数据进行重排序，使高度相似的基因排列在同一特定区域，由此定位一个类；将找到的类从当前数据集中排除，对剩下的数据重复以上处理，直至满足所提出的停止的条件为止。本方法算法简单，易于实现；数据输入顺序不敏感，且抗噪能力强。结果中的类质量好，不相似的元素不会划分在一个类中，而且根据需要具有延展性。

Description

基于主集分割的基因芯片表达谱数据聚类方法

技术领域

本发明属于计算机算法、数据挖掘技术领域，具体涉及基于主集分割的基因芯片表达谱数据聚类方法。

背景技术

基因芯片技术是90年代发展起来的高效、快速的生物实验技术，具有并行化、集约化和多参量等特点。将基因芯片技术与大规模数据统计以及数据挖掘方法相结合，是当今进行生物研究的重要方法之一。扫描芯片实验后带荧光信号的芯片图像而得到的大量的数据被称为基因芯片数据(microarray data)，也称为基因表达谱数据(gene expression profiles)或简称为表达谱数据，它表示了大量基因在不同实验条件下的表达谱，具有典型的高噪、高变异、高维、高相关、维数差异大等特点，通过对表达谱数据进行统计分析或数据挖掘分析，可以找到具有生物特征意义的信息。

聚类分析作为数据挖掘的一种技术，被广泛的应用于基因芯片数据分析。早期生物实验证明表达相似的基因具有相似功能，基于这一假设，聚类被用来寻找相似功能的基因/样本，进行基因功能预测等等。现在有很多种算法都较为广泛地运用于表达谱数据聚类分析中，并且解决了一些问题，例如Kmeans聚类、层次聚类、自组织映射神经网络(Self-Organizing Maps)，Self-Organizing Tree Algorithm，Quality-based Clustering等已被成功的应用于高维基因芯片数据的分析中。大多数聚类算法起源于非生物相关的研究领域，在实际应用中，存在一些不足之处。如Kmeans方法在处理较大数据量时的扩展性和效率都较理想，但需要预先定义类的个数，且可能陷入局部最优，结果受噪声和异常值的影响较大；自组织图(SOM)可以处理部分数据及含缺失值的数据，算法稳健，结果易于可视化，且可用于模糊聚类，但要预先确定的参数太多，这些参数具有一定经验性，改变参数数值往往会极大影响聚类结果；层级聚类可以得到基因之间层级性的关系，结果用树状图表示，易于可视化，但结果不稳定，受初始变量和所选参数的影响较大，且不能直接得到相似基因组成的类。且在绝大多数方法中，所有元素都被划分到某个类中，即便是某个元素与类中其他元素的相似性并不高，但由于与其它类相似性更低，因此只能被归并到这类中，造成了该类的″噪声″。一个类的总体表达效果将受噪声影响，这样的结果也将不利于进一步的分析。

发明内容

本发明目的在于提出一种所用参数比较少，结果受噪声影响小的基因芯片表达谱数据聚类方法。

本发明提出的基因芯片表达谱数据聚类方法，是一种迭代的聚类算法，我们称之为主集分割法(Ds_Clust)。从本质上讲，该算法是一个启发式的两步方法，依次地确定各个类(类的个数是预先未知的，不是算法的参数)。其具体步骤如下：

首先，利用主集方法对基因或样本(原始数据)进行重排序，使表达谱高度相似的基因排列在同一特定区域，由此可以定位一个类；类的分割界线通常难于确定，本发明提出了一种定位标准，根据类内元素间的平均相似性远高于类外元素间的平均相似性的性质，从排序后的基因或样本序列划分出一个具有高相似表达的基因或样本的类；将找到的类从当前数据集中排除，对剩下的数据重复以上处理，直到满足所提出的循环停止条件为止。

将基因或样本看作是高维空间中的点，或是图像中的像素的分割。将基因或样本之间的相似性关系类比于图像中像素之间的相似性关系，可以借用图像分割中的算法解决基因或样本聚类/分类问题。

本发明中提出的主集分割法(Ds_Clust)，它主要有以下几个优点：

(1)算法过程简单，易于实现。

(2)不用预先输入类的个数，减少了对先验知识的以及专业知识的依赖。

(3)结果中的类质量好。可以找到高度相似的元素，不会将不相似的元素划分在一个类中。

(4)对数据的输入顺序不敏感。

(5)结果形式简单，易于使用。

(6)结果类根据需要具有延展性

(7)抗噪能力强。

附图说明

图1为对原始序列进行重新排序图示。

图2为整个算法的实现过程图示。其中，C_i，R_i分别是第i个类的中心和半径，i＝1，2…。

具体实施方式

下面进一步描述本发明的具体实施方式：

(1)用主集方法对原基因或样本序列进行重排序。

主集的定义与类的定义具有相似性，必须满足两个条件：主集内部的元素具有相同性质，主集内部的元素与其外的元素性质不同。可以用二次型法，通过使用迭代的方法求解单纯型的极大值的过程来求解主集。

对于一个带约束条件的极值问题：

最大化：f(x)＝x′Wx    并满足：x∈Δ   (1)

这里：Δ＝{x∈R″：x_i≥0 for all i∈V and e′x＝1}，

这个约束条件称作R″上的标准单纯形，e是全1列向量。

通过求解迭代方程，可以求得二次型的极大值。使用离散动态方程：

$x_i(t+1)=x_i\left(t\right)\frac{{\left(\mathrm{WX}\left(t\right)\right)}_i}{X^T\left(t\right)\mathrm{WX}\left(t\right)}$ (2)

其中x_i(t)是向量X(t)的分量，t是迭代步数；向量x的支持集可以定义为向量x中非零元素的下标集合，也就是：

           σ(x)＝{i∈V:x_i＞0}         (3)

通过求解二次型问题(1)的局部最大值获得向量x的支持集，支持集所对应的元素所组成的类为最初意义上的主集。在应用过程中，我们注意到：二次型(1)的局部最优解x实际上给我们提供了类的隶属度信息，x的分量大小表示分量对应的点属于当前类的可能性。而使用向量x的支持集作为分割标准在高噪音、高相关的基因表达谱数据聚类中并不适用。因此，可以利用向量x中所包含的类的隶属度信息，对原基因或样本进行重排序，通过重排序使得具有高相似表达的基因/样本排列在整个序列的某个特定位置。

用G表示基因表达谱数据矩阵，G的每一行代表一种基因，每一列代表一种实验条件，G(i，k)表示基因i在条件k下的表达值，G中基因之间的相互关系可以用一个关系矩阵A＝(a_ij)表示，A是一个对称矩阵，其元素a_ij表示对应两个元素i和j的相关性，其值越大，表示两元素之间的相关性越大，相反，元素间的相关性就越小：

$a_{\mathrm{ij}}=\exp (-\frac{\left|\right|g_i-g_j\left|\right|}{\mathrm{\δ}})$ (4)

其中δ＞0是缩放因子，起到调节作用，以控制聚类敏感度。g_i和g_j分别表示第i和j个元素的数据向量，在对基因表达谱进行聚类时，g_i和g_j就是G的第i行和j行基因表达向量。还可以通过多种方式计算元素间的相似性，例如欧式距离，Pearson系数等等。

令W＝A，通过公式(2)计算向量X(t)，其每一分量x_i(t)对应原数据集G中的一条基因或样本向量，数值越大的分量对应的基因或样本间相似性越高。按照x_i(t)的值从大到小对原始基因或样本进行排序，并将排序后对应的基因表达向量依次记为g₁′g₂′，…，g_n′。此时，相似性高的元素排列在序列前段，该区域元素密度高，以此作为聚类的中心，是合理的推测。这样一来，使我们在下一步计算中有了明确的方向，极大的简化了计算。

(2)根据类内元素的平均相似性远高于类外元素的平均相似性的标准，从重新排序之后的基因/样本序列中划分出一个具有高相似性的基因/样本的类。在排序后的基因序列g₁′，g₂′，…，g_n′中，属于同一个主集的基因由于对应的x_i(t)较大会排列在序列尾部(若按从小到大排列)，从该序列中划分出一个类，具体步骤如下：

a.以欧式距离作为相似性的判断，先判断g₁′和g₂′之间的距离。假如g₁′和g₂′之间的距离大于一个初始设定的距离值(该值的确定具有经验性)，则停止分割，并记录分割失败的累计次数。若是在第(2)步中未找到合理聚类的连续次数大于某一数值，则停止整个算法的迭代过程；否则，将g₁′和g₂′从剩余数据集中去除掉，退出步骤(2)返回步骤(1)；若g₁′和g₂′之间的距离小于上述初始设定的距离值，则将g₁′和g₂′作为初始类，并用g₁′和g₂′的平均值作为初始类中心。

b.接着判断g₃′是否属于这个类，若不属于，则停止分割；若属于，则将g₃′加入类中，并用相同的方法，反复迭代步骤b依次判断g₄′，g₅′等是否属于这一类，直到找到一个聚类，并输出该类元素为止。

那么，假如{g₁′，g₂′…，g_k′}属于同一个类，如何判断g_k+1′是否属于这个类呢？令D＝{g₁′，g₂′…，g_k′}为类内元素的集合，U＝{g_k+2′，g_k+3′…，g_n′}为类外元素集合，则：

$d_D=\frac{1}{k}\underset{g_i\∈D}{\mathrm{\Σ}}{d}_{i,k+1}$ (5)

$d_U=\frac{1}{n-k-1}\underset{g_i\∈U}{\mathrm{\Σ}}{d}_{i,k+1}$ (6)

其中d_i，k+1表示g_i′和g_k+1′之间的距离，本文中采用欧氏空间距离。

根据类的性质，若

                   d_U＞d_D                              (7)

即g_k+1′到类{g₁′，g₂′…，g_k′}的距离平均值小于g_k+1′到其他非主集基因的距离平均值，那么认为g_k+1′是属于这个类。此判断对于类的性质的规定，将影响最终的聚类结果。由于在对基因表达谱进行聚类时，目的是寻找表达谱数据高度相似的基因，而条件(7)过于宽松，我们采用的判断条件是d_U≥3d_D，这样规定将有助于减少结果类中的噪音。若放松此约束条件，则结果中类的半径增大、个数减小。

(3)将第(2)步中找到的类从剩余的数据中去除掉，如果数据集中存在未被分类的元素，则返回第(1)步继续执行；若整个数据集已被完全分类，则停止迭代。

应用实例

附图1显示了通过上述算法中的第(1)步对原始序列进行重新排序，并对排序后的元素序号重新标记的结果，重排序后的序列用g₁′，g₂′，…，g_n′表示。

该算法的整个过程可以用图2来表示。每次循环过程中找到一个类，然后从原始数据集中去除当前找到的类，并继续处理剩余数据。结果形式为各个找到的类，这些类的半径、中心由原始数据集本身的性质决定，并且互不相关。

Claims (2)

1、一种基于主集分割的基因芯片表达谱数据聚类方法，其特征在于具体步骤如下；

首先，利用主集方法对基因或样本进行重新排序，使表达谱高度相似的基因排列在同一特定区域，由此确定一个类；其中，类的分割的定位标准是，根据类内元素间平均相信性高于类外元素间的平均相似性的性质，从排序后的基因或样本序列中划出具有高相似表达的基因或样本的类；然后，将找到的类从当前数据集中排除，对剩下的数据重复以上过程，直至满足所提出的循环条件为上。

2、根据权利要求1所述的聚类方法，其特征在于所述的利用主集方法对基因或样本进行重新排序的步骤如下：

(1)用G表示基因表达谱数据矩阵，G的每一行代表一种基因，每一列代表一种实验条件，G(i，k)表示基因i在条件k下的表达值，g_i和g_j分别表示第i和j个元素的数据向量，用a_ij表示对应两个元素i和j的相关性：

$a_{\mathrm{ij}}=\exp (-\frac{|g_i-g_j|}{\mathrm{\δ}})---\left(4\right)$

δ＞0，为缩放因子，记A＝(a_ij)表示G中基因之间相互关系；

对于一个带约束条件的极值问题，用离散动态方程表示为：

$x_i(t+1)=x_i\left(t\right)\frac{{\left(\mathrm{WX}\left(t\right)\right)}_i}{X^T\left(t\right)\mathrm{WX}\left(t\right)}---\left(2\right)$

其中x_i(t)是向量X(t)的分量，t是迭代步数；向量x的支持集可以定义为向量x中非零元素的下标集合，也就是：

               σ(x)＝{i∈V：x_i＞0}                          (3)

令公式(2)中的W＝A，求解该二次型问题，计算得向量x(t)，其每一分量x_i(t)对应原数据集G中的一条基因或样本向量，按照x_i(t)的值，从大到小对原始基因或样本进行排序，将排序后对应的基因表达向量依次记为g₁′，g₂′，...，g_n′；

(2)从重新排序之后的基因序列g₁′，g₂′，...，g_n′中，划分出一个类，具体步骤如下：

a.以欧式距离作为相似性的判断，先判断g₁′和g₂′之间的距离，假如g₁′和g₂′之间的距离大于一个初始设定的距离值，则停止分割，并记录分割失败的累计次数；若是在第(2)步中未找到合理聚类的连续次数大于某一数值，则停止整个算法的迭代过程；否则，将g₁′和g2′从剩余数据集中去初掉，退出步骤(2)返回步骤(1)；若g₁′和g₂′之间的距离小于上述初始设定的距离值，则将g₁′和g₂′作为初始类，并用g₁′和g₂′的平均值作为初始类中心；

b.接着判断g₃′是否属于这个类，若不属于，则停止分割；若属于，则将g₃′加入类中，并用相同的方法，反复迭代步骤b依次判断g₄′，g₅′等是否属于这一类，直到找到一个聚类，并输出该类元素为止；

(3)将第(2)步中找到的类从剩余的数据中去除掉，如果数据集中存在未被分类的元素，则返回第(1)步继续执行；若整个数据集已被完全分类，则停止迭代。

Images