本发明涉及生物系统相变临界点检查的技术领域，具体涉及一种基于相对熵指标(relativeentropyscore,记为res)检测复杂生物系统相变临界点的方法。
背景技术：
糖尿病和癌症等复杂疾病的进展通常是一个非线性过程，具有三个阶段，即正常状态、疾病前状态和疾病状态，其中疾病前状态是疾病状态之前的临界状态或临界点。传统的生物标记物旨在通过利用观察到的分子的差异表达信息来识别疾病状态，但是由于正常和疾病前状态之间通常没有显著差异，所以可能无法检测到疾病前状态。因此，对疾病前状态进行信号通知是一个挑战，这实际上意味着疾病预测。
下面介绍计算方法的理论推导如下：
在临界相变之前和附近不同的动态特性：
复杂疾病进展的动力学可以通过以下非线性离散时间动态系统来表示：
z(t)＝f(z(t-1)；p),(1)
这里z(t)＝(z1(t),z2(t),…,zn(t))是n维状态向量或在时间t＝1,2,…的变量，p＝(p1,…,ps)是一个代表缓慢变化因素的参数向量或者驱动因素，例如，遗传因素(snp,cnv等)，表观遗传因素(methylation，acetylation等)或者环境因素。f:rn×rs×rn是一个非线性函数。对于这样一个非线性系统，该系统在
对于z附近的系统(1)，在p到达pc之前，系统应该保持稳定的平衡
·这组中的变量zi(t)之间的相关性增加；
·这组变量zi(t)与其他组的变量zj(t)之间的相关性下降；
·这组变量zi(t)的标准差增加。
因此，正常状态与疾病前状态之间的动力学存在显著差异。正常状态是一个具有高回弹的稳定状态，对参数扰动不敏感，因此可以将其建模为平稳的马尔科夫过程。当系统处于正常状态时，z(t)和z(t-1)的分布之间没有显著的变化，即随着时间的推移，概率分布几乎保持不变。相比之下，有着低回弹性的疾病前状态对参数变化敏感，它的动力学或概率分布随着时间的推移而变化。这样，疾病前状态被建模为一个时变的马尔科夫过程。当系统处于疾病前状态时，z(t)的分布和z(t-1)的分布有显著的差别。基于这些动态特性，可以识别从正常状态到疾病前状态的切换时间。
大多数生物分子通过与功能模块或模块之间的其它生物分子的相互作用来执行其功能。这种模块间和模块内的互连性表明，特定遗传异常的影响不仅影响携带它的基因产物的活性，而且可以沿着由生物分子组成的网络的链接延伸，并改变其它基因产物的活性。因此，了解生物分子的相互作用网络环境对于确定影响生物分子的缺陷的表型影响至关重要。
技术实现要素：
本发明的目的是通过利用正常状态和疾病状态之间的不同特性，提出一种基于相对熵指标(relativeentropyscore)检测复杂生物系统相变临界点的方法，在复杂疾病的生物过程中，在临界点到达之前识别出疾病前状态。特别地，识别疾病前状态就相当于检测两个网络发生差异的切换点。
为了研究网络系统的演化，本发明使用了一个差异网络，该差异网络集合了差异边，即量化差异网络中每条差异边的统计重要性(即相对熵指标，res)。
本发明的目的可以通过采取如下技术方案达到：
一种基于相对熵指标检测复杂生物系统相变临界点的方法，所述的方法包括如下步骤：
s1、把连续时间观测数据序列ot＝{o1,o2,o3,…,ot}转换为时序差异网络序列{dn2,dn3,…,dnt-1,dnt}；
先建立相关网络，通过将相关性映射到现有功能网络，即string网络，在每个采样时间点对观测序列{o1,o2,o3,…,ot}构造相关网络序列{n1,…,nt}，其中，nt表示t时刻的相关网络，连接两个节点的每条边表示两个生物分子之间的相关性，而仅连接一个节点的每条边表示生物分子的自调节或变化，随后，选择一个参数α，使得pearson相关系数pcc满足下式：|pcc|≥α，其中参数α是基于特定真实数据的待确定参数，把满足以上条件的pcc的相关网络的边保留下来，不满足以上条件的边除去，得到相关网络；
s2、准备参照样本，以正常时期提取的样本为参照样本。对于真实的数据集，我们通常选择来自正常组织的样本作为参照样本；
s3、根据参照样本，拟合生物分子的分布，具体如下：
对于生物分子gi，基于参照样本{s1，s2，…，sk}中的表达水平拟合高斯分布；然后，得到一个k维向量(area(dgi(s1)),area(dgi(s2)),…,area(dgi(sk)))，其中，area(dgi(sk))代表第k个样本中的生物分子gi被高斯分布所确定的累积面积；
s4、按照如下公式构造参照分布p
其中，
s5、计算相对熵指标(即，relativeentropyscore)，记为res
其中，res_n代表通过正常样本求得的相对熵指标，
其中，
h<u,v>(xv,xul)代表研究对象v与研究对象u的第l个正常样本之间的边特征，h<v,u>(xu,xvp)代表研究对象u与研究对象v的第p个正常样本之间的边特征，xv代表研究对象v的所有正常样本，xu代表研究对象u的所有正常样本，xul代表研究对象u的第l个正常样本，xvs代表研究对象v的第s个正常样本，而
其中，p(xv1)代表研究对象v的第1个正常样本的分布，p(xv2)代表研究对象v的第2个正常样本的分布，…，p(xvm)代表研究对象v的第m个正常样本的分布，p(xul)代表研究对象u的第l个正常样本的分布；
同理可得
其中，res_d代表通过疾病样本求得的相对熵指标，h<u,v>(yv,yul)代表研究对象v与研究对象u的第l个疾病样本之间的边特征，h<v,u>(yu,yvs)代表研究对象u与研究对象v的第s个疾病样本之间的边特征，yv代表研究对象v的所有疾病样本，yu代表研究对象u的所有疾病样本，yul代表研究对象u的第l个疾病样本，yvs代表研究对象v的第s个疾病样本，p代表离散概率分布，分布p满足
进一步地，该检测复杂生物系统相变临界点的方法需要至少3个样本。
进一步地，所述的相对熵指标(res)在不同状态下有不同的特征，疾病状态下相对熵指标(res)的值较正常状态下的值要小。
进一步地，所述的参数α的选取原则是使得正常状态下的差异网络存在尽量少的差异边，从而凸显出存在一定数量差异边的疾病前状态。
本发明相对于现有技术具有如下的优点及效果：
本发明给出了一个基于相对熵指标(res)的计算方法，用来识别即将到来的临界转变，它被真实数据集证明是有效的。值得注意的是，本发明的目的是检测从正常状态(或疾病前状态)产生的早期预警信号，而不是找到发生质变的疾病状态(或疾病前状态)的迹象。其中本发明的创新之处有：
1、传统的方法只能判断个体是处于健康状态还是疾病状态，而对于健康状态的极限状态也即临界转变关键期不能有效地察觉，本发明采用时间差异网络的计算方法，能够准确地反映出复杂疾病发展过程中的前疾病时期或预测复杂疾病恶化的发生；
2、现有技术中单变量或少数变量受噪声影响太大，临界点信号不明显，本发明方法可以克服；
3、本发明方法采用非监督学习以及前向算法，使高通量数据实现了实际运算；
4、本发明方法中非常巧妙的模型设计将连续的基因表达数据转换为观测数据，这是难点和关键点所在。
附图说明
图1是本发明公开的一种基于相对熵指标检测复杂生物系统相变临界点的方法流程示意图；
图2(a)是肺鳞状细胞癌(lusc)数据集上在第一种情况下，即：ⅰα期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图2(b)是肺鳞状细胞癌(lusc)数据集上在第二种情况下，即：ⅰα与ⅰβ期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图2(c)是肺鳞状细胞癌(lusc)数据集上在第三种情况下，即：ⅰα，ⅰβ与ⅱα期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图2(d)是肺鳞状细胞癌(lusc)数据集上在第四种情况下，即：ⅰα，ⅰβ，ⅱα与ⅱβ期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图2(e)是肺鳞状细胞癌(lusc)数据集上在第五种情况下，即：ⅰα，ⅰβ，ⅱα，ⅱβ与ⅲα期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图2(f)是肺鳞状细胞癌(lusc)数据集上基因表达和相对熵指标(res)结果的比较示意图；
图2(g)是肺鳞状细胞癌(lusc)数据集上由相对熵指标(res)构成的网络的动力学演化示意图；
图3(a)是肺腺癌(luad)数据集上在第一种情况下，即：ⅰ期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图3(b)是肺腺癌(luad)数据集上在第二种情况下，即：ⅰ与ⅰα期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图3(c)是肺腺癌(luad)数据集上在第三种情况下，即：ⅰ，ⅰα与ⅰβ期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图3(d)是肺腺癌(luad)数据集上在第四种情况下，即：ⅰ，ⅰα，ⅰβ与ⅱα期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图3(e)是肺腺癌(luad)数据集上在第五种情况下，即：ⅰ，ⅰα，ⅰβ，ⅱα与ⅱβ期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图3(f)是肺腺癌(luad)数据集上在第六种情况下，即：ⅰ，ⅰα，ⅰβ，ⅱα，ⅱβ与ⅲα期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图3(g)是肺腺癌(luad)数据集上在第七种情况下，即：ⅰ，ⅰα，ⅰβ，ⅱα，ⅱβ，ⅲα与ⅲβ期样本作为对照样本，剩余样本作为实验组，基因表达和相对熵指标(res)生存分析结果的比较示意图，左侧代表数据集上基因表达的生存分析结果，右侧为数据集上相对熵指标(res)的生存分析结果；
图3(h)是肺腺癌(luad)数据集上基因表达和相对熵指标(res)结果的比较示意图；
图3(i)是肺腺癌(luad)数据集上由相对熵指标(res)构成的网络的动力学演化示意图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
实施例
如图1所示，本发明所实施的案例公开了一种基于相对熵指标检测复杂生物系统相变前的临界状态的方法。
以下给出计算相对熵指标(res)的节点特征和边特征的数据矩阵示意。
首先是节点特征，用x,y来区分正常样本与疾病样本。研究对象u的第m个正常样本记为xum，同理研究对象u的第n个疾病样本记为yun。此外，研究对象v的第m个正常样本记为xvm，同理研究对象v的第n个疾病样本记为yvn。
接下来是边特征。
h<u,v>(xv,xul)代表研究对象v的正常样本与研究对象u的第l个正常样本之间的边特征；
h<v,u>(xu,xvp)代表研究对象u的正常样本与研究对象v的第p个正常样本之间的边特征；
h<u,v>(yv,yul)代表研究对象v的疾病样本与研究对象u的第l个疾病样本之间的边特征；
h<v,u>(yu,yvp)代表研究对象u的疾病样本与研究对象v的第p个疾病样本之间的边特征；
则通过正常样本求得的相对熵指标res_n：
其中，
xv代表研究对象v的所有正常样本，xu代表研究对象u的所有正常样本，xul代表研究对象u的第l个正常样本，xvs代表研究对象v的第s个正常样本，p(xv1)代表研究对象v的第1个正常样本的分布，p(xv2)代表研究对象v的第2个正常样本的分布，…，p(xvm)代表研究对象v的第m个正常样本的分布。p(xul)代表研究对象u的第l个正常样本的分布。
同理可得通过疾病样本求得的相对熵指标res_d
其中，h<u,v>(yv,yul)代表研究对象v与研究对象u的第l个疾病样本之间的边特征，h<v,u>(yu,yvp)代表研究对象u与研究对象v的第p个疾病样本之间的边特征，yv代表研究对象v的所有疾病样本，yu代表研究对象u的所有疾病样本，yul代表研究对象u的第l个疾病样本，yvs代表研究对象v的第s个疾病样本。
根据图1公开的流程示意图。
该实施例得到的结果如下：
1、预测真实数据集的临界点
本实施例将基于相对熵指标的方法应用于两个个真正的实验数据集即肺鳞状细胞癌(lusc)和肺腺癌(luad)。
2、基于相对熵指标在2个肿瘤数据集的应用
为了进一步证明该方法的有效性，将其应用于2个肿瘤数据集：肺鳞状细胞癌、肺腺癌，所有这些数据都来自于tcga癌基因图谱，由肿瘤和肿瘤邻近样本组成。根据tcga癌基因图谱相应的临床资料，将肿瘤分为不同的分期。肺鳞状细胞癌、肺腺癌可以分为7个时期。在所有2个数据集中，根据相对熵指标(res)算法计算每一对基因的相对熵指标(res)。最后，观察每对基因的相对熵指标(res)值的变化确定肿瘤的临界时期。
相对熵指标(res)成功地识别了所有两种癌症恶化前的关键阶段。为了验证所识别的临界时期，需对临界转化前和转换临界转化后的样本进行kaplan-meier(log-rank)生存分析进行比较(图2(a)-图2(e),图3(a)-图3(g))。临界转化前的样本的预后寿命通常比临界转化后的样本高。具体地说，对于肺鳞状细胞癌，从图2(c)可以看出，临界时期前样本(ia-iia期的样本)的存活时间比临界时期后样本(iib-iv阶段的样本)的存活时间长得多，两组样本生存曲线之间存在显著差异(显著值p＝0.042)。肺腺癌的ⅱb期前后样本的生存曲线有显著性差异(p＝0.015，图3(g))，临界前样本(ia-iib期的样本)的存活时间比临界时期后样本(iiia-iv期的样本)长得多。这些结果表明所确定的临界阶段是准确的，并且与预后密切相关。
综上所述，本发明利用观测到的正常和疾病前状态分子间的差异关联信息，提出了一种基于时间差异网络的计算方法，能够准确地反映疾病前状态或预测严重疾病的发生。这种差异网络与现有方法的区别在于，技术人员研究基因或者蛋白质的差异关联(或者相关性)，而不是基因或者蛋白质的差异表达。这一工作的理论基础是利用动态网络生物标志物对临界状态进行量化。
上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。
1.一种基于相对熵指标检测复杂生物系统相变临界点的方法，其特征在于，所述的方法包括如下步骤：
s1、把连续时间观测数据序列ot＝{o1，o2，o3，...，ot}转换为时序差异网络序列{dn2，dn3，...，dnt-1，dnt}；
先建立相关网络，通过将相关性映射到现有功能网络，即string网络，在每个采样时间点对观测序列{o1，o2，o3，...，ot}构造相关网络序列{n1，...，nt}，其中，nt表示t时刻的相关网络，连接两个节点的每条边表示两个生物分子之间的相关性，而仅连接一个节点的每条边表示生物分子的自调节或变化，随后，选择一个参数α，使得pearson相关系数pcc满足下式：|pcc|≥α，其中参数α是基于特定真实数据的待确定参数，把满足以上条件的pcc的相关网络的边保留下来，不满足以上条件的边除去，得到相关网络；
s2、准备参照样本，以正常时期提取的样本为参照样本。对于真实的数据集，我们通常选择来自正常组织的样本作为参照样本；
s3、根据参照样本，拟合生物分子的分布，具体如下：
对于生物分子gi，基于参照样本{s1，s2，...，sk}中的表达水平拟合高斯分布；然后，得到一个k维向量
s4、按照如下公式构造参照分布p
其中，
s5、计算相对熵指标，记为res
其中，res_n代表通过正常样本求得的相对熵指标，
其中，
h<u，v>(xv，xul)代表研究对象v与研究对象u的第l个正常样本之间的边特征，h<v，u>(xu，xvp)代表研究对象u与研究对象v的第p个正常样本之间的边特征，xv代表研究对象v的所有正常样本，xu代表研究对象u的所有正常样本，xul代表研究对象u的第l个正常样本，xvs代表研究对象v的第s个正常样本，而
其中，p(xv1)代表研究对象v的第1个正常样本的分布，p(xv2)代表研究对象v的第2个正常样本的分布，...，p(xvm)代表研究对象v的第m个正常样本的分布，p(xul)代表研究对象u的第l个正常样本的分布；
同理可得
其中，res_d代表通过疾病样本求得的相对熵指标，h<u，v>(yv，yul)代表研究对象v与研究对象u的第l个疾病样本之间的边特征，h<v，u>(yu，yvs)代表研究对象u与研究对象v的第s个疾病样本之间的边特征，yv代表研究对象v的所有疾病样本，yu代表研究对象u的所有疾病样本，yul代表研究对象u的第l个疾病样本，yvs代表研究对象v的第s个疾病样本，p代表离散概率分布，分布p满足
2.根据权利要求1所述的一种基于相对熵指标检测复杂生物系统相变临界点的方法，其特征在于，该检测复杂生物系统相变临界点的方法需要至少3个样本。
3.根据权利要求2所述的一种基于相对熵指标检测复杂生物系统相变临界点的方法，其特征在于，所述的相对熵指标在不同状态下有不同的特征，疾病状态下相对熵指标的值较正常状态下的值要小。
4.根据权利要求1所述的一种基于相对熵指标检测复杂生物系统相变临界点的方法，其特征在于，所述的参数α的选取原则是使得正常状态下的差异网络存在尽量少的差异边，从而凸显出存在一定数量差异边的疾病前状态。
技术总结