CN106202984A - 一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法 - Google Patents

Abstract

本发明提供的一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法，将随机森林模型和复杂网络中的节点介数结合在一起，用最小网格聚类方法验证，提供一种新的视角分析发现肿瘤致病因子和诊断标记物。通过生物信息学，数学统计分析，建立不同miRNA网络数据的相关性，用平均边覆盖和度分布比较miRNA表达量网络、miRNA结构网络和miRNA序列网络，从而分析整合了与肿瘤相关的不同miRNA的数据网络，在各种不同类型miRNA数据层之间建立起聚类特征的内部连结，以形成一个复杂的网络。这种不同种类生物分子之间的高度整合提供了一种更加简便，准确度更高的筛选方法，从而为癌症诊断和药物的发现提供更有价值的参考。

Description

一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法

技术领域

本发明涉及肿瘤标志物技术领域，具体涉及一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法。

背景技术

癌症是严重威胁人类生存和社会发展的重大疾病和严重的公共卫生问题之一，癌症控制已成为世界各国政府的卫生战略重点。MicroRNAs(miRNAs)是在真核生物中发现的一类内源性的具有调控功能的非编码RNA，其大小长约20～25个核苷酸。近些年，对miRNA研究越来越多，miRNA的表达水平与癌症的类型、分期及病人的其它临床数据相关，在癌症生物学的几乎所有方面都扮演角色，如增殖、凋亡、入侵、转移和血管生成。

在选择肿瘤标志物时，仅仅选择一种血清miRNA作为肿瘤标志物往往特异性较低，若测定多种miRNA组合或miRNA表达谱，可提高诊断的准确性。但是目前遇到的问题是癌症的复杂性。尽管系统分子生物学技术极大地提高了发现导致肿瘤恶化的miRNA表达谱异常的可能性，但产生的海量分子生物学数据很少被系统分析和利用。

尽管现有技术中有研究对癌症大数据和癌症miRNA表达谱进行数据挖掘，找出数据中蕴含的信息，例如，Wu等人综合比较了线性判别方法、二次判别分析，k近邻，bagging和boosting分类树(装袋和提高分类树)，支持向量机和随机森林方法在卵巢癌分类中的应用，随机森林方法优于其他方法(Comparison ofstatistical methodsforclassification ofovarian cancerusing mass spectrometry data.BaolinWu,Bioinformatics,2003,19(13):1636–1643.)，但是该方法使用比较复杂。另外，McDermott等人利用LuminimalA-like乳腺癌患者和对照组的血液样本，提取RNA，逆转录并进行微阵列分析确定76个差异表达miRNA，进一步通过神经网络数据挖掘算法识别10个差异表达miRNA，其中组合miR-29、miR-181a和miR-652经过二进制逻辑回归可较好区分癌症组和对照组，准确率0.8，经RQ-PCR检测这三个miRNA在乳腺癌中低表达(IdentificationandValidationofOncologic miRNABiomarkers for LuminalA-like BreastCancer.McDermottAM,PLoS ONE,2014,9(1):e87032.doi:10.1371/journal.pone.0087032)。但是人工神经网络需要大量的参数，如网络结构、权值和阈值的初始值。由于不能观察中间的学习过程，输出结果难以解释，会影响到结果的可信度和可接受程度，并且人工神经网络算法学习时间长，有时甚至可能达不到学习的目的。

发明内容

有鉴于此，本发明的目的在于通过对癌症大数据和癌症miRNA表达谱的数据挖掘，找出数据中蕴含的隐藏信息，找到miRNA序列、二级结构和表达量的关系，提供一种新的视角分析发现肿瘤诊断标记物。

为了实现上述发明目的，本发明提供以下技术方案：

本发明提供了一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法，包括以下步骤：

1)提供肿瘤病人的miRNA数据，所述miRNA数据包括肿瘤病人正常组织的miRNA数据和肿瘤组织的miRNA数据；

2)将所述步骤1)得到的正常组织的miRNA数据和肿瘤组织的miRNA数据分别进行过滤；

3)将所述步骤2)过滤后的miRNA数据分别做归一化处理，得到标准化的正常组织和肿瘤组织的miRNA数据；

4)将所述步骤3)得到的标准化的miRNA数据与下载的miRNA成熟序列进行序列比对，分别得到miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵；

5)将所述步骤4)得到的三个矩阵作为网络中节点间的距离，构建miRNA表达量网络、miRNA结构网络和miRNA序列网络；

6)比较所述步骤5)中的三个网络的平均边覆盖和度分布，得到相似二级结构的miRNA表达量相似；

7)采用随机森林模型对所述步骤2)中得到相同病人的正常组织的miRNA数据和肿瘤组织的miRNA的表达量选择最佳子集，计算最佳子集的miRNA的MIC值(MaximalInformation coefficient,MIC)；

8)基于所述步骤7)的最佳子集miRNA的MIC值分别构建正常组织和肿瘤组织的miRNA网络；

