CN109841265B - 使用片段化模式确定血浆游离核酸分子组织来源的方法和***及应用 - Google Patents

Abstract

本发明提供了一种使用片段化模式确定血浆游离核酸分子组织来源的方法和***及应用。该方法通过将核酸分子的测序结果与多个参考序列进行比对处理，分别确定与所述参考序列对应的已知基因的频谱参数；然后基于预定组织中多个已知基因的表达量信息，以及多个已知基因的所述频谱参数，确定所述核酸分子与所述预定组织的相关性；其中，所述频谱参数是通过确定参考序列上位点的核小体保护分数，确定多个核小体周期作为第一参数以及与所述多个核小体周期对应的傅里叶强度作为第二参数，从而确定所述参考序列对应所述已知基因的所述频谱参数。将本发明所提供的方法应用于癌症的检测中，能够更加灵敏和准确。

Description

使用片段化模式确定血浆游离核酸分子组织来源的方法和系 统及应用

技术领域

本发明涉及生物信息领域，具体涉及一种使用片段化模式确定血浆游离核酸分子组织来源的方法和***及应用。

背景技术

癌症，即恶性肿瘤和血液癌症的统称，是细胞发生不可控增殖、且伴随转移的恶性疾病。癌症的发病率居高不下，晚期癌症的生存率较低，使得癌症的早发现、早诊断、早治疗尤为重要。

液体活检技术，不同于术中取固体肿瘤组织进行化验检测，是通过抽血即可安全、便捷、准确地检测癌症的一种方法。其主要原理为：血浆游离DNA(cfDNA)来源于人体各组织器官凋亡后释放的片段化基因组。当癌症发生时，cfDNA中即含有癌细胞的基因组片段。另外，cfDNA的半衰期约为十几分钟到几个小时，这使得抽血监测癌症成为可能。

但抽血进行癌症监测，预测组织病变的准确性还需要进一步提高。

发明内容

本发明旨在至少在一定程度上解决相关技术中的技术问题之一。为此，本发明的一个目的在于提出一种确定核酸分子与预定组织相关性的方法以及确定核酸分子与预定组织相关的***、确定血浆游离核酸组织来源的方法。

本发明是基于如下发现所获得的：

细胞凋亡时的基因组片段化模式一般被认为与其染色质结构有关。片段化模式，顾名思义，是指细胞内的核酸酶将DNA打断为长度不一的片段时产生的片段端点模式。染色质排列开放的基因组区域比紧密的区域更易结合核酸酶从而被打断；不结合任何蛋白质的裸露DNA更容易被核酸酶打断，而被核小体、转录因子等保护的区域更不容易被打断。因此，由基因组的片段化模式可以推得各基因区域的核小体排布。核小体在基因组上的排布具有极强的周期性，根据现有MNase-seq的研究，相邻核小体的间距约为185bp，包含约145bp的核心缠绕区域和约40bp的连接区域。

利用公共数据集中各组织细胞的基因表达数据可推得基因启动子区域的核小体排布。具体来说，基因表达越高，其启动子区域的核小体排布约紧密，反之，基因表达越低，其启动子区域的核小体排布约松散。利用cfDNA的片段化模式筛查癌症。目前在癌症液体活检领域中，多数研究关注cfDNA中的癌症特异突变、拷贝数变异和甲基化模式，2015年以来，开始涌现利用cfDNA片段化模式检测癌症的研究方法。因为不同组织细胞的染色质结构不同，在发生凋亡时基因组的片段化模式不同，所以我们通过检测cfDNA中的片段化模式推得混合核小体排布，结合公共数据集中各组织细胞的基因表达数据推得的核小体排布，可解耦出cfDNA的组织来源分布，从而对受试者是否含有癌组织进行筛查。

为此，本发明提供了如下技术方案：

根据本发明的第一方面，本发明提供了一种确定核酸分子与预定组织相关性的方法，包括：将所述核酸分子的测序结果与多个参考序列进行比对处理，其中，所述测序结果由多个测序读段构成，所述多个参考序列的每一个分别是基于多个已知基因之一的序列确定的；基于所述比对处理的结果，针对所述多个参考序列的每一个，分别确定与所述参考序列对应的所述已知基因的频谱参数；基于所述预定组织中多个所述已知基因的表达量信息，以及多个所述已知基因的所述频谱参数，确定所述核酸分子与所述预定组织的相关性；其中，所述频谱参数是针对所述参考序列的每一个，分别通过下列步骤确定的：针对所述参考序列的至少一部分位点，分别确定每个位点的核小体保护分数；基于所述每个位点的核小体保护分数，确定第一参数和第二参数，所述第一参数表征多个核小体周期，所述第二参数表征分别与所述多个核小体周期对应的傅里叶强度；基于所述第一参数和所述第二参数，确定所述参考序列对应所述已知基因的所述频谱参数。

通过本发明提供的确定核酸分子与预定组织相关性的方法，能够判断核酸分子与预定组织的相关性。通过研究不同个体中相同核酸分子与预定组织的相关性，能够判断不同的状态。例如，能够用于人群的癌症早期筛查，还可用于癌症病人的复发风险预测与监控，干预癌症病人的诊疗。而且不仅可以用于检测多种癌症，可检测的癌症类型取决于搜集到哪些组织的表达谱数据。

根据本发明的实施例，以上所述确定核酸分子与预定组织相关性的方法可以进一步包括如下技术特征：

在本发明的一些实施例中，所述核酸分子为血浆游离核酸分子。

在本发明的一些实施例中，所述参考序列长度为5000～11000bp。参考序列的长度应大于25个核小体周期，长度过短会造成第二参数估计不准确，而参考序列过长则会引入与基因表达不相关的区域，从而降低方法的准确性。

在本发明的一些实施例中，所述参考序列包括所述已知基因的序列。

在本发明的一些实施例中，所述参考序列是基于所述已知基因的转录起始位点确定的。基因的转录起始位点附近富含启动子和转录因子结合位点，其核小体特征与基因的表达量密切相关，因此选取转录起始位点附近区域作为参考序列有利于确定核酸分子与预定组织的相关性。

在本发明的一些实施例中，所述测序是通过高通量测序进行的。

