专利名称：基因聚类程序、基因聚类方法及基因聚类分析装置的制作方法
技术领域：
本发明涉及基因聚类程序、基因聚类方法及基因聚类分析装置。更具体地，本发明涉及能够基于基因表达量随时间变化的相似性将每个基因归类为特定簇的基因聚类程序等。
背景技术：
在系统生物学的领域中，已经开始尝试基于基因表达量、基因定位和基因活性随时间变化的测量数据对由基因形成的细胞内信号网络进行解释。细胞内信号网络由动态变化的层次网络体系结构构成。近来已经有人提出了将 “蝴蝶结信号网络”作为构成细胞内信号网络的一种基本网络体系结构(非专利文献1和非专利文献2)。蝴蝶结信号网络(以下简称为“蝴蝶结网络”)具有好比蝴蝶结的网络体系结构， 蝴蝶结的结被想象成作为分类器的核分子，用于调节细胞对刺激物的免疫反应。也就是说， 在蝴蝶结网络中，细胞内和细胞间信号转导的各种输入聚集到置于结中的核分子中。于是， 核分子的细胞内浓度随着输入而发生改变，以根据浓度而激活位于信号下游的特定基因簇，从而表达特定的输出。已经报道了蝴蝶结网络可用于免疫细胞之间的信号转导、代谢信号转导(非专利文献1)、toll样受体信号转导(非专利文献i)和上皮生长因子信号转导(非专利文献 3)。已经研究了，蝴蝶结网络是一种卓越的网络体系结构，其坚固且具有用于进化的灵活性 (非专利文献4和非专利文献5)。在蝴蝶结网络中，位于信号下游的基因基于预定的核分子的浓度聚类成基因簇。 为了基于基因表达量、基因定位和基因活性随时间变化的测量数据来辨别各个基因所属于的簇并分析蝴蝶结网络，需要一种出色几何工具来解释整个网络体系结构从而预测簇间的关系。迄今，已基于k_均值法(非专利文献6)、层次聚类(非专利文献7)和自组织映射 (非专利文献8)研究了这样的工具。然而，这几种工具都有只能以一个步骤进行算术处理的缺点。S卩，由于簇被重叠以形成各数据要素的层级，所以层次聚类只作成不可改变的树状图。此外，层次聚类基于一对一的相似性来聚类基因，因此最终归到一个簇的基因可能互相不具有生物学关联。基于自组织映射(SOM)的工具(例如“基因簇(GENECLUSTER) ”)在数据的初步分析方面尤其出色，但需要为簇数量的预测初始值提前设定网格大小。传统的k-均值法同样需要提前设定簇数量，而且可能提供在生物学上没有意义的结果，因为聚类结果取决于所设定的数量。“GENEI^attern”(非专利文献9)通过横向整合这些传统工具而获得，是目前可用的最有效的工具。但是，它还不具备足够的性能来基于例如用来解释蝴蝶结网络的基因表达量随时间变化的数据正确地对每个基因进行聚类。
__专禾1J文献 1 :"The Edinburghhuman metabolic network reconstruction and its functional analysis，，，Molecular System Biology, 2007 ;3 :135。非专利文献 2:“A comprehensive map of the toll-like receptor signaling network",Molecular System Biology, 2006 ；2 :2006. 0015。非专利文献 3:“A comprehensive pathway map of epidermal growth factor receptor signaling，，，Molecular System Biology, 2005 ； 1 :2005. 0010。非专利文献 4 :"Bow ties, metabolism and disease，，，Trends in Biotechnology, 2004 ；22 (9) 446-50 # # ^lJ i K 5 :"Biological robustness", Nature Reviews Genetics, 2004 ； 5(11) :826-37o非专利文献 6 "‘Systematic determination of genetic network architecture", Nature Genetics,1999 ；22(3) :281_285。__专禾U文献 7 !"Cluster analysis and display of genome-wide expression patterns，，, Proceeding of National Academy of Sciences, 1998 ；95 (25) : 14863-14868。非专禾丨J 文献 8 ‘‘ Interpreting patterns of gene expression with self-organizing maps :Methods and application to hematopoietic differentiation", Proceeding of National Academy of Sciences,1999 ；96 (6) 2907-2912。非专利文献 9 :"GenePattern 2. 0”，Nature Genetics, 2006 ；38 :500_501。
发明内容
本发明要解决的问题因此，本发明的主要目的是提供基因聚类工具，其无需先验数据预测，就能够基于基因表达量随时间变化的数据进行高精度的基因聚类。解决问题的方法鉴于上述问题，本发明提供了一种基因聚类程序，用于至少执行以下步骤步骤 (1)，基于表示基因表达量随时间变化的数据计算反映数据间相似性的特征值；步骤0)， 基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量；步骤(3)，在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N ；以及步骤(4)，基于布尔矩阵N聚类数据。在该基因聚类程序中，在步骤中，通过线性回归分析或小波变换基于数据计算特征值。在步骤O)中，用核方法或余弦相似性基于特征值计算本征向量。此外，在步骤(3)中，用对称最近邻滤波(FSNN)算法将相似性矩阵M转换成布尔矩阵N。而且，在步骤(3)中，在用FSNN算法进行转换后，通过图形拉普拉斯(graph Laplacian)、马尔科夫链(Markov chain)、双随机近似(DSA)算法或双随机尺度(DSQ算法中的任一种将矩阵标准化。在该基因聚类程序中，在步骤中，通过最大期望(EM)算法和完全正分解(CP)算法进行软聚类。而且，在步骤(4)中，在软聚类后用布莱格曼-亚瑟-瓦斯尔维斯基初始化(Breg man-Arthur-Vassilvitskiiinitialization, BAV)算法进行硬聚类。本发明还提供了一种记录计算机可读的基因聚类程序的记录介质。本发明还提供了一种基因聚类方法，至少包括以下步骤步骤(1)，基于表示基因表达量随时间变化的数据计算反映数据间相似性的特征值；步骤O)，基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量；步骤(3)，在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N ；以及步骤(4)，基于布尔矩阵N聚类数据。此外，本发明还提供了一种基因聚类分析装置，至少包括装置(1)，用于基于表示基因表达量随时间变化的数据计算反映数据间相似性的特征值；装置O)，用于基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量；装置(3)，用于在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N ；以及装置(4)，用于基于布尔矩阵N聚类数据。本发明的效果本发明提供了基因聚类工具，其无需先验数据预测，就能够基于基因表达量随时间变化的数据进行高精度的基因聚类。
图1是示出了根据本发明的基因聚类程序中的处理步骤的流程图。图2是示出了表示基因表达量随时间变化的数据实例的示图，该数据由根据本发明的基因聚类程序进行处理。图3是示出了通过小波变换进行数据处理的概念图。图4是示出了用于创建基因表达量随时间变化的直方图的方法的概念图。图5是示出了计算特征值的步骤前后数据维数变化的概念图。图6是示出了基因i的对称最近邻的概念图。图7是示出了从相似性矩阵M到布尔矩阵N的转换过程的概念图。图8是示出了布尔矩阵和DSS矩阵的概念图。图9是示出了直到在根据本发明的基因聚类程序中获得最终聚类结果的过程中的数据处理的概念图。图10是示出了根据本发明的基因聚类分析装置的构造实例的框图。
具体实施例方式根据本发明的基因聚类方法至少执行以下步骤步骤(1)，基于表示基因表达量随时间变化的数据计算反映数据间相似性的特征值；步骤O)，基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量；步骤(3)，在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N;以及步骤G)，基于布尔矩阵N聚类数据。以下将对每个步骤逐一进行描述。1.特征值的计算该步骤相当于“基于表示基因表达量随时间变化的数据计算反映数据间相似性的特征值”的步骤⑴(见图1中的Si)首先，通过线性回归分析或小波变换(哈尔小波变换(Haar wavelet transform) 或多贝西小波变换(Daubechies wavelet transform))，利用D4-20尺度函数系数基于表示基因表达量随时间变化的数据计算反映数据间相似性的特征值。图2示出了表示基因表达量随时间变化的数据的实例。所示数据是对三个基因a、b、c在四个时间点1、2、3、4测量到的表达量。线性回归分析是用于对表示表达量变化的变化曲线进行比较的简单方法。而小波变换能够收集一定时间内变化曲线的所有信息。因此，小波变换甚至能够分析只在某个时间点提供表达数据的基因，而这在传统分析方法中会因为测量数据不完整从分析中被排除。图3示出了通过小波(哈尔小波变换)变换进行数据处理的概念图。在该小波变换中，利用直方图代替变化曲线对基因表达量随时间的变化数据(在这里，数据随时间从9、7、3到幻进行处理，该直方图被分解成例如一组四个的哈尔小波分量(见图3A)。数据用四维形式的平均数[9，7，3，5]、二维形式的平均数[8，4]和系数[1，_1]、一维形式的平均数[6]和系数[2]来表示。因此，通过一维小波转换数据被处理为[6(基)， 2，1，_1(系数)](见图；3)。小波转换利用直方图以这种方式对基因表达量的变化数据进行处理，从而能够用与使用变化曲线处理相比明显数量更少的系数进行最恰当的拟合。图4示出了基于图1所示表示基因表达量随时间变化的数据创建直方图的方法的概念图。图4A中实线或虚线所示的表达量变化可以转换成图4B所示直方图。在本步骤中，基因表达量随时间的变化数据被处理为以该方式转换成的直方图， 而特征值被计算为如上所述的一组系数以减少数据维数。图5示出了本步骤前后数据维数变化的示意图。2.相似性矩阵本征向量的计算接着，基于所算出的特征值通过核(热核)方法或余弦相似性对所有基因组合计算相似性矩阵M(半正定矩阵M)的本征向量。以下，相似性矩阵M简称为“矩阵M”。本步骤相当于根据本发明的基因聚类程序中“基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量”的步骤O)(见图1中S2)。(2-1)基于核方法的矩阵M当两个基因为i和j (i和j均是1以上的整数)，通过核方法在矩阵M中的行i和列j的输入被定义为式(1)。输入表示基因i和基因j之间的相似性。[数学式1]
权利要求
1.一种基因聚类程序，用于至少执行以下步骤步骤(1)，基于表示基因表达量随时间的变化的数据计算反映数据间相似性的特征值；步骤O)，基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量； 步骤( ，在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N ；以及步骤G)，基于布尔矩阵N聚类数据。
2.根据权利要求1所述的基因聚类程序，其中，在所述步骤(1)中，通过线性回归分析或小波变换基于所述数据计算所述特征值。
3.根据权利要求2所述的基因聚类程序，其中，在所述步骤O)中，利用核方法或余弦相似性基于所述特征值计算所述本征向量。
4.根据权利要求3所述的基因聚类程序，其中，在所述步骤(3)中，通过对称最近邻滤波(FSNN)算法将相似性矩阵M转换成布尔矩阵N。
5.根据权利要求4所述的基因聚类程序，其中，在所述步骤(3)中，在用对称最近邻滤波算法进行转换后，通过图形拉普拉斯、马尔科夫链、双随机近似(DSA)算法或双随机尺度 (DSS)算法中的任一种将矩阵标准化。
6.根据权利要求5所述的基因聚类程序中，其中，在所述步骤中，通过最大期望 (EM)算法和完全正分解(CP)算法进行软聚类。
7.根据权利要求6所述的基因聚类程序中，其中，在所述步骤(4)中，在软聚类后通过布莱格曼-亚瑟-瓦斯尔维斯基初始化(BAV)算法进行硬聚类。
8.—种记录介质，记录计算机可读的、根据权利要求1所述的基因聚类程序。
9.一种基因聚类方法，至少包括以下步骤步骤(1)，基于表示基因表达量随时间的变化的数据计算反映数据间相似性的特征值；步骤( ，基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量； 步骤( ，在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N ；以及步骤G)，基于布尔矩阵N聚类数据。
10.一种基因聚类分析装置，至少包括装置(1)，用于基于表示基因表达量随时间的变化的数据计算反映数据间相似性的特征值；装置O)，用于基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量； 装置(3)，用于在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N ；以及装置G)，用于基于布尔矩阵N聚类数据。
全文摘要
本发明要提供一种基因聚类工具，其无需先验数据预测，就能够基于基因表达量随时间变化的数据进行高精度的基因聚类。本发明提供了一种基因聚类程序，用于至少执行以下步骤步骤S1，基于表示基因表达量随时间的变化的数据计算反映数据间相似性的特征值；步骤S2，基于算出的特征值对所有基因组合计算相似性矩阵M的本征向量；步骤S3，在保持本征向量的本征值的同时将相似性矩阵M转换成布尔矩阵N；以及步骤S4，基于布尔矩阵N聚类数据。
文档编号G06F19/20GK102227731SQ200980147398
公开日2011年10月26日 申请日期2009年12月1日 优先权日2008年12月2日
发明者北野宏明, 弗兰克·尼尔森, 理查德·诺克, 纳塔利娅·波卢利亚赫 申请人:索尼公司