9)计算在正常组织和肿瘤组织的miRNA网络中同一节点的节点介数，比较并统计正常组织和肿瘤组织的网络中节点介数发生较大变化的点，从而筛选出肿瘤miRNA标志物的种类。

10)采用最小网格聚类方法对所述步骤4)得到的miRNA序列比对值矩阵进行聚类，将序列相似的miRNA聚为一类，根据已知的miRNA的功能推知具有相同序列的miRNA的功能，从而验证筛选出肿瘤miRNA标志物的种类是否准确。

优选的，所述步骤2)中过滤具体为挑选得到相同病人的正常组织的miRNA数据和肿瘤组织的miRNA数据去除miRNA表达量低于5％的数据。

优选的，所述步骤3)中归一化处理是采用Z标准化处理；Z标准化公式如式I所示为：

优选的，所述步骤5)中构建miRNA表达量网络、miRNA结构网络和miRNA序列网络的方法具体为：miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵，作为网络中miRNA结点间的距离，miRNA表达量网络分为三个部分构建，第一部分所选阈值为25，第二部分所选阈值为33，第三部分所选阈值为161；所述miRNA结构网络分为三部分构建，第一部分所选阈值为16，第二部分所选阈值为27，第三部分所选阈值为49；所述miRNA序列网络一次构建完成。

优选的，所述步骤6)中平均边覆盖是采用以下方法进行计算：将复杂网络之间的平均边覆盖<o>定义为：

$<o>=\frac{1}{k}\underset{i,j}{\&Sigma;}{o}_{l,j},o_{i,j}=\frac{1}{M}\underset{\mathrm{\&alpha;}}{\&Sigma;}{a}_{i,j}^{\&lsqb;\mathrm{\&alpha;}\&rsqb;}$

其中，k为至少在一个网络层中出现的结点对数；a_i,j＝1表示结点和结点之间有连边，a_i,j＝0表示结点和结点之间无连边；<0>＝1时表示当且仅当所有网络连边都相同。

优选的，所述步骤7)中选择最佳子集后还包括：十折交叉对数据进行验证。

优选的，所述步骤7)中选择最佳子集是采用随机森林模型分类器对miRNA表达量进行过滤和封装，对正常组织和肿瘤组织的miRNA表达量进行特征选择。

优选的，所述步骤7)中最佳子集miRNA的MIC值采用如下计算方法得到：

首先定义两个联合随机变量(X,Y)特征矩阵的分布规律的性质，然后使用G(k,l)表示；k,l均为正整数；联合随机变量(X,Y)分布在[0,1]×[0,1]中，定义网格G，使得(X,Y)|G＝(col_G(X),row_G(Y)),这里的col_G(X)表示网格G的X轴上的网格列数，row_G(Y)表示网格G的Y轴上的网格行数；

$I(X,Y)=\underset{y\&Element;Y}{\&Sigma;}\underset{x\&Element;X}{\&Sigma;}p(x,y)log\left(\frac{p(x,y)}{p\left(x\right)p\left(y\right)}\right)$

其中X,Y表示随机变量，p(x,y)表示联合概率分布，p(x)p(y)表示边际概率分布；

最大信息系数思想：把散点图中的数据用网格进行分区，封存在网格里，计算变量间的MIC值转化成网格的最佳分区；

用M(X,Y)表示(X,Y)的群体特征矩阵如式Ⅱ所示：

I^*((X,Y),s,t)＝maxI((X,Y)|G)

$I^{*}\left(\right(X,Y),s,t)=H\left(x\right)+H\left(y\right)-H(x,y)=\underset{y\&Element;Y}{\&Sigma;}\underset{x\&Element;X}{\&Sigma;}p(x,y)\log \left(\frac{p(x,y)}{p\left(x\right)p\left(y\right)}\right)$

s·t＜B(n)＝n^0.6

其中：G∈G(k,l)，I(X,Y)，代表的交互信息；n表示数据的数量，s,t分别表示轴和轴分区的数量；

所述的MIC值计算公式为MIC＝max_s.t＜B(n)M(X,Y)_s,t，。

优选的，所述计算节点介数的方法是采用Matlab程序运算。

优选的，所述最小网格聚类方法具体为：

定义：x设为任一miRNA，Y为miRNA集合，则D(x,Y)＝min_y∈YD(x,y)为到的距离；记[x]为包含miRNA：x的类，Γ为所有miRNA的集合；最小网格聚类方法如下：

(1)求(x₁,x₂)＝argmax_x,yD(x,y)，设x₁,y₁为两个初始类，CL₁＝[x₁]，CL₂＝[x₂]。C₁＝{CL₁,CL₂}，R₁＝Γ\C₁，集合R₁为剩下的miRNA集合；

(2)n步以后，可得C_n＝{CL₁,CL₂,…,CL_n+1}，R_n＝Γ\C_n，其中CL_i为单点类，i＝1,2,…,n+1，在n+1步，设为剩余miRNA与集合距离，若MD>ε，ε为设定的最小网格阈值，需找设C_n+1＝{C_n,x_MD}，R_n+1＝Γ\C_n+1，直到满足MD＜ε；