在本发明的一些实施例中，所述已知基因包括Ensembl第75版中GRCh37的所有基因。

在本发明的一些实施例中，针对每个位点，所述核小体保护分数是通过下列步骤确定的：

(1)基于所述位点，确定比对窗口，所述比对窗***有所述位点；

(2)确定与所述比对窗口匹配的所述测序读段；

(3)将所述测序读段区分为第一类型测序读段和第二类型测序读段，所述第一类型测序读段完全覆盖所述比对窗口，所述第二类型测序读段不完全覆盖所述比对窗口；

(4)基于所述第一类型测序读段的数目和所述第二类型测序读段的数目，确定所述核小体保护分数。

在本发明的一些实施例中，所述比对窗口长度为100～140。窗口大小应与核小体核心组蛋白结合DNA的长度一致，一般选择为120bp，过短会将转录因子的保护位点误认为是核小体保护位点，过长则会丢失部分排布紧密的核小***点。

在本发明的一些实施例中，所述位点位于所述比对窗口的中部。

在本发明的一些实施例中，在步骤(4)中，所述核小体保护分数与所述第一类型测序读段的数目呈正相关，与所述第二类型测序读段的数目呈负相关。

在本发明的一些实施例中，在步骤(4)中，所述核小体保护分数是通过公式第一类型测序读段的数目-所述第二类型测序读段确定的。

在本发明的一些实施例中，基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数，确定所述参考序列对应的多个核小体周期作为第一参数，以及所述多个核小体周期对应的傅里叶强度作为第二参数。

在本发明的一些实施例中，所述第一参数和第二参数是通过傅里叶变换确定的。

在本发明的一些实施例中，所述第一参数和第二参数是通过下列步骤确定的：

(a)基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数进行平滑矫正后作图得到第一图，其中，x轴为位点与基因转录起始位点的距离，y轴为矫正后的核小体保护分数；

(b)将步骤步骤(a)中得到的所述第一图进行傅里叶变换，得到第二图，所述第二图的x轴为核小体周期，y轴为与核小体周期对应的傅里叶强度。

在本发明的一些实施例中，所述频谱参数是基于下列公式确定的:

其中，T代表核小体周期，F(T)代表核小体周期T对应的傅里叶强度，a代表核小体周期的下限，b代表核小体周期的上限。

在本发明的一些实施例中，a＝168，b＝208。该区间是通过最大化健康人核酸分子与血液组织的相关性得到的，代表着与基因表达最相关的核小体周期范围，其他可行参数包括：a∈[120,170]，b∈[205,212]。在本发明的一些实施例中，在确定参数a、b后，所述频谱参数仅与所述第一、第二参数有关。

在本发明的一些实施例中，所述相关性是通过下列公式确定的：

其中，ρ_i代表所述核酸分子与预定组织i的皮尔逊相关系数，

代表所述核酸分子在多个已知基因的所述频谱参数，A_i代表预定组织i在所述多个已知基因的表达量，Var代表方差，Cov代表协方差。

在本发明的一些实施例中，所述预定组织选自：肝、胆、脾、肺、肾、膀胱、食管、胃、小肠、结肠、直肠、十二指肠、阑尾、胰岛、唾液腺、扁桃体、甲状腺、甲状旁腺、肾上腺、乳腺、卵巢、输卵管、子宫、子宫内膜、宫颈、***、***、睾丸、附睾、精囊、皮肤、脂肪组织、大脑皮层、眼、心肌、骨骼肌、平滑肌、***、血液、骨髓中的至少一种。

根据本发明的第二方面，本发明提供了一种确定核酸分子与预定组织相关性的***，包括：比对装置，所述比对装置用于将所述核酸分子的测序结果与多个参考序列进行比对处理，其中，所述测序结果由多个测序读段构成，所述多个参考序列的每一个分别是基于多个已知基因之一的序列确定的；频谱参数确定装置，所述频谱参数确定装置与所述比对装置相连，所述频谱参数确定装置针对所述多个参考序列的每一个，分别确定与所述参考序列对应的所述已知基因的频谱参数；相关性确定装置，所述相关性确定装置与所述频谱参数确定装置相连，所述相关性确定装置基于所述预定组织中多个所述已知基因的表达量信息，以及多个所述已知基因的所述频谱参数，确定所述核酸分子与所述预定组织的相关性；其中，所述频谱参数确定装置针对所述参考序列的每一个，分别通过如下单元确定的：核小体保护分数确定单元，所述核小体保护分数确定单元针对所述参考序列的至少一部分位点，分别确定每个位点的核小体保护分数；参数确定单元，所述参数确定单元与所述核小体保护分数确定单元相连，所述参数确定单元基于所述每个位点的核小体保护分数，确定第一参数和第二参数，所述第一参数表征多个核小体周期，所述第二参数表征分别与所述多个核小体周期对应的傅里叶强度；频谱参数计算单元，所述频谱参数计算单元与所述参数确定单元相连，所述频谱参数计算单元基于所述第一参数和所述第二参数，确定所述参考序列对应所述已知基因的所述频谱参数。

根据本发明的实施例，以上所述确定核酸分子与预定组织相关性的***可以进一步包括如下技术特征：

在本发明的一些实施例中，以上***中，所述核酸分子为血浆游离核酸分子；

在本发明的一些实施例中，以上***中，所述参考序列长度为5000～11000bp。

在本发明的一些实施例中，以上***中，所述参考序列包括所述已知基因的序列。

在本发明的一些实施例中，以上***中，所述参考序列是基于所述已知基因的转录起始位点确定的。

在本发明的一些实施例中，以上***中，所述测序是通过高通量测序进行的。在本发明的一些实施例中，以上***中，所述第一参数和第二参数是通过傅里叶变换确定的。

在本发明的一些实施例中，以上***中，所述第一参数和第二参数是通过下列步骤确定的：

(a)基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数进行平滑矫正后作图得到第一图，其中，x轴为位点与基因转录起始位点的距离，y轴为矫正后的核小体保护分数；

(b)将步骤步骤(a)中得到的所述第一图进行傅里叶变换，得到第二图，所述第二图的x轴为核小体周期，y轴为与核小体周期对应的傅里叶强度。

在本发明的一些实施例中，所述频谱参数是基于下列公式确定的:

其中，T代表核小体周期，F(T)代表核小体周期T对应的傅里叶强度，a代表核小体周期的下限，b代表核小体周期的上限；在本发明的一些实施例中，以上***中，a＝168，b＝208。

在本发明的一些实施例中，以上***中，所述相关性是通过下列公式确定的：

其中，ρ_i代表所述核酸分子与预定组织i的皮尔逊相关系数，

代表所述核酸分子在多个已知基因的所述频谱参数，A_i代表预定组织i在所述多个已知基因的表达量，Var代表方差，Cov代表协方差。

在本发明的一些实施例中，以上***中，针对每个位点，所述核小体保护分数是通过下列步骤确定的：

(1)基于所述位点，确定比对窗口，所述比对窗***有所述位点；

(2)确定与所述比对窗口匹配的所述测序读段；

(3)将所述测序读段区分为第一类型测序读段和第二类型测序读段，所述第一类型测序读段完全覆盖所述比对窗口，所述第二类型测序读段不完全覆盖所述比对窗口；

(4)基于所述第一类型测序读段的数目和所述第二类型测序读段的数目，确定所述核小体保护分数。

在本发明的一些实施例中，以上***中，所述比对窗口长度为100～140bp。

在本发明的一些实施例中，以上***中，所述位点位于所述比对窗口的中部。

在本发明的一些实施例中，以上***中，在步骤(4)中，所述核小体保护分数与所述第一类型测序读段的数目呈正相关，与所述第二类型测序读段的数目呈负相关。

在本发明的一些实施例中，以上***中，在步骤(4)中，所述核小体保护分数是通过公式第一类型测序读段的数目-所述第二类型测序读段的数目确定的。

在本发明的一些实施例中，以上***中，基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数，确定所述参考序列对应的多个核小体周期作为第一参数，以及所述核小体周期对应的傅里叶强度作为第二参数。

根据本发明的第三方面，本发明提供了一种确定血液无细胞核酸组织来源的方法，包括：(i)根据本发明第一方面任一实施例所述的方法，确定所述血液无细胞核酸与预定组织的相关性；(ii)将所述相关性进行排序，以便确定所述血浆游离核酸的组织来源。

根据本发明的实施例，以上确定血液无细胞核酸组织来源的方法可以进一步包括如下技术特征：

在本发明的一些实施例中，将步骤(ii)中得到的排名与对照排名进行比较，以便确定排名变化的组织来源。

在本发明的一些实施例中，所述对照排名是基因正常个体或者具有已知状态的个体确定的。

在本发明的一些实施例中，所述排序是按照由高到低进行的，并且选择排序升高的组织作为重点监控组织。

本发明所取得的有益效果为：与现有利用cfDNA片段化模式筛查癌症的技术相比，本发明的泛癌种筛查的液体活检分析方法能更充分利用cfDNA片段化模式信息，与基因表达的相关性更大，在测序较深时能更灵敏、准确地检测多种癌症。

附图说明

图1为根据本发明的实施例提供的技术原理图，展示了本发明首创的频谱加权核小体周期特征与现有技术的关键差异。

图2为根据本发明的实施例提供的相关性散点图，表明了本发明使用的频谱加权核小体周期特征与基因表达的相关性比现有技术更高。

图3为根据本发明的实施例提供的测序深度与相关性的关系图，证明了本发明使用的频谱加权核小体周期特征与基因表达的相关性随着测序深度的增加而增大，现有技术则无此现象。

图4为根据本发明的实施例提供的使用本发明与现有技术检测49个癌症样本的结果。当测序较浅时，两种方法均不足以用于筛查；当测序较深时，使用本发明能更灵敏地检测出真实癌症组织。

图5为根据本发明的实施例提供的使用本发明与现有技术进行下采样分析的结果。

图6为根据本发明的实施例提供的确定核酸分子与预定组织相关性的***的结构示意图。

具体实施方式

下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的，旨在用于解释本发明，而不能理解为对本发明的限制。

确定核酸分子与预定组织相关性的方法以及确定核酸分子与预定组织相关性的 ***

根据本发明的一个方面，本发明提供了一种确定核酸分子与预定组织相关性的方法，包括：将所述核酸分子的测序结果与多个参考序列进行比对处理，其中，所述测序结果由多个测序读段构成，所述多个参考序列的每一个分别是基于多个已知基因之一的序列确定的；基于所述比对处理的结果，针对所述多个参考序列的每一个，分别确定与所述参考序列对应的所述已知基因的频谱参数；基于所述预定组织中多个所述已知基因的表达量信息，以及多个所述已知基因的所述频谱参数，确定所述核酸分子与所述预定组织的相关性；其中，所述频谱参数是针对所述参考序列的每一个，分别通过下列步骤确定的：针对所述参考序列的至少一部分位点，分别确定每个位点的核小体保护分数；基于所述每个位点的核小体保护分数，确定第一参数和第二参数，所述第一参数表征多个核小体周期，所述第二参数表征分别与所述多个核小体周期对应的傅里叶强度；基于所述第一参数和所述第二参数，确定所述参考序列对应所述已知基因的所述频谱参数。

本文中，确定核酸分子与预定组织相关性指的是：通过对核酸分子的片段化模式进行处理分析，获得所述核酸分子的核小体特征与不同组织基因表达的相关性的大小，然后通过比较这些相关性的大小，可以用来判断核酸分子与这些组织的关系。例如，可以判断这些核酸分子来自于这些组织的可能性；例如，可以通过比较这些核酸分子与不同组织的相关性，判断来自于哪些组织的成分更多。

本文中，术语“参考序列”指的是能够通过比对来表征所述核酸分子的基因或者基因组来源，长度等各方面信息的核酸序列。根据本发明的实施例，所述参考序列可以是核酸序列已知的基因或者基因组。通过将所述核酸分子的测序结果对多个参考序列进行比对处理，借助于参考序列的信息，来获得所述核酸分子的基因组位置、长度等片段化模式信息，这些信息可以用来指示不同的核酸分子。这些参考序列可以是以基因组或者基因等形式存在。例如，参考序列可以是已知基因的序列。参考序列的长度可以根据基因的不同，长度可以是5000～11000bp，例如可以是10000bp。参考序列可以是从已知基因的转录起始位点开始，向下游延伸到5000～11000bp左右。这些已知基因可以是Ensembl第75版中GRCh37的所有基因。