(3)假设上述步骤在N步停止，有C_N＝{CL₁,CL₂,…,CL_N+1}，对于任意x∈R_N的，寻找与其距离最近的类，若y_k∈CL_k则是CL_k类；

定义：设CL_i为任一个聚类，则，为聚类CL_i的半径。

本发明提供的一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法，通过构建正常组织和肿瘤组织的miRNA网络，得到两个网络结构差异较大的节点介数，从而筛选出肿瘤miRNA标志物的种类。将随机森林模型和复杂网络中结点重要性度量—介数结合在一起，提供一种新的视角分析发现致病因子或诊断标记物。筛选微小RNA分子组学等系统生物组学数据，构建了由多变量组成的癌症分子信息数据结构网络。通过生物信息学，数学统计，建立不同分子网络数据的相关性，从而整合了与肿瘤相关的不同数据层的数据，在各种不同类型生物学数据层之间建立起高度的内部连结，从而形成一个复杂的网络。这种不同种类生物分子之间的高度整合，从而提供了一种方法更加简便，准确度高的筛选方法，从而为癌症诊断和药物的发现提供有价值的参考。

说明书附图

图1为本发明实施例1中miRNA表达量网络、miRNA二级结构网络和miRNA比对序列网络；

图2为本发明实施例2中miRNA表达量网络和二级结构网络度分布情况；

图3为本发明实施例2中乳腺癌病人正常组织和肿瘤组织miRNA表达量的MIC网络；

图4为本发明实施例3中前列腺病人正常组织和肿瘤组织miRNA表达量的MIC网络。

具体实施方式

本发明提供了一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法，包括以下步骤：

1)提供肿瘤病人的miRNA数据，所述miRNA数据包括肿瘤病人正常组织的miRNA数据和肿瘤组织的miRNA数据；

2)将所述步骤1)得到的正常组织的miRNA数据和肿瘤组织的miRNA数据分别进行过滤；

3)将所述步骤2)过滤后的肿瘤组织的miRNA数据分别做归一化处理，得到标准化的肿瘤组织的miRNA数据；

4)将所述步骤3)得到的标准化的miRNA数据与下载的miRNA成熟序列进行序列比对，分别得到miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵；

5)将所述步骤4)得到的三个矩阵作为网络中节点间的距离，构建miRNA表达量网络、miRNA结构网络和miRNA序列网络；

6)比较所述步骤5)中的三个网络的平均边覆盖和度分布，发现相似二级结构的miRNA，它们的表达量也相似；

7)采用随机森林模型对所述步骤2)中得到相同病人的正常组织的miRNA数据和肿瘤组织的miRNA的表达量选择最佳子集，计算最佳子集的miRNA的MIC值；

8)基于所述步骤7)的最佳子集miRNA的MIC值分别构建正常组织和肿瘤组织的miRNA网络；

9)计算在正常组织和肿瘤组织的miRNA网络中的同一节点的节点介数，比较并统计正常组织和肿瘤组织的网络中节点介数发生较大变化的点，从而筛选出肿瘤miRNA标志物的种类。

10)采用最小网格聚类方法对所述步骤4)得到的miRNA序列比对值矩阵进行聚类，将序列相似的miRNA聚为一类，根据已知的miRNA的功能推知具有相同序列的miRNA的功能，从而验证筛选出肿瘤miRNA标志物的种类是否准确。

本发明提供的一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法，将随机森林模型和复杂网络中结点重要性度量—介数结合在一起，提供一种新的视角分析发现致病因子或诊断标记物。

首先，提供肿瘤病人的miRNA数据，所述miRNA数据包括肿瘤病人正常组织的miRNA数据和肿瘤组织的miRNA数据。