其中，术语“片段化模式”指的是核酸分子所呈现出的片段化模式，例如在哪些位置发生断裂从而片段化等等。例如，细胞内的核酸酶能够将DNA打断成长度不一的片段，从而产生不同的片段化模式。例如，染色质排列开放的基因组区域比紧密的区域更易结合核酸酶从而被打断。再例如，不结合任何蛋白质的裸露DNA更容易被核酸酶打断，而被核小体、转录因子等保护的区域更不容易被打断。通过对这些片段化模式信息进行分析，可以用来指示不同的核酸分子。

根据本发明的实施例，核酸分子可以是游离于细胞外的核酸分子。例如可以为血浆游离DNA(cfDNA)。血浆游离DNA来源于人体各组织器官凋亡后释放的片段化基因组。通过对这些血浆游离DNA中片段化基因组的片段化模式信息分析，就可以用来指示不同的核酸分子。这些血浆游离DNA中的核酸分子，来自于人体的多种组织。通过对不同核酸分子的片段化模式信息分析，可以确定核酸分子与不同组织的相关性，从而可以分析血浆游离DNA中核酸分子与组织的相关性的强弱。

在获取核酸分子的测序结果时，可以通过高通量测序来进行。在将核酸分子进行测序获得测序结果时，可以通过对待测样本的核酸序列进行测序文库制备、上机测序获得。根据本发明的实施例，在制备测序文库时可以根据所选择的测序方法的要求进行，测序方法依据所选的测序平台的不同，可选择但不限于Illumina公司的Hisq2000/2500测序平台、Life Technologies公司的Ion Torrent平台、BGI的BGISEQ平台，测序方式为双末端测序，获得的下机数据是测读出来的片段，称为读段(reads)。在制备测序文库时，可以利用DNA单链全基因组建库、常规的DNA双链全基因组建库、pre-BS的全基因组甲基化建库等等。

本文中，“将所述核酸分子的测序结果与多个参考序列进行比对处理”中所称的“比对”即匹配的含义。具体比对时，可以利用已知比对软件进行，例如bowtie2、SOAP、BWA和TeraMap等，本申请对此不作限制。在比对过程中，根据比对参数的设置，一对或一条读段最多允许有n个碱基错配(mismatch)，例如设置n为1或2，若读段中有超过n个碱基发生错配，则视为该对读段无法比对到参考序列，或者，若错配的n个碱基全部位于读段对中的一个读段，则视为该读段对中的该读段无法比对到参考序列。同时可以提供samtools、bedtools等软件提取DNA片段端点、长度等信息。

本文中，术语“核小体保护分数”即WPS(window protection score)，表征基因组各个位点被核小体保护的程度，这个值越大表示该位点存在核小体的可能性越大。WPS值受到该位点比对上的核酸分子数量的影响，随机效应较大，可通过平滑、减均值两步处理进行矫正，再使用矫正后WPS做后续分析。

本文中，术语“频谱参数”即频谱加权核小体周期SWP(spectral weightedperiod)，是由核小体保护分数的波形图通过傅里叶变换转换到频谱空间做周期分解得到的，表示一段基因组区域内的核小体平均间距。这个值越大表示区域内的核小体排列越疏松，反之，SWP越小表示核小体排列越紧密。

根据本发明的实施例，每个位点的核小体保护分数可以通过下列步骤确定：(1)基于所述位点，确定比对窗口，所述比对窗***有所述位点；(2)确定与所述比对窗口匹配的所述测序读段；(3)将所述测序读段区分为第一类型测序读段和第二类型测序读段，所述第一类型测序读段完全覆盖所述比对窗口，所述第二类型测序读段不完全覆盖所述比对窗口；(4)基于所述第一类型测序读段的数目和所述第二类型测序读段的数目，确定所述核小体保护分数。

针对每个位点，比对窗口的长度可以为100～140bp。例如，可以是100bp～130bp，110bp～120bp等等。每个位点位于比对窗口的中部。例如，当选择所述比对窗口大约是120bp时，可以选择位点上游60bp和位点下游60bp作为相应的比对窗口。而针对每个位点的比对窗口，确定能够完全覆盖所述比对窗口的测序读段作为第一类型的测序读段，确定不完全覆盖所述比对窗口的测序读段作为第二类型的测序读段。所谓覆盖，指的是测序读段的5’末端和3’末端与比对窗口对齐或者在比对窗口以外。完全覆盖即测序读段的5’末端和3’末端在比对窗口之外或者与比对窗口对齐。不完全覆盖即测序读段的5’末端落在比对窗口内，或者测序读段的3’末端落在比对窗口内。根据本发明的实施例，每个位点的核小体保护分数与所述第一类型测序读段的数目呈正相关，与所述第二类型测序读段的数目呈负相关。在至少一些实施方式中，可以利用第一类型测序读段的数目与第二类型测序读段的数目的差值来作为每个位点的核小体保护分数。

在一种具体实施方式中，利用如下方法计算位点的核小体保护分数：对基因组每任一位点i，考虑以i为中心(2d+1)长的窗口[i-d,i+d]内，假设某片段比对到基因组后的端点位置为[lo,hi]，若[lo,hi]与[i-d,i+d]的交集为[i-d,i+d]，则认为片段完全覆盖该窗口；反之，若[lo,hi]与[i-d,i+d]的交集不等于[i-d,i+d]且不为空，则认为片段不完全覆盖该窗口。若完全覆盖该窗口的片段有m条，不完全覆盖该窗口的片段有n条，则WPS＝m-n。依此类推，可计算基因组任一位置的WPS值。

进一步地，可以根据参考序列的位点信息以及所述位点信息对应的核小体保护分数，确定所述参考序列对应的多个核小体周期作为第一参数，所述核小体周期对应的傅里叶强度作为第二参数。本文中，每个位点对应一个核小体保护分数，如上所述每个位点对应的核小体保护分数可能相同，也可能不同，在一段基因组区域内呈现为高高低低的波形图。将核小体保护分数波形图中的峰值视为核小***点，则相邻核小体的间距即为核小体周期。在一段基因组区域内，将核小体保护分数波形图做傅里叶变换，可理解为将疏密不均的多个核小体进行周期分解，由此可得到多个核小体周期及其对应的周期占比，术语“第一参数”用来表征这里的多个核小体周期。术语“第二参数”用来表征与所述多个核小体周期对应的傅里叶强度，即周期占比。