本发明中，所述肿瘤病人的miRNA数据的来源为癌症基因信息数据库The CancerGenomeAtlas/TCG(https://tcga-data.nci.nih.gov/tcga/dataAccessMatrix.h tm)。

得到的相同病人的正常组织的miRNA数据和肿瘤组织的miRNA数据后，本发明对得到的正常组织的miRNA数据和肿瘤组织的miRNA数据分别进行过滤。

本发明中，所述过滤具体优选为挑选得到相同病人的正常组织的miRNA数据和肿瘤组织的miRNA数据中去除miRNA表达量低于5％的数据。

得到的过滤后的miRNA数据后，本发明为减少后续数据的差距，对所述过滤后的肿瘤组织的miRNA数据分别做归一化处理，得到标准化的肿瘤组织的miRNA数据。

本发明中，所述的归一化处理优选采用Z标准化处理进行；所述的Z标准化公式如式I所示为：

得到的标准化的肿瘤组织的miRNA数据后，本发明对得到的标准化肿瘤组织的miRNA数据与下载的miRNA成熟序列进行序列比对，分别得到miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵。本发明中，所述的miRNA成熟序列的来源优选为miRBase(http://www.mirbase.org)数据库下载得到。

本发明中，将标准化的肿瘤组织的miRNA数据按照下载的miRNA成熟序列的数量调整，使两者的数据数量相同。

本发明中，所述miRNA表达水平之间的欧式距离矩阵优选将肿瘤组织的miRNA数据与miRNA成熟序列表达水平之间的欧氏距离形成欧式距离矩阵。

本发明中，所述二级结构间的差异矩阵优选将标准化的肿瘤组织的miRNA序列与miRNA成熟序列的二级结构间的差异数据形成矩阵。

本发明中，所述miRNA的序列比对值矩阵优选是将标准化的肿瘤组织的miRNA序列与miRNA成熟序列的miRNA的序列比对值形成矩阵。miRNA序列比对优选用Needlemen-Wunsch算法求得的matlab函数(nwalign)。

得到miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵后，本发明根据所述得到三个矩阵分别构建miRNA表达量网络、miRNA结构网络和miRNA序列网络。

本发明中，所述miRNA表达水平之间的欧式距离矩阵的数值点作为网络中节点间的距离，构建miRNA表达量网络。本发明中，所述miRNA表达量网络分为三个部分构建，第一部分所选阈值为25，第二部分所选阈值为33，第三部分所选阈值为161。

本发明中，所述二级结构间的差异矩阵的数值点作为网络中节点间的距离，构建miRNA结构网络。本发明中，所述miRNA结构网络分为三部分构建，第一部分所选阈值为16，第二部分所选阈值为27，第三部分所选阈值为49。

本发明中，所述miRNA序列比对值矩阵的数值点作为网络中节点间的距离，构建miRNA序列网络。本发明中，所述miRNA序列网络优选一次构建完成，所选阈值为31.2。

得到miRNA表达量网络、miRNA结构网络和miRNA序列网络后，本发明对所述miRNA表达量网络、miRNA结构网络和miRNA序列网络通过计算两者平均边覆盖和度分布的趋势，具有相似二级结构的miRNA，它们的表达量也相似。

本发明中，所述平均边覆盖具体技术方法优选为将复杂网络之间的平均边覆盖(average edge overlap)<o>定义为：

$<o>=\frac{1}{k}\underset{i,j}{\&Sigma;}{o}_{l,j},o_{i,j}=\frac{1}{M}\underset{\mathrm{\&alpha;}}{\&Sigma;}{a}_{i,j}^{\&lsqb;\mathrm{\&alpha;}\&rsqb;}$

其中，k为至少在一个网络层中出现的结点对数；a_i,j＝1表示结点和结点之间有连边，a_i,j＝0表示结点和结点之间无连边；<0>＝1时表示当且仅当所有网络连边都相同。

本发明中，所述miRNA表达量网络和结构网络的平均边覆盖和度分布相似，得到具有相似二级结构的miRNA它的表达量也相似的结论。

基于具有相似二级结构的miRNA它的表达量也相似的结论，采用随机森林模型对所述过滤后的正常组织的miRNA数据和肿瘤组织的miRNA的表达量选择最佳子集，计算正常组织的miRNA数据和肿瘤组织miRNA数据的最佳子集的miRNA数据间的MIC值。

本发明中，所述随机森林模型具体优选为采用随机森林模型分类器进行选择。所述选择最佳子集是采用随机森林模型分类器对miRNA表达量进行过滤和封装，对正常组织和肿瘤组织的miRNA表达量进行特征选择。所述过滤和封装具体优选采用R语言caret包中sbf和rfe方法进行特征选择，以降低维度。

本发明中，所述选择最佳子集后优选还包括采用十折交叉对数据进行验证，以检测最佳子集的选择的准确率。所述最佳子集的选择经十折交叉验证准确率≥98.50％时才能进入后续分析。

所述十折交叉的具体方法没有特殊限制，采用本领域技术人员所熟知的十折交叉的技术方法即可。

本发明中，所述计算最佳子集miRNA的MIC值具体是最佳子集中每个miRNA之间的MIC值。

本发明中，所述最佳子集miRNA的MIC值优选采用如下计算方法得到：

首先定义两个联合随机变量(X,Y)特征矩阵的分布规律的性质，然后使用G(k,l)表示；k,l均为正整数；联合随机变量(X,Y)分布在[0,1]×[0,1]中，定义网格G，使得(X,Y)|G＝(col_G(X),row_G(Y)),这里的col_G(X)表示网格G的X轴上的网格列数，row_G(Y)表示网格G的Y轴上的网格行数；

$I(X,Y)=\underset{y\&Element;Y}{\&Sigma;}\underset{x\&Element;X}{\&Sigma;}p(x,y)\log \left(\frac{p(x,y)}{p\left(x\right)p\left(y\right)}\right)$

其中X,Y表示随机变量，p(x,y)表示联合概率分布，p(x)p(y)表示边际概率分布；

最大信息系数思想：把散点图中的数据用网格进行分区，封存在网格里，计算变量间的MIC值转化成网格的最佳分区；

用M(X,Y)表示(X,Y)的群体特征矩阵如式Ⅱ所示：

I^*((X,Y),s,t)＝maxI((X,Y)|G)

$I^{*}\left(\right(X,Y),s,t)=H\left(x\right)+H\left(y\right)-H(x,y)=\underset{y\&Element;Y}{\&Sigma;}\underset{x\&Element;X}{\&Sigma;}p(x,y)log\left(\frac{p(x,y)}{p\left(x\right)p\left(y\right)}\right)$

s·t＜B(n)＝n^0.6

其中：G∈G(k,l)，I(X,Y)代表的交互信息；n表示数据的数量，s,t分别表示轴和轴分区的数量；

所述的MIC值计算公式为MIC＝max_s.t＜B(n)M(X,Y)_s,t，可知：0≤MIC≤1。对于B(n)＝n^α中的α＝0.6，n表示数据的数量，涉及探索最大上限网格B(n)的划分。B(n)设置太大将导致随机数据的MIC值非零，每个点将落到自身的小空格里；B(n)设置过小将导致MIC只能检测出明显和简单的函数关系。所述B(n)的大小依赖n的大小。MIC-P²:反映非线性关系的程度。其中P为pearson相关系数，由于MIC代表一般性的相关关系，MIC-P²值越大，非线性程度越高。

基于所述正常组织miRNA的最佳子集的MIC值和肿瘤组织miRNA的最佳子集的MIC值分别构建正常组织miRNA网络和肿瘤组织的miRNA网络。

本发明中，所述构建正常组织miRNA网络和肿瘤组织的miRNA网络的具体方法优选为每两个miRNA之间的MIC值作为网络中两个结点边的权重，两个网络阈值都选0.35，互信息大于0.35的两个节点之间有连边。

计算正常组织和肿瘤组织的miRNA网络中同一节点的节点介数，比较并统计正常组织和肿瘤组织的网络中节点介数发生较大变化的点，从而筛选出肿瘤miRNA标志物的种类。

本发明中，所述计算节点介数的方法优选采用Matlab程序运算。

采用最小网格聚类方法对所述步骤4)得到的miRNA序列比对值矩阵进行聚类，将序列相似的miRNA聚为一类，根据已知的miRNA的功能推知具有相同序列的miRNA的功能，从而验证筛选出肿瘤miRNA标志物的种类是否准确。

本发明中，所述步骤10)中最小网格聚类方法具体优选为：