在至少一些实施方式中，可以通过傅里叶变换确定相应的第一参数和第二参数。在至少一些实施方式中，所述第一参数和第二参数可以通过下列步骤确定：(a)基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数进行作图得到第一图，其中，x轴为位点与基因转录起始位点的距离，y轴为矫正后的核小体保护分数；(b)将步骤步骤(a)中得到的所述第一图进行傅里叶变换，得到第二图，所述第二图的x轴为核小体周期，y轴为核小体周期对应的傅里叶强度。

同时，本发明还提供了一种确定核酸分子与预定组织相关性的***，如图6所示包括：比对装置、频谱参数确定装置、相关性确定装置，其中所述频谱参数确定装置与所述比对装置相连，所述相关性确定装置与所述频谱参数确定装置相连，所述比对装置用于将所述核酸分子的测序结果与多个参考序列进行比对处理，其中，所述测序结果由多个测序读段构成，所述多个参考序列的每一个分别是基于多个已知基因之一的序列确定的；所述频谱参数确定装置针对所述多个参考序列的每一个，分别确定与所述参考序列对应的所述已知基因的频谱参数；所述频谱参数确定装置针对所述参考序列的每一个，分别通过如下单元确定所述频谱参数：核小体保护分数确定单元，所述核小体保护分数确定单元针对所述参考序列的至少一部分位点，分别确定每个位点的核小体保护分数；参数确定单元，所述参数确定单元与所述核小体保护分数确定单元相连，所述参数确定单元基于所述每个位点的核小体保护分数，确定第一参数和第二参数，所述第一参数表征多个核小体周期，所述第二参数表征分别与所述多个核小体周期对应的傅里叶强度；频谱参数计算单元，所述频谱参数计算单元与所述参数确定单元相连，所述频谱参数计算单元基于所述第一参数和所述第二参数，确定所述参考序列对应所述已知基因的所述频谱参数。

上述对于确定核酸分子与预定组织相关性的方法中的优点以及技术特征的描述，同样适用于本发明所提供的确定核酸分子与预定组织相关性的***，在此不再赘述。

利用本发明所提供的确定核酸分子与预定组织相关性的方法或者确定核酸分子与预定组织相关性的***，来确定核酸分子与预定组织的相关性，可以用于分析预定组织的状态，例如，预定组织是否发生癌变等异常。

在一些实施方式中，本发明提供了一种确定血浆游离核酸组织来源的方法，包括：(i)根据本发明确定核酸分子与预定组织相关性的方法，确定所述血浆游离核酸与预定组织的相关性；(ii)将所述相关性进行排序，来确定所述血浆游离核酸的组织来源。通过对正常个体或者具有已知状态的个体进行分析，得到核酸分子与预定组织相关性的排名，作为对照排名。然后将上述步骤(ii)中所获得的排名与对照排名进行比较，确定排名变化的组织来源。在一些优选实施方式中，所述排序是按照相关性从高到低进行的，并且选择排序升高的组织作为重点监控组织。

以cfDNA为例，与现有利用cfDNA片段化模式筛查癌症的技术相比，利用本发明提供的方法来检测组织的异常，本发明所提供的方法能更充分利用cfDNA片段化模式信息，而且通过核小体周期作为第一参数，与多个核小体周期对应的傅里叶强度作为第二参数，关联第一参数和第二参数，获得核酸分子对应已知基因的频谱参数，与基因表达的相关系数更大，在测序较深时能更灵敏、准确地检测多种癌症。

下面将结合实施例对本发明的方案进行解释。本领域技术人员将会理解，下面的实施例仅用于说明本发明，而不应视为限定本发明的范围。实施例中未注明具体技术或条件的，按照本领域内的文献所描述的技术或条件或者按照产品说明书进行。所用试剂或仪器未注明生产厂商者，均为可以通过市购获得的常规产品。

实施例1

实施例1以cfDNA为例，提供了一种通过抽血获取cfDNA中有关频谱参数的检测方法，具体如下：

对检测受试者而言，其只需要抽取约10ml外周血。将受试者的外周血经过分离血浆、提取cfDNA、全基因组建库和双端测序之后得到cfDNA序列，再将所得cfDNA序列匹配到基因组，获得cfDNA片段来源的基因组位置、长度等片段化模式信息，进而计算全基因组单碱基分辨率的核小体保护分数WPS，最后得到cfDNA中的各个基因的频谱加权核小体周期SWP(即频谱参数)。具体包括如下步骤：

1、分离血浆：从检测受试者体内获得外周血血液，首先4度、1600g离心10分钟使血液分层，然后将上层血浆转移至空管，再4度、16000g离心血浆10分钟以完全除去细胞组分，最后取两步离心后的上清液进行后续实验，或置于-80度冷冻保存。分离血浆不彻底会直接导致提取到的cfDNA中含有血细胞的基因组污染，故一般以两步离心为宜，另外，血浆样本冻存时间过长也会导致cfDNA出现基因组污染现象。

2、提取cfDNA：利用QIAGEN公司的

Circulating Nucleic Acid Kit提取cfDNA。

3、建库和测序：利用ThruPLEX-FD和ThruPLEX DNA-seq试剂盒进行双链建库，并利用测序平台HiSeq 2000和NextSeq 500进行双端测序。

4、利用比对软件bowtie2将测序获得的cfDNA与GRCh37基因组进行比对，并利用samtools软件提取片段端点、长度等信息。

5、计算核小体保护分数(WPS)：对基因组每任一位点i，考虑以i为中心(2d+1)长的窗口[i-d,i+d]内，假设某片段比对到基因组后的端点位置为[lo,hi]，若[lo,hi]与[i-d,i+d]的交集为[i-d,i+d]，则认为片段完全覆盖该窗口；反之，若[lo,hi]与[i-d,i+d]的交集不等于[i-d,i+d]且不为空，则认为片段不完全覆盖该窗口。若完全覆盖该窗口的片段有m条，不完全覆盖该窗口的片段有n条，则WPS＝m-n。依此类推，可计算基因组任一位置的WPS值。在该实施例中，取d＝60以匹配单个核小体的大小，并使用长度为[120,180]的片段计算WPS。