定义：x设为任一miRNA，Y为miRNA集合，则D(x,Y)＝min_y∈YD(x,y)为到的距离；记[x]为包含miRNA：x的类，Γ为所有miRNA的集合；最小网格聚类方法如下：

(1)求(x₁,x₂)＝argmax_x,yD(x,y)，设x₁,y₁为两个初始类，CL₁＝[x₁]，CL₂＝[x₂]。C₁＝{CL₁,CL₂}，R₁＝Γ\C₁，集合R₁为剩下的miRNA集合；

(2)n步以后，可得C_n＝{CL₁,CL₂,…,CL_n+1}，R_n＝Γ\C_n，其中CL_i为单点类，i＝1,2,…,n+1，在n+1步，设为剩余miRNA与集合距离，若MD>ε，ε为设定的最小网格阈值，需找设C_n+1＝{C_n,x_MD}，R_n+1＝Γ\C_n+1，直到满足MD＜ε；

(3)假设上述步骤在N步停止，有C_N＝{CL₁,CL₂,…,CL_N+1}，对于任意x∈R_N的，寻找与其距离最近的类，若y_k∈CL_k则是CL_k类；

定义：设CL_i为任一个聚类，则，为聚类CL_i的半径。

下面结合实施例对本发明提供的一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法进行详细的说明，但是不能把它们理解为对本发明保护范围的限定。

实施例1

研究数据的来源为癌症基因信息数据库(The Cancer Genome Atlas/TCGA)(https://tcga-data.nci.nih.gov/tcga/dataAccessMatrix.htm)。选取浸润性乳癌病人miRNASeq数据。其中，miRNA的数据为1034个数，来自1206个病人。其中在miRNA表达数据中有101个为乳腺癌病人的正常组织miRNA数据，其余为乳腺癌病人肿瘤组织miRNA数据。在正常组织和肿瘤组织的miRNA数据中，有很多miRNA不表达或表达率低，去除miRNA表达量低于5％的个体得到正常组织和肿瘤组织的miRNA数据。从miRBase(http://www.mirbase.org)数据库下载收集540个miRNA的成熟序列。

对肿瘤组织的miRNA表达水平做归一化处理，对肿瘤组织的miRNA数据进行Z标准化，从而缩小构建矩阵的差距。将选择540个Z标准化的肿瘤组织的miRNA数据与540个miRNA的成熟序列进行比对，分别得到miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵。将三个矩阵的数值点作为网络中节点间的距离分别构建miRNA表达量网络、miRNA结构网络和miRNA序列网络。

用平均边覆盖和度分布比较miRNA表达量网络、miRNA结构网络和miRNA序列网络。如附图1所示为miRNA表达量、miRNA二级结构和miRNA序列网络，其中A)为miRNA表达量网络，分为三个部分共有3025条边，6个孤立点；B)为miRNA二级结构网络，三部分共有2798条边，2个孤立点；C)miRNA序列网络，共有2897条边，4个孤立点。miRNA表达水平网络与结构网络和序列网络的平均边覆盖分别为0.507和0.5314，三个网络的总边数大致相同，结构网络有2798个总边数，表达量网络有3025个总边数，比对序列值网络有2897个总边数，三个网络的边仅约占完全图边数的0.02。所占比例较小的情况下，表达水平网络与结构网络的边覆盖已经超过50％。统计两个网络的重要属性，结果如表1所示，这两个网络的度分布函数大致服从幂分布，如图2A-B所示，用Matlab进行曲线拟合，得到的曲线如图2C-D所示，得到miRNA表达量网络的度分布函数为f(x)＝0.3618x^-0.08919-0.2588，SSE为0.002626，R-square为0.9013。miRNA结构网络的度分布函数f(x)＝0.199x^-0.2244-0.08338，SSE为0.001386，R-square为0.9481。miRNA表达量网络和结构网络的度分布相似，大致服从幂律分布。miRNA的表达量网络和序列网络边覆盖为0.5314，与结构网络的边覆盖为0.507，有理由相信在癌症网络研究中miRNA的表达量和二级结构是有关联的。