6、计算频谱参数(SWP)：取基因的转录起始位点(Transcription Start Site，TSS)下游10kb区域，即：若基因位于基因组正链，则取[TSS,TSS+10000)；若基因位于基因组负链上，则取(TSS-10000,TSS]。然后，将基因区域内各位置的WPS(记为

)进行快速傅里叶变换提取核小体的频域信息，得到各个周期T下的傅里叶强度

则该基因的频谱参数为：

其中，为了使SWP与基因表达的相关性最大，取a＝168，b＝208。同时，本实施例提供了基于前述SWP特征解耦cfDNA组织来源分布的算法。由于基因的表达与SWP存在负相关关系，因而本发明借助公共数据集中大量的基因表达数据，对检测受试者的cfDNA加权核小体频谱与各组织细胞的基因表达计算相关系数并排序，即可得到该检测受试者cfDNA的来源组织排序。

在计算cfDNA组织来源分布的方法中，所使用的表达数据的组织类型决定了可溯源的cfDNA来源组织类型，进而决定了可检测的癌症类型。本发明中，使用THE HUMANPROTEIN ATLAS数据库中的transcript_rna_tissue.tsv.zip和transcript_rna_celline.tsv.zip作为各组织细胞的基因表达数据集，再计算最大的转录本表达值作为基因的表达值，并将表达值矩阵中的零值替换为0.04取对数作为对数表达矩阵。

进一步地，记cfDNA的频谱参数为

M为基因的个数。记各组织细胞中的对数表达矩阵为A，A中第i行第j列的A_ij代表第i个组织中基因j的表达值，1≤i≤N，1≤j≤M，N为样本个数，将A中第i行记为A_i＝(A_i1,A_i2,…,A_iM)。那么，cfDNA与任一组织i的皮尔逊相关系数为

则

即为组织i在cfDNA来源组织排序中的排名，而cfDNA的来源组织排序为r＝(r₁,r₂,…r_N)。

为了详细阐述频谱参数相比于现有方法的区别，本实施例提供附图1、2、3详细阐述改进后核小体特征与基因表达具有更强的相关性，且该相关性随着测序深度的增大而增强。

附图1展示了高、中、低表达的三个基因的核小体保护分数、改进前后的核小体特征及其与对数表达值的皮尔逊相关系数。其中，现有的核小体特征的计算方式为傅里叶强度在193-199bp(中间图的虚线)的平均值。低表达基因ENSG00000066468的FFT强度主峰在199bp，核小体排布疏松，频谱参数为197.1，现有的核小体特征为1114992，对数表达值为-4.82；中表达基因ENSG00000001617的FFT强度主峰在193bp，相对于低表达基因的核小体排布更紧密，频谱参数为193.6，现有的核小体特征为3767714，对数表达值为1.74；高表达基因ENSG00000003756的FFT强度有多个主峰，核小体排布时而疏松时而密集，频谱参数为185.0，现有的核小体特征为205691，对数表达值为5.04。使用这三个基因改进前后的核小体特征与对数表达值算得皮尔逊相关系数分别为-0.06、-0.91，表明了改进后的核小体特征(频谱参数)与基因表达的相关性更强，更能保持基因表达的高低次序信息。

附图2展示了三个健康人cfDNA测序样本(BH01、IH01、IH02，来源于GEO公共数据集GSE71378)的改进前后核小体特征及其与血细胞系U-937的基因表达间的相关关系。其中，左侧三图的横坐标为现有核小体特征，即193-199bp的平均傅里叶强度，右侧三图的横坐标为改进后核小体特征，即频谱参数SWP，六图的纵坐标均为U-937细胞系的对数表达值，每个点代表一个基因。A、B图为BH01的结果，改进前后的线性回归R分别为0.36、0.46；C、D图为IH01的结果，改进前后的线性回归R分别为0.25、0.42；E、F图为IH02的结果，改进前后的线性回归R分别为0.15、0.37。可见改进后的频谱加权核小体特征与基因表达更相关。

附图3展示了在改进前后测序深度对核小体特征与基因表达的相关性的影响。四图的横坐标均为cfDNA的测序深度，纵坐标均为核小体特征与基因表达的负相关系数的最小值，左侧两图为现有方法，即使用现有核小体特征，右侧两图为改进方法，即使用改进后的核小体特征。具体实施方式为：从公共数据集GSE71378和SRA438908中取65个cfDNA测序样本，计算其测序深度和与基因表达的相关系数绘制为A、B图中的真实数据(空心菱形)。同时，对每一个cfDNA样本，随机置换所有基因的核小体特征100次，再与各组织的基因表达做相关分析求负相关系数的最小值，绘制为A、B图中的置换数据(实心圆点)。因为置换数据打乱了核小体排布与基因表达的相关性，所以其相关系数在两种方法中均不超过-0.05，而真实数据中，改进方法得到的相关系数与测序深度呈现负相关关系。

为了进一步验证这种负相关关系，从65个cfDNA样本中选取11个测序较深的样本(BH01、IC15、IC17、IC20、IC35、IC37、SRR3819936、SRR3819937、SRR3819938、SRR3819939、SRR3819940)进行降采样分析，将原始测序读段按步长为5％，从5％至95％的比例随机采样，每种比例重复50次绘制于图3的C、D图中，同一样本的降采样结果用实线连接。结果证明，改进的核小体特征与基因表达的相关性会随着测序深度的增大而增强，而现有核小体特征与基因表达的相关性大多无此现象，少数样本甚至随着测序深度的增加有减弱的趋势。随着测序成本的不断下降，测序深度的不断加深，改进方法在癌症检测领域将具备更大的潜力。

实施例2

实施例提供一种利用前述cfDNA组织来源分布筛查潜在癌症的检测方法。比较健康人、检测受试者cfDNA的组织来源排序，找出排序异常变化的组织即为检测受试者潜在的癌变组织。从将受试者的组织来源排序与健康人做比较，找出排名变化大的组织即为潜在的癌变组织。

将公共数据集GSE71378中健康人的cfDNA全基因组测序样本(BH01)按照实施例1的方法进行处理，得到健康人cfDNA的来源组织排序r_healthy＝(r₁,r₂,…,r_N)。同样，对于未知cfDNA样本进行同样步骤得到

计算未知cfDNA样本与健康人cfDNA样本的差异

最终，未知cfDNA样本的潜在癌症组织排序为rank(Δr)，将同一种组织的最大排名变化视为该组织的最终排名。本实施例中共使用了40种组织，具体包括肝、胆、脾、肺、肾、膀胱、食管、胃、小肠、结肠、直肠、十二指肠、阑尾、胰岛、唾液腺、扁桃体、甲状腺、甲状旁腺、肾上腺、乳腺、卵巢、输卵管、子宫、子宫内膜、宫颈、***、***、睾丸、附睾、精囊、皮肤、脂肪组织、大脑皮层、眼、心肌、骨骼肌、平滑肌、***、血液、骨髓。

在公共数据集GSE71378和SRA438908中49个癌症病人的cfDNA WGS数据中实施本发明，具体包括膀胱癌、乳腺癌、结直肠癌、食管癌、头颈癌、肾癌、肝癌、肺癌、卵巢癌、胰腺癌、***癌、皮肤癌、睾丸癌、子宫癌、胃癌这15大类癌症，其中，有10例受试者数据测序较深(大于8*基因组)，剩余39例样本测序深度不超过3*基因组。分别使用改进方法(三角)和现有方法(圆点)得到的真实癌症组织在预测组织中的排名见图4，同一样本在两种方法的排名用实线连接。在10例高深度cfDNA样本中，有9例样本使用改进方法的结果优于现有方法，且改进方法将前3名准确率从10％提高到了50％，将前10名准确率从50％提高到了90％。而在39例低深度样本中，由于信息量不足，两种方法总体效果都不足以用于癌症检测。为了进一步探索适合应用本发明的cfDNA测序深度，将9例高深度cfDNA样本(除去用两种方法均预测较差的SRR3819938样本)进行降采样分析，所得结果绘制于图5中，黑色实线为改进方法，灰色虚线为现有方法，黑色虚线为排名为10的基准线。子图标题为样本编号及测序深度，横坐标为降采样比例，纵坐标为预测排名，同样将原始测序读段按步长为5％，从5％至95％的比例随机采样，每种比例重复50次，方差线的最大、最小值分别代表25％、75％分位数。与图3结果类似，随着测序深度的加深，改进方法趋于稳定，且当测序深度高于8*时，除了IC20的预测排名有轻微扰动外，其他样本的真实癌组织始终排名前10。另一方面，现有方法无法随着测序深度的加深而变好，反而在IC17、IC37中会随着测序深度的加深而变差，该现象与图3结论一致。

从以上结果可以看出，本发明使用了更全面的核小体排布信息，首创频谱加权核小体周期特征，检测多种癌种的准确率更高。

在本发明的描述中，术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、“第二”的特征可以明示或者隐含地包括至少一个该特征。在本发明的描述中，“多个”的含义是至少两个，例如两个，三个等，除非另有明确具体的限定。

在本发明中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”、“固定”等术语应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或成一体；可以是机械连接，也可以是电连接或彼此可通讯；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通或两个元件的相互作用关系，除非另有明确的限定。对于本领域的普通技术人员而言，可以根据具体情况理解上述术语在本发明中的具体含义。

在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不必针对的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外，在不相互矛盾的情况下，本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。

尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。

Claims (9)

1.一种确定核酸分子与预定组织相关性的方法，其特征在于，包括：

将所述核酸分子的测序结果与多个参考序列进行比对处理，其中，所述测序结果由多个测序读段构成，所述多个参考序列的每一个分别是基于多个已知基因之一的序列确定的；

基于所述比对处理的结果，针对所述多个参考序列的每一个，分别确定所述核酸分子在所述多个已知基因的频谱参数；

基于所述预定组织中多个所述已知基因的表达量信息，以及所述核酸分子在所述多个已知基因的频谱参数，确定所述核酸分子与所述预定组织的相关性；

其中，所述频谱参数是针对所述参考序列的每一个，分别通过下列步骤确定的：

针对所述参考序列的至少一部分位点，分别确定每个位点的核小体保护分数；

基于所述每个位点的核小体保护分数，确定第一参数和第二参数，所述第一参数表征多个相邻核小体的间距，即多个核小体周期，所述第二参数表征分别与所述多个核小体周期对应的傅里叶强度；

基于所述第一参数和所述第二参数，确定所述参考序列对应所述已知基因的所述频谱参数，

所述第一参数和第二参数是通过下列步骤确定的：

(a)基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数进行平滑矫正后作图得到第一图，其中，x轴为位点与基因转录起始位点的距离，y轴为矫正后的核小体保护分数；

(b)将步骤(a)中得到的所述第一图进行傅里叶变换，得到第二图，所述第二图的x轴为核小体周期，y轴为与核小体周期对应的傅里叶强度；

所述频谱参数是基于下列公式确定的:

其中，T代表核小体周期，F(T)代表核小体周期T对应的傅里叶强度，a代表核小体周期的下限，b代表核小体周期的上限；

所述相关性是通过下列公式确定的：

其中，ρ_i代表所述核酸分子与预定组织i的皮尔逊相关系数，

代表所述核酸分子在多个已知基因的所述频谱参数，A_i代表预定组织i在所述多个已知基因的表达量，Var代表方差，Cov代表协方差。

2.根据权利要求1所述的方法，其特征在于，所述核酸分子为游离于细胞外的核酸分子；

所述参考序列长度为5000～11000bp；

所述参考序列包括所述已知基因的序列；

所述参考序列是基于所述已知基因的转录起始位点确定的；

所述测序是通过高通量测序进行的；

所述已知基因包括Ensembl第75版中GRCh37的所有基因。

3.根据权利要求1所述的方法，其特征在于，针对每个位点，所述核小体保护分数是通过下列步骤确定的：