表1miRNA表达量网络和miRNA结构网络属性

实施例2

将实施例1下载得到的乳腺癌病人的正常组织和肿瘤组织的miRNA数据，提取乳癌病人正常组织和肿瘤组织的miRNA，病人数101个，miRNA数据为622个。采用随机森林模型对乳腺癌病人正常组织的miRNA数据和肿瘤组织的miRNA进行封装和过滤选择筛选择最佳子集。为了选择具有最小基因数量的并保持最高分类准确率的miRNA子集，采用十折交叉验证进行试验以评估分类器模型，得到miRNA分类结果如表2所示。对于乳腺癌数据集，当miRNA子集数量为50时，得到交叉验证准确率为98.50％。筛选出这50个miRNA，分别计算病人正常组织和肿瘤组织的50个miRNA之间的互信息(MIC)值，根据所得到的值构建MIC网络，附图3为正常组织和肿瘤组织的MIC网络。每两个miRNA之间的互信息作为网络中两个结点边的权重，两个网络阈值都选0.35，互信息大于0.35的两个节点之间有连边。一个miRNA在两个网络中，节点介数差别大，说明此miRNA在正常组织和肿瘤组织间发生了较大变化，也许在癌症发生发展中起重要作用。

提取10个节点介数差异较大的miRNA，hsa-mir-101-2、hsa-mir-10b、hsa-mir-130b、hsa-mir-190b、hsa-mir-193a、hsa-mir-204、hsa-mir-28、hsa-mir-365-2、hsa-mir-375、hsa-mir-192。查阅文献资料可知，其中hsa-mir-101-2侧翼区rs462480和rs1053872遗传突变可增加患乳腺癌的风险。抑制hsa-mir-10的表达可部分逆转由于TGF-b1诱导的乳腺癌细胞的EMT过程。hsa-mir-130b被验证在三阴性乳腺癌中可直接抑制细胞周期蛋白G2基因(CCNG2--一种调节细胞周期的关键基因)。BMP-6(骨形态发生蛋白)已被确认为乳腺癌分化和转移相关的肿瘤抑制剂，miR-192表达量升高可抑制细胞生长,在乳腺癌组中miR-192表达量明显下降，与BMP-6表达量成正比，可调控miR-192表达发挥BMP-6对细胞增殖的抑制作用。hsa-mir-193a作为新的肿瘤抑制剂作用于EGFR驱动细胞周期，抑制乳腺肿瘤细胞周期进程和增值、hsa-mir-204与乳腺癌的TMN分期和转移相关，mir-204低表达的患者比高表达的患者存活时间短。众所周知，BCL2蛋白是的调节线粒体生理和细胞死亡的基因，miR-195，miR-24-2和miR-365-2通过直接结合BCL2基因的3'-UTR结合位点负调控BCL2。miR-195，miR-24-2和miR-365-2的异位表达可导致BCL2蛋白水平的显著减少。hsa-mir-375在乳腺小叶增生差异表达并促进乳腺腺泡极性的损失。其中hsa-mir-101-2、hsa-mir-10b、hsa-mir-193a、hsa-mir-204、hsa-mir-28、hsa-mir-365-2在乳腺癌表达水平中下调。hsa-mir-130b、hsa-mir-190b、hsa-mir-375、hsa-mir-192腺癌表达水平中上调。

表2采用随机森林分类器所获得的miRNA十折交叉验证准确率

从miRBase查询miRNA成熟序列并做序列比对，得到一个540×540的距离矩阵，利用最小网格聚类，将这些miRNA做聚类分析。这个聚类算法把序列相似的miRNA聚为一类，具有相似结构的miRNA，功能也相似，其中如表3所示的第七类含有23个miRNA。其中hsa-let-7家族与乳腺癌相关，hsa-let-7i，hsa-mir-34a，hsa-mir-34b，hsa-mir-520b，hsa-mir-708，hsa-mir-548l与癌症转移相关，hsa-mir-28，hsa-mir-376a-2在正常组织与肿瘤组织的网络中，前后介数相差较大，在这两个网络中节点重要性发生很大变化，可推测也与乳癌转移相关。