(1)基于所述位点，确定比对窗口，所述比对窗***有所述位点；

(2)确定与所述比对窗口匹配的所述测序读段；

(3)将所述测序读段区分为第一类型测序读段和第二类型测序读段，所述第一类型测序读段完全覆盖所述比对窗口，所述第二类型测序读段不完全覆盖所述比对窗口；

(4)基于所述第一类型测序读段的数目和所述第二类型测序读段的数目，确定所述核小体保护分数；

所述比对窗口长度为100～140bp；

所述位点位于所述比对窗口的中部；

在步骤(4)中，所述核小体保护分数与所述第一类型测序读段的数目呈正相关，与所述第二类型测序读段的数目呈负相关；

在步骤(4)中，所述核小体保护分数是通过公式第一类型测序读段的数目-所述第二类型测序读段的数目确定的；

基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数，确定所述参考序列对应的多个核小体周期作为第一参数，以及所述多个核小体周期对应的傅里叶强度作为第二参数。

4.根据权利要求1所述的方法，其特征在于，

a＝168，b＝208；

在确定参数a、b后，所述频谱参数仅与所述第一、第二参数有关。

5.根据权利要求1所述的方法，其特征在于，

所述预定组织选自：肝、胆、脾、肺、肾、膀胱、食管、胃、小肠、结肠、直肠、十二指肠、阑尾、胰岛、唾液腺、扁桃体、甲状腺、甲状旁腺、肾上腺、乳腺、卵巢、输卵管、子宫、子宫内膜、宫颈、***、***、睾丸、附睾、精囊、皮肤、脂肪组织、大脑皮层、眼、心肌、骨骼肌、平滑肌、***、血液、骨髓中的至少一种组织。

6.一种确定核酸分子与预定组织相关性的***，其特征在于，包括：

比对装置，所述比对装置用于将所述核酸分子的测序结果与多个参考序列进行比对处理，其中，所述测序结果由多个测序读段构成，所述多个参考序列的每一个分别是基于多个已知基因之一的序列确定的；

频谱参数确定装置，所述频谱参数确定装置与所述比对装置相连，所述频谱参数确定装置针对所述多个参考序列的每一个，分别确定所述核酸分子在所述多个已知基因的频谱参数；

相关性确定装置，所述相关性确定装置与所述频谱参数确定装置相连，所述相关性确定装置基于所述预定组织中多个所述已知基因的表达量信息，以及所述核酸分子在所述多个已知基因的频谱参数，确定所述核酸分子与所述预定组织的相关性；

其中，所述频谱参数确定装置针对所述参考序列的每一个，分别通过如下单元确定所述频谱参数：

核小体保护分数确定单元，所述核小体保护分数确定单元针对所述参考序列的至少一部分位点，分别确定每个位点的核小体保护分数；

参数确定单元，所述参数确定单元与所述核小体保护分数确定单元相连，所述参数确定单元基于所述每个位点的核小体保护分数，确定第一参数和第二参数，所述第一参数表征多个核小体周期，所述第二参数表征分别与所述多个核小体周期对应的傅里叶强度；

频谱参数计算单元，所述频谱参数计算单元与所述参数确定单元相连，所述频谱参数计算单元基于所述第一参数和所述第二参数，确定所述参考序列对应所述已知基因的所述频谱参数，

所述第一参数和第二参数是通过下列步骤确定的：

(a)基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数进行平滑矫正后作图得到第一图，其中，x轴为位点与基因转录起始位点的距离，y轴为矫正后的核小体保护分数；

(b)将步骤步骤(a)中得到的所述第一图进行傅里叶变换，得到第二图，所述第二图的x轴为核小体周期，y轴为与核小体周期对应的傅里叶强度；

所述频谱参数是基于下列公式确定的:

其中，T代表核小体周期，F(T)代表核小体周期T对应的傅里叶强度，a代表核小体周期的下限，b代表核小体周期的上限；

所述相关性是通过下列公式确定的：

其中，ρ_i代表所述核酸分子与预定组织i的皮尔逊相关系数，

代表所述核酸分子在多个已知基因的所述频谱参数，A_i代表预定组织i在所述多个已知基因的表达量，Var代表方差，Cov代表协方差。

7.根据权利要求6所述的***，其特征在于，所述核酸分子为游离于细胞外的核酸分子；

所述参考序列长度为5000～11000bp；

所述参考序列包括所述已知基因的序列；

所述参考序列是基于所述已知基因的转录起始位点确定的；所述测序是通过高通量测序进行的；

a＝168，b＝208。

8.根据权利要求6所述的***，其特征在于，针对每个位点，所述核小体保护分数是通过下列步骤确定的：

(1)基于所述位点，确定比对窗口，所述比对窗***有所述位点；

(2)确定与所述比对窗口匹配的所述测序读段；

(3)将所述测序读段区分为第一类型测序读段和第二类型测序读段，所述第一类型测序读段完全覆盖所述比对窗口，所述第二类型测序读段不完全覆盖所述比对窗口；

(4)基于所述第一类型测序读段的数目和所述第二类型测序读段的数目，确定所述核小体保护分数；

所述比对窗口长度为100～140bp；

所述位点位于所述比对窗口的中部；

在步骤(4)中，所述核小体保护分数与所述第一类型测序读段的数目呈正相关，与所述第二类型测序读段的数目呈负相关；

在步骤(4)中，所述核小体保护分数是通过公式第一类型测序读段的数目-所述第二类型测序读段的数目确定的；

基于所述参考序列的位点信息以及所述位点信息对应的核小体保护分数，确定所述多个核小体周期作为第一参数，以及与所述多个核小体周期对应的傅里叶强度作为第二参数。

9.一种确定血浆游离核酸组织来源的方法，其特征在于，包括：

(i)根据权利要求1～5中任一项所述的方法，确定所述血浆游离核酸与预定组织的相关性；

(ii)将所述相关性进行排序，以便确定所述血浆游离核酸的组织来源；

将步骤(ii)中得到的排名与对照排名进行比较，以便确定排名变化的组织来源；

所述对照排名是基于正常个体或者具有已知状态的个体确定的；

所述排序是按照由高到低进行的，并且选择排序升高的组织作为重点监控组织。

Images