类52中含有三个miRNA，hsa-mir-210，hsa-mir-331，hsa-mir-940，其中hsa-mir-210，hsa-mir-940都与细胞侵袭相关。hsa-mir-331在正常组织与肿瘤组织的网络中，前后介数相差较大，可推测也与细胞侵袭相关。

如表4所示的类14中含有14个元素，其中hsa-let-7家族与乳腺癌相关，hsa-mir-1256抑制前列腺癌细胞的增值和扩散，hsa-mir-146a可以诱导细胞凋亡,hsa-mir-195，hsa-mir-31可以抑制癌细胞的增值，扩散和转移。hsa-mir-1323参与人肺癌细胞的抗辐射。miR-192的过表达显著诱导膀胱癌细胞、食管癌细胞has-mir-1307在正常组织与肿瘤组织的网络中，前后介数相差较大，可推测也与细胞扩散和增值相关。

由此得之，随机森林模型和复杂网络分析方法选取癌症的重要生物分子方法是可靠的，可为癌症的预测、治疗和药物的选择与研发提供有价值的参考。

表3第7类miRNA的种类

表4第14类miRNA的种类

实施例3

选取与乳腺癌不相关的另种癌症--前列腺癌，采用实施例1和实施2的方法下载序列并处理，对随机森林模型和复杂网络分析方法选取乳癌发生重要生物分子方法进行验证。当miRNA子集数量为30时，得到95.27％的交叉验证准确率。构建前列腺病人正常组织和肿瘤组织miRNA表达量的MIC网络，见附图4，图4-A)图为前列腺癌病人正常组织miRNA表达量的MIC网络，节点越大，颜色越深表明节点介数越大；图4-B)前列腺癌病人肿瘤组织miRNA表达量的MIC网络，节点越大，颜色越深表明节点介数越大。

计算MIC网络中节点的节点介数，选取8个节点介数差异较大的miRNA，hsa-mir-183、hsa-mir-19a、hsa-mir-221、hsa-mir-25、hsa-mir-375、hsa-mir-381、hsa-mir-92a-1、hsa-mir-96。其中七个都与前列腺癌相关，具体信息见表5，除了hsa-mir-381。

表5与前列腺癌相关miRNA功能

由以上实施例可知，本发明提供的一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法，采用随机森林模型和复杂网络分析方法选取癌症的重要生物分子方法是可靠的，可为癌症的预测、治疗和药物的选择与研发提供有价值的参考。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (10)

1.一种基于多层复杂网络对肿瘤miRNA标志物的筛选方法，其特征在于，包括以下步骤：

1)提供肿瘤病人的miRNA数据，所述miRNA数据包括肿瘤病人正常组织的miRNA数据和肿瘤组织的miRNA数据；

2)将所述步骤1)得到的正常组织的miRNA数据和肿瘤组织的miRNA数据分别进行过滤；

3)将所述步骤2)过滤后的肿瘤组织的miRNA数据分别做归一化处理，得到标准化的肿瘤组织的miRNA数据；

4)将所述步骤3)得到的标准化的miRNA数据与下载的miRNA成熟序列进行序列比对，分别得到miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵；

5)将所述步骤4)得到的三个矩阵分别作为网络中节点间的距离，分别构建miRNA表达量网络、miRNA结构网络和miRNA序列网络；

6)比较所述步骤5)中的三个网络的平均边覆盖和度分布，得到具有相似二级结构的miRNA，所述相似二级结构的miRNA的表达量也相似；

7)采用随机森林模型对所述步骤2)中得到相同病人的正常组织的miRNA数据和肿瘤组织的miRNA的表达量选择最佳子集，计算最佳子集的miRNA的MIC值；

8)基于所述步骤7)的最佳子集miRNA的MIC值分别构建正常组织和肿瘤组织的miRNA网络；

9)计算在正常组织和肿瘤组织的miRNA网络中同一节点的节点介数，比较并统计正常组织和肿瘤组织的网络中节点介数发生较大变化的点，从而筛选出肿瘤miRNA标志物的种类。

10)采用最小网格聚类方法对所述步骤4)得到的miRNA序列比对值矩阵进行聚类，将序列相似的miRNA聚为一类，根据已知的miRNA的功能推知具有相同序列的miRNA的功能，从而验证筛选出肿瘤miRNA标志物的种类是否准确。

2.根据权利要求1所述的方法，其特征在于，所述步骤2)中过滤具体为从挑选的相同病人的正常组织的miRNA数据和肿瘤组织的miRNA数据中去除miRNA表达量低于5％的数据。

3.根据权利要求1所述的方法，其特征在于，所述步骤3)中归一化处理是采用Z标准化处理；Z标准化公式如式I所示：

4.根据权利要求1所述的方法，其特征在于，所述步骤5)中构建miRNA表达量网络、miRNA结构网络和miRNA序列网络的方法具体为：miRNA表达水平之间的欧式距离矩阵、二级结构间的差异矩阵和miRNA的序列比对值矩阵，作为网络中miRNA结点间的距离，miRNA表达量网络分为三个部分构建，第一部分所选阈值为25，第二部分所选阈值为33，第三部分所选阈值为161；所述miRNA结构网络分为三部分构建，第一部分所选阈值为16，第二部分所选阈值为27，第三部分所选阈值为49；所述miRNA序列网络一次构建完成。

5.根据权利要求1所述的方法，其特征在于，所述步骤6)中平均边覆盖是采用以下方法进行计算：将复杂网络之间的平均边覆盖<o>定义为：

$<o>=\frac{1}{k}\underset{i,j}{\&Sigma;}{o}_{l,j},o_{i,j}=\frac{1}{M}\underset{\mathrm{\&alpha;}}{\&Sigma;}{a}_{i,j}^{\&lsqb;\mathrm{\&alpha;}\&rsqb;}$

其中，k为至少在一个网络层中出现的结点对数；a_i,j＝1表示结点和结点之间有连边，a_i,j＝0表示结点和结点之间无连边；<0>＝1时表示当且仅当所有网络连边都相同。

6.根据权利要求1所述的方法，其特征在于，所述步骤7)中选择最佳子集后还包括：十折交叉验证对数据进行验证。

7.根据权利要求1所述的方法，其特征在于，所述步骤7)中选择最佳子集是采用随机森林模型分类器对miRNA表达量进行过滤和封装，对正常组织和肿瘤组织的miRNA表达量进行特征选择。

8.根据权利要求1所述的方法，其特征在于，所述步骤7)中最佳子集miRNA的MIC值采用如下计算方法得到：

首先定义两个联合随机变量(X,Y)特征矩阵的分布规律的性质，然后使用G(k,l)表示；k,l均为正整数；联合随机变量(X,Y)分布在[0,1]×[0,1]中，定义网格G，使得(X,Y)|G＝(col_G(X),row_G(Y)),所述col_G(X)表示网格G的X轴上的网格列数，row_G(Y)表示网格G的Y轴上的网格行数；

$I(X,Y)=\underset{y\&Element;Y}{\&Sigma;}\underset{x\&Element;X}{\&Sigma;}p(x,y)log\left(\frac{p(x,y)}{p\left(x\right)p\left(y\right)}\right)$

其中X,Y表示随机变量，p(x,y)表示联合概率分布，p(x)p(y)表示边际概率分布；

最大信息系数思想：把散点图中的数据用网格进行分区，封存在网格里，计算变量间的MIC值转化成网格的最佳分区；

用M(X,Y)表示(X,Y)的群体特征矩阵如式Ⅱ所示：

I^*((X,Y),s,t)＝maxI((X,Y)|G)

$I^{*}\left(\right(X,Y),s,t)=H\left(x\right)+H\left(y\right)-H(x,y)=\underset{y\&Element;Y}{\&Sigma;}\underset{x\&Element;X}{\&Sigma;}p(x,y)log\left(\frac{p(x,y)}{p\left(x\right)p\left(y\right)}\right)$

s·t＜B(n)＝n^0.6

其中：G∈G(k,l)，I(X,Y)，代表的交互信息；n表示数据的数量，s,t分别表示轴和轴分区的数量；

所述的MIC值计算公式为MIC＝max_s.t＜B(n)M(X,Y)_s,t，。

9.根据权利要求1所述的方法，其特征在于，所述计算节点介数的方法是采用Matlab程序运算。

10.根据权利要求1所述的方法，其特征在于，所述步骤10)中最小网格聚类方法具体为：

定义：x设为任一miRNA，Y为miRNA集合，则D(x,Y)＝min_y∈YD(x,y)为到的距离；记[x]为包含miRNA：x的类，Γ为所有miRNA的集合；最小网格聚类方法如下：

(1)求(x₁,x₂)＝argmax_x,yD(x,y)，设x₁,y₁为两个初始类，CL₁＝[x₁]，CL₂＝[x₂]。C₁＝{CL₁,CL₂}，R₁＝Γ\C₁，集合R₁为剩下的miRNA集合；

(2)n步以后，可得C_n＝{CL₁,CL₂,…,CL_n+1}，R_n＝Γ\C_n，其中CL_i为单点类，i＝1,2,…,n+1，在n+1步，设为剩余miRNA与集合距离，若MD>ε，ε为设定的最小网格阈值，需找设C_n+1＝{C_n,x_MD}，R_n+1＝Γ\C_n+1，直到满足MD＜ε；

(3)假设上述步骤在N步停止，有C_N＝{CL₁,CL₂,…,CL_N+1}，对于任意x∈R_N的，寻找与其距离最近的类，若y_k∈CL_k则是CL_k类；

定义：设CL_i为任一个聚类，则，为聚类CL_i的半径。