用于衡量基因组不稳定性的检测方法、设备、终端设备和计算机可读存储介质

技术领域

本发明属于生物技术领域，具体的，本发明涉及用于衡量基因组不稳定性的检测方法、 设备、终端设备和计算机可读存储介质。

背景技术

大多数癌细胞都具有基因组不稳定的特征，这是细胞分裂过程中基因组改变的趋势。 癌症通常是由多个控制细胞分裂的基因和抑癌基因的失活导致的。研究显示，基因组的完 整性受到多种监测机制的密切监测，包括DNA损伤检查点、DNA修复机制和有丝分裂检查 点等。这些机制中任何一种调控上的缺陷都常常导致基因组不稳定，从而使细胞易于发生 恶性转化。例如，组蛋白尾部的翻译后修饰与细胞周期和染色质结构的调控密切相关，DNA 甲基化状态与基因组完整性有关。基因组的不稳定性为个体提供了更短的细胞周期和/或绕 过细胞内和免疫控制系统的优势，从而使癌细胞具有生长优势并被选为恶性转化细胞。基 因组不稳定性包括碱基对突变频率增加、微卫星不稳定性(MSI)等小的结构变异，以及染色 体数目或结构改变等显著的结构变异。目前，关于基因组不稳定起源的机制尚未明确。

基因组不稳定导致的染色体结构异常是癌细胞的一个关键特征。癌基因或抑癌基因的 单核苷酸多态性或拷贝数变异，是肿瘤发生与恶化的直接驱动因素。基因组不稳定可引起 基因组范围的遗传变异，这些变异会赋予肿瘤细胞克隆性生长与遗传进化的优势，最终导 致肿瘤的耐药以及肿瘤的复发。因此，在细胞内存在多个进化上保守的通路，它们通过启 动DNA修复过程或启动细胞凋亡来响应这些错误。目前，已知有几种DNA修复途径在DNA 损伤后被激活，一般可分为核苷酸切除修复(NER)、碱基切除修复(BER)、错配修复(MMR)、 DNA双链断裂修复(DSBR)。

在DNA双链损伤中，DNA双链断裂主要通过同源重组(homologous recombinationrepair，HRR)与非同源末端连接(non-homologous end-joining,，NHEJ)两条途径进行 修复。NHEJ不要求断裂末端的序列同源，修复贯穿于整个细胞周期中频繁且复杂性较低的DSB；HRR频率远低于NHEJ，主要介导诸如复制叉断裂等复杂性或危险度较高的DSB修复。DNA损伤修复机制对于维持基因组的稳定性具有十分重要的作用。

近期研究表明，杂合性缺失(LOH)，大区段转换(LST)和跨端粒等位基因不平衡(TAI) 三种基于DNA的基因组不稳定性的评估方法，能够反映肿瘤潜在的同源重组DNA修复缺陷 的程度。

然而，目前对于基因组不稳定性的分析仍有待改进。

发明内容

本发明旨在至少在一定程度上解决相关技术中的技术问题之一。由此，本发明提出了 可用于基因组不稳定性例如拷贝数异常检测的测序数据分析方法、测序数据分析设备、终 端设备和计算机可读存储介质。

根据本发明的第一方面，本发明提供了一种测序数据分析方法或者用于衡量基因组不 稳定性的检测方法，根据本发明的实施例，该方法包括：(1)针对待测样本，获取来自多 个捕获区域的测序数据，所述多个捕获区域的每一个均含有至少一个SNP位点；(2)针对所述多个捕获区域的每一个，确定所述捕获区域的BAF数值和LRR数值，所述BAF数值表 征所述捕获区域中中位SNP位点突变基因型频率，所述LRR数值是通过公式Log2(待测样 本的测序深度/对照样本的测序深度)确定的；(3)基于所述捕获区域的BAF数值和所述 LRR数值，对所述多个捕获区域进行离群点去除，以便获得经过过滤处理的所述多个捕获 区域；(4)基于所述捕获区域的所述BAF数值和所述LRR数值，对至少一条染色体进行至 少一轮分段处理，以便分别获得BAF一级分选片段和LRR一级分选片段；(5)针对所述 BAF一级分选片段和LRR一级分选片段，进行断点集合合并，并基于所述断点的组合，确 定二级分选片段；(6)基于预定的合并阈值，对所述二级分选片段进行迭代式合并处理， 以便获得三级分选片段；(7)基于预定片段长度阈值对所述三级分选片段进行过滤和切割 处理，以便获得多个四级分选片段，所述四级分选片段的长度为3～11Mb，可选的，所述四 级分选片段的长度为4.5～5.5Mb；(8)基于所述四级分选片段中所包含捕获区域的BAF数 值，将所述四级分选片段与对照样本中相应区域的BAF数值进行同分布检验，并基于预定 的检验阈值和各所述四级分选片段中BAF离群点去除前后的密度比值，将各所述四级分选 片段归类为平衡片段、非平衡片段、纯合片段或者非纯合片段，其中，所述四级分选片段 中BAF离群点去除前后的密度比值小于预定的纯合阈值，则所述四级分选片段为纯合片段， 所述四级分选片段中BAF离群点去除前后的密度比值不小于所述预定的纯合阈值，则所述 四级分选片段为非纯合片段；(9)基于所述四级分选片段中所包含捕获区域的BAF数值和 LRR数值，针对所述多个四级分选片段的每一个，分别确定BAF均值、拷贝数均值以及所 述四级分选片段的BAF均值在全部四级分选片段中的排序，所述四级分选片段的拷贝数均 值在全部四级分选片段中的排序；(10)针对所述四级分选片段的BAF均值、拷贝数均值 百分比序数进行平面坐标轴聚类作图，以便确定所述待测样本的肿瘤纯度；(11)基于所 述肿瘤纯度，将所述平衡片段的拷贝数均值及拷贝数百分比序数进行聚类分析，以便确定 所述待测样本的倍性校正因子；(12)基于所述待测样本的肿瘤纯度和所述倍性校正因子， 对所述多个三级分选片段中的每一个分别进行BAF数值校正和拷贝数校正；(13)基于校 正后的BAF数值和拷贝数数值，确定A基因型参数值和B基因型参数值；(14)基于所述 A基因型参数值和B基因型参数值以及所述三级分选片段的长度，确定所述三级分选片段 的变异类型，所述变异类型包括：杂合性缺失、长片段拷贝数状态转变和跨端粒等位基因 不平衡。

根据本发明的实施例，利用上述方法，能够有效地通过对基因组进行测序所得到的测 序数据进行分析，得到相应的基因组不稳定性例如拷贝数异常等，并且能够确定变异类型， 诸如变异类型包括：杂合性缺失、长片段拷贝数状态转变和跨端粒等位基因不平衡。

在本发明的第二方面，本发明还提出了一种测序数据分析设备或者用于衡量基因组不 稳定性的检测设备，根据本发明的实施例，该测序数据分析设备包括：测序数据获取模块， 用于针对待测样本，获取来自多个捕获区域的测序数据，所述多个捕获区域的每一个均含 有至少一个SNP位点；BAF-LRR数值确定模块，用于针对所述多个捕获区域的每一个，确 定所述捕获区域的BAF数值和LRR数值，所述BAF数值表征所述捕获区域中中位SNP位点 突变基因型频率，所述LRR数值是通过公式Log2(待测样本的测序深度/对照样本的测序深度)确定的；离群点去除模块，用于基于所述捕获区域的BAF数值和所述LRR数值，对 所述多个捕获区域进行离群点去除，以便获得经过过滤处理的所述多个捕获区域；分段处 理模块，用于基于所述捕获区域的所述BAF数值和所述LRR数值，对至少一条染色体进行 至少一轮分段处理，以便分别获得BAF一级分选片段和LRR一级分选片段；断点集合合并 模块，用于针对所述BAF一级分选片段和LRR一级分选片段，进行断点集合合并，并基于 所述断点的组合，确定二级分选片段；二级分选片段合并模块，用于基于预定的合并阈值， 对所述二级分选片段进行迭代式合并处理，以便获得三级分选片段；四级分选片段获取模 块，用于基于预定片段长度阈值对所述三级分选片段进行过滤和切割处理，以便获得多个 四级分选片段，所述四级分选片段的长度为3～11Mb，可选的，所述四级分选片段的长度为 4.5～5.5Mb；四级分选片段归类模块，用于基于所述四级分选片段中所包含捕获区域的BAF 数值，将所述四级分选片段与对照样本中相应区域的BAF数值进行同分布检验，并基于预 定的检验阈值和各所述四级分选片段中BAF离群点去除前后的密度比值，将各所述四级分 选片段归类为平衡片段、非平衡片段、纯合片段或者非纯合片段，其中，所述四级分选片 段中BAF离群点去除前后的密度比值小于预定的纯合阈值，则所述四级分选片段为纯合片 段，所述四级分选片段中BAF离群点去除前后的密度比值不小于所述预定的纯合阈值，则 所述四级分选片段为非纯合片段；四级分选片段参数确定模块，用于基于所述四级分选片 段中所包含捕获区域的BAF数值和LRR数值，针对所述多个四级分选片段的每一个，分别 确定BAF均值、拷贝数均值以及所述四级分选片段的BAF均值在全部四级分选片段中的排 序，所述四级分选片段的拷贝数均值在全部四级分选片段中的排序；肿瘤纯度确定模块， 用于针对所述四级分选片段的BAF均值、拷贝数均值百分比序数进行平面坐标轴聚类作图， 以便确定所述待测样本的肿瘤纯度；倍性校正因子确定模块，用于基于所述肿瘤纯度，将 所述平衡片段的拷贝数均值及拷贝数百分比序数进行聚类分析，以便确定所述待测样本的 倍性校正因子；三级分选片段校正模块，用于基于所述待测样本的肿瘤纯度和所述倍性校 正因子，对所述多个三级分选片段中的每一个分别进行BAF数值校正和拷贝数校正；基因 型参数值确定模块，用于基于校正后的BAF数值和拷贝数数值，确定A基因型参数值和B 基因型参数值；变异类型确定模块，用于基于所述A基因型参数值和B基因型参数值以及 所述三级分选片段的长度，确定所述三级分选片段的变异类型，所述变异类型包括：杂合 性缺失、长片段拷贝数状态转变和跨端粒等位基因不平衡。

根据本发明的实施例，利用该测序数据分析设备，能够有效地实施上述测序数据分析方法， 能够有效地通过对基因组进行测序所得到的测序数据进行分析，得到相应的基因组不稳定 性例如拷贝数异常等，并且能够确定变异类型，诸如变异类型包括：杂合性缺失、长片段 拷贝数状态转变和跨端粒等位基因不平衡。

在本发明的第三方面，本发明还提出了一种终端设备，根据本发明的实施例，该终端 设备包括：存储器和处理器；所述存储器用于存储程序代码；所述处理器，调用所述程序代码，当程序代码被执行时，用于执行如前面所述的方法。

在本发明的第四方面，本发明还提出了一种计算机可读存储介质，根据本发明的实施 例，该计算机可存储介质包括指令，当其在计算机上运行时，使得计算机执行如前面所述 的方法。

需要说明的是，本领域技术人员能够理解，在本文中针对用于衡量基因组不稳定性的 检测方法或者用于衡量基因组不稳定性的检测设备所描述的特征和优点，同样适用于其他 方面，在此不再赘述。

本发明的附加方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明 显，或通过本发明的实践了解到。

附图说明

本发明的上述和/或附加的方面和优点从结合下面附图对实施例的描述中将变得明显 和容易理解，其中：

图1是根据本发明一个实施例的测序数据分析方法的流程示意图；

图2是根据本发明另一个实施例的测序数据分析设备的结构示意图；

图3是根据本发明另一个实施例的去除离群点的示意图；

图4是根据本发明另一个实施例的分段方法示意图和倍性校正的参考示意图；

图5是根据本发明另一个实施例的采用公开软件PureCN的分段结果；

图6是根据本发明另一个实施例的S1样本利用PureCN(左图)和本发明的用于衡量基 因组不稳定的检测方法(在本文中有时也称为“GSA算法”)(右图)样本1号染色体分 段对比图；

图7是根据本发明另一个实施例的区域LRR分布对比，依次为：1号染色体0-109M区域(左上)，以及0-3.7M区域(右上)、16-22M区域(左下)，以及74-78M区域(右下)；

图8是根据本发明另一个实施例的区域BAF分布比较，左为1号染色体0-109M区域，右为1号染色体74-78M区域；

图9是根据本发明另一个实施例的S1样本利用PureCN(左图)和本发明GSA算法(右图)样本2号染色体分段对比图；

图10是根据本发明另一个实施例的是根据本发明另一个实施例的S1样本利用PureCN (左图)和本发明GSA算法(右图)样本3号染色体分段对比图；

图11是根据本发明另一个实施例的是根据本发明另一个实施例的S1样本利用PureCN (左图)和本发明GSA算法(右图)样本3号染色体末段对比图；；

图12～16是根据本发明另一个实施例的S1样本利用PureCN(左图)和本发明GSA算法(右图)4-22号染色体末段的分段效果对比；

图17～22是根据本发明另一个实施例的S2样本利用PureCN(左图)和本发明GSA算法 (右图)对不同染色体的分段效果对比；

图23是根据本发明另一个实施例的利用PureCN和本发明GSA算法计算所得的纯度与 理论纯度的对比；

图24是根据本发明另一个实施例的用PureCN和本发明GSA算法计算所得的倍性的对 比。

具体实施方式

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

在本发明的第一方面，本发明提供了一种测序数据分析方法，根据本发明的实施例， 参考图1，该方法包括：

S10：获取测序数据

在该步骤中，针对待测样本，获取来自多个捕获区域的测序数据，所述多个捕获区域 的每一个均含有至少一个SNP位点。

根据本发明的实施例，捕获区域可以通过对已知SNP位点的位置进行分析获得，从而 实现效率和准确性之间的平衡。根据本发明的实施例，可以采用发明人提出的新型确定捕 获区域的方法来对参考基因组划分若干个捕获区域(在本文中有时也称为“捕获区域”)。 从而，针对相应的捕获区域设计捕获探针，以便构建来自这些捕获区域的测序数据。

进一步，本领域技术人员可以通过常规SNP分析方法对针对上述针对特定区域的下机 数据(测序数据)进行信息学分析，可以获得针对待测样本，在已知的SNP位点的测序结 果。例如，可以在对测序数据进行质控，过滤后采用已知的算法进行处理，例如包括但不限于：SOAPsnp、AtlasSNP2、GATK和SAMtools call。在此不再赘述。

根据本发明的实施例，上述捕获区域的测序数据是可以通过采用本发明的试剂盒获得， 该用于衡量基因组不稳定性的试剂盒含有探针集合，所述探针集合是通过下列步骤确定的：

(1)基于预先确定的预期区间数目和预期间隔的至少之一，将参考基因组序列划分为多个一级区域，所述一级区域含有至少一个已知的SNP位点；

(2)针对所述多个一级区域的每一个，进行SNP过滤，以便获得一级高频SNP， 并将所述一级区域划分为含有所述一级高频SNP的一级区域和不含所述一级高频SNP 的一级区域，所述一级高频SNP为所述一级区域中的频率最高SNP并且所述一级高频 SNP的频率不低于预定频率；

(3)选择相邻所述一级高频SNP之间的距离超过预定阈值的区域作为缺口区域，所述预定阈值不低于所述预期间隔的1.5倍；

(4)基于所述预期间隔，将所述缺口区域划分为至少一个二级区域；

(5)针对所述至少一个二级区域的每一个，分别进行二级高频SNP搜寻，所述二 级高频SNP的频率不低于所述预定频率：

(a)由所述二级区域的中心点向两侧至少一次延伸处理，所述延伸处理由所述中心点两侧延伸预定长度；和

(b)在所述至少一次延伸预定长度之后，在所得到的区域中寻找频率最高的SNP，如果该SNP的频率不低于10％，则选择该SNP作为所述二级高频SNP，其中，在获得所 述二级高频SNP时停止所述延伸处理，或者在经过所述延伸处理后的区域长度超过所 述预定阈值时停止所述延伸处理，

(6)基于所述一级高频SNP和所述二级高频SNP汇总作为起始点SNP，确定三级 区域，所述三级区域是在所述起始点SNP的两侧延伸预定长度而确定的；和

(7)基于所述三级区域，构建特异性识别所述三级区域的探针，所述探针适于确定基因组不稳定性。

根据本发明的实施例，采用上述方法得到的探针能够有效地获取相应的捕获区域的数据，并且能够控制捕获区域的数目，提高分析的效率。

S20：确定BAF数值和LRR数值

在该步骤中，针对所述多个捕获区域的每一个，确定所述捕获区域的BAF数值和LRR 数值，所述BAF数值表征所述捕获区域中中位SNP位点突变基因型频率，所述LRR数值是通过公式Log2(待测样本的测序深度/对照样本的测序深度)确定的。

在本文中，BAF是指将每个捕获区域bed区域内的中位SNP位点突变基因型频率作为 该区域的BAF。LRR是通过下列步骤确定的：统计每个捕获区域(有时也称为bed区域)内的平均深度，并进行GC含量和捕获区域校正，根据样本平均测序深度进行归一化，并根据log2(肿瘤组织样本深度/对照样本深度)计算LRR。其中对照样本深度是使用的对照样本(在本文中所采用的术语“对照样本”应做广义理解，既可以是与测序样本相同来源的个 体的其他组织细胞，例如血细胞，也可以来自其他已知病理状态的个体的生物样本，例如 正常人的样本，具体的，可以采用历史血细胞样本的测序数据)经过矫正归一化后的平均 深度。

根据本发明的实施例，针对多个所述捕获区域的每一个，可以通过下列步骤确定所述 BAF数值：

(2-1)确定所述捕获区域中的每一个所述SNP位点的次等位基因频率，以便获得捕获 区域内所有位点的次等位基因频率；

(2-2)确定所述每个捕获区域中所有位点的次等位基因频率的中位数值，作为所述捕 获区域的BAF数值。

根据本发明的实施例，本领域技术人员可以采用常规的方式获得各位点的次等位基因 频率。优选的，针对给定的SNP，可以通过统计支持各SNP类型的测序读段(read)数目来判断次等位基因频率，举例而言，针对给定的SNP位点，共有100条read支持该SNP位 点(与该位点唯一比对)，其中，支持碱基A、T、G和C的测序读段数目分别为50、40、 5和5个，则T为次等位基因，其次等位基因频率为40％。进一步，将捕获区域中所涉及的 全部SNP位点的次等位基因频率进行统计，确定这些次等位基因频率的中位数，从而可以 得到该捕获区域的BAF数值。

根据本发明的实施例，针对每一个所述捕获区域，可以通过下列步骤确定所述LRR数 值：

(2-i)针对所述捕获区域，确定所述待测样本的归一化后的测序深度，以及对照样本 归一化后的测序深度；和

(2-ii)按照下列公式确定所述捕获区域的LRR数值：

Log2(待测样本的测序深度/对照样本的测序深度)。

S30：进行离群点去除

在该步骤中，基于所述捕获区域的BAF数值和所述LRR数值，对所述多个捕获区域进 行离群点去除，以便获得经过过滤处理的所述多个捕获区域。

根据本发明的实施例，可以通过做图，例如按照捕获区域在染色体上的位置变异与相 应的BAF数值或者LRR数值进行做图，从而可以基于所得到的的绘图结果可视化地去除明 显偏离的离群点。通过人工操作既可以肉眼判断明显偏离的数据点。另外，根据本发明的 实施例，还可以通过下列步骤进行离群点的自动识别和处理。

首先，进行所述离群点去除之前，预先使用临近线性插补的方法对异常的BAF数值和 异常的LRR数值进行校正处理，所述校正处理包括：

确定BAF数值或者LRR数值异常(例如，根据本发明的实施例，显示的数值数量级明显与其他数据不相符合的数值，或者被显示为无穷大inf的数据)的捕获区域，在所述捕 获区域的前后各延伸不小于5个捕获区域，并将延伸后所得到区域的平均BAF数值或LRR 数值作为异常数值的校正数值。本领域技术人员能够理解的是，延伸捕获区域的数目足以 使延伸后的校正数值不再为异常数值。需要说明的是，BAF和LRR均以捕获区域为单位， 一个捕获区域一个数值。

参考图3，另外，根据本发明的实施例，所述BAF数值离群点是通过下列步骤确定的： 选择BAF数值小于0.05或者大于0.95的数值点作为BAF数值离群点。根据本发明的实施例，所述LRR数值离群点是通过下列步骤确定的：针对每个捕获区域，分别构建以该捕获 区域为中心的LRR数值-位置索引二元正态分布；针对预定捕获区域，以所述预定捕获区域 为中心，前后延伸预定数目的相邻捕获区域，优选延伸100个所述相邻捕获区域，以获得LRR数值分析区域；基于所述LRR数值分析区域中所包含所述捕获区域的对应LRR数值-位 置索引二元正态分布，确定每个所述LRR数值-位置索引二元正态分布在所述预定捕获区域 位置的概率密度；叠加计算全部在所述预定捕获区域位置的所述概率密度，以便确定所述 预定捕获区域的概率密度总数值；如果所述概率密度总数值低于预定离群点阈值，则确定 所述预定捕获区域为LRR数值离群点。

根据本发明的实施例，首先会给每个点赋予一个二维正态分布的密度函数，该密度函 数协方差均为0，x向方差则参考段落总长度。对于每个点，叠加计算相邻点的概率密度， 对于靠近段落边缘的点，则需要将无对称位置的点的密度做翻倍处理。所有点中密度总值 低于阈值或处于某一分位数(例如30％)以下的，将被去除。

由此，通过上述处理可以得到经过过滤的捕获区域数据。

S40：分段处理，获得一级分选片段

在该步骤中，基于所述捕获区域(本领域技术人员能够理解的是，这里所处理的捕获 区域是指经过步骤S30过滤得到的捕获区域)的BAF数值和所述LRR数值，对至少一条染色体的至少一部分进行至少一轮分段处理，以便分别获得BAF一级分选片段和LRR一级分选片段。

参考图4，

在步骤(4)中，选择完整染色体序列的至少一部分(根据本发明的实施例，优选采用 完整的染色体全长序列，例如Y染色体以外的其他染色体)作为起始分析区域，并基于所述起始分析区域作为根字节，分别基于LRR数值和BAF数值进行所述分段处理，

其中，所述分段处理进一步包括：

(4-1)，并计算待分段处理区域的参数平均值

(4-2)针对所述待分段处理区域中的多个所述捕获区域，确定每个所述捕获区域的参 数与所述参数平均值的差异

(4-3)针对所述待分段处理区域，沿着所述待分段处理区域的延伸方向，由第一个所 述捕获区域开始，针对所述捕获区域的每一个进行所述差异数值的一次累积加和

(4-4)在所述待分段处理区域中，确定所述累积加和的极值点，选择与所述待分段处 理区域端点之一距离均超过预定预定距离阈值的极值点作为分割点，利用所述分割点将所 述起始分析区域进行分割为多个子区域；

(4-5)针对所述多个子区域重复步骤(4-1)～(4-4)至少一次，直到新生成段落不再小于预定长度阈值，以便获得多个含有至少50所述捕获区域的段落，以便获得所述分选片段。

所述预定长度阈值足以包括不小于50个所述捕获区域。

根据本发明的实施例，还可以进一步包括：

(4-6)对长度小于预定长度阈值(例如50个所述捕获区域的长度，诸如3Mb)或SNP数量小于阈值(例如平均50个捕获区域所包含的SNP数目)的临近段落进行合并，并将合 并得到的片段作为所述分选片段，

所述合并进一步包括：

相邻的分选片段分别基于BAF，LRR数值进行同分布的循环统计检验，并将分布偏差小 于同分布检验阈值的相邻段落进行合并。

根据本发明的实施例，在本步骤中，通过递归的方式建立并遍历段落树状结构：首先 将整个段落作为根节点，其所有分割的子段落作为子节点，子节点再次分割后段落作为该 子节点的子节点。

满足下述条件时，需要发起递归判断，方法为：判断该段落是否为根节点，如是则结 束递归调用过程，如否，则回溯该节点的上一个平级子节点递归调用分割过程。如果没有 未处理的平级子节点，则回溯根节点递归调用分割过程。

对于某一待分割段落，分割过程为：首先判断段落是否被分割过，如被分割过，则发 起递归判断。如果段落未被分割过，但所包含数据点数小于等于阈值，则发起递归判断。如果该节点未被分割过，且数据点数大于阈值，标记该段落为已分割，并计算数据偏离均值的游程

如果极大值和极小值点(即断点)分别处于段落两端，或者该点与段落两端之间的差 值小于阈值，则不记录这2个点，并发起递归判断，从两个点中间部分开始进入下一个循环，生成一个子节点。如果其中一个点到段落两端的距离大于阈值，则记录这个点，并将 这整个序列分成两份(生成2个子节点),并记录两个子节点各自的起点和终点。如果两个 点距离段落两端的距离大于阈值，则记录三个分段(生成3个子节点)，并记录三个子节 点各自的起点和终点。按照上述方法直至处理最后一个子节点。

如果分割的子节点长度小于阈值，则向上回溯，处理倒数第二个子节点，否则重新分 段，继续生成新的子节点，并处理最后一个子节点，如果所有子节点处理完，则处理其父节点的平级节点；直至最终回溯到根节点，不再有新的子节点，则结束。

对经过粗分的各个段落进行合并。合并包括两个步骤：首先合并小段落，可以单纯根 据段落长度不满足阈值，或长度不满足阈值且点数不满足阈值进行合并，阈值需参考预期 分辨率进行调整。随后，相邻的片段需进行同分布的统计检验。检验指标为以下三种：

前段与后段点分布的概率密度函数(通常使用核密度估计获取)的差值之和

前段与后段点分布的积累分布函数(对核密度曲线进行数值积分获取)的差值之和 两个段落归并后总的标准差与两个段落各自标准差之加权平均的差值。

通过循环的方式进行合并，首先合并根据某一检测指标判断最接近同分布的两个相邻 片段，并重新计算新合并的段落和其前后段落的指标，如此循环直到所有指标均达到阈值。

S50：进行断点集合合并，确定二级分选片段

在该步骤中，针对所述BAF一级分选片段和LRR一级分选片段，进行断点集合合并，并基于所述断点的组合，确定二级分选片段。

例如，作为示例，BAF一级分选片段的片段位点分别为1～10,20～30，LRR一级分选片 段为5～15,25～35，则断点为1,5,10,15,20,25,30和35，二级分选片段为 1～5,5～10,10～15,15～20,20～25,25～30,30～35。该例子仅做说明使用，而不以任何方式限制 本发明。

S60：对二级分选片段进行合并，获得三级分选片段

在该步骤中，基于预定的合并阈值(例如片段长度不超过1M或片段长度不超过2M且 SNP数量小于50或SNP数量小于5的二级分选片段)，对所述二级分选片段进行合并处理，以便获得三级分选片段。

S70：对三级分选片段进行过滤和切割处理，得到四级分选片段

在该步骤中，基于长度对所述三级分选片段进行过滤和切割处理，以便获得多个四级 分选片段，所述四级分选片段的长度为3～11Mb。

S80：对四级分选片段进行归类

在该步骤中，基于所述四级分选片段中所包含捕获区域的BAF数值，将所述四级分选 片段与对照样本中相应区域的BAF数值进行同分布检验，并基于预定的检验阈值和各所述 四级分选片段含有的BAF离群点去除前后的密度比值，将各所述四级分选片段归类为平衡 片段、非平衡片段、非纯合片段或者纯合片段，其中，所述四级分选片段中所包含的BAF离群点去除前后的密度比值小于预定的纯合阈值，则所述四级分选片段为纯合片段，不小于所述预定的纯合阈值，则所述四级分选片段为非纯合片段。

根据本发明的实施例，段落分割完毕后，选取所有长度超过阈值(如10MB)的段落。将这些段落按照阈值长度进行切割，舍弃掉两端长度不足的段落。然后对于余下所有段落，计算段落信息，包括该段落的BAF均值、BAF均值在所有这些段落中的百分比排序、拷贝 数的均值(2的LRR次方乘以2)，拷贝数在所有段落中的百分比排序。然后将该段落的 BAF分布与对照参比(如前所述，可以为历史血细胞参考基线或者待测对象的血细胞测序 数据)进行同分布检测，如果指标低于平衡性阈值，则将该段落标记为平衡，如果高于不 平衡阈值，则将该段标记为不平衡。如果段落中BAF位点经过离群点去除前后的密度比值 低于纯合性阈值，则将该段落标记为纯合。

根据本发明的实施例，在该步骤中，进一步包括：

(8-1)针对多个所述四级分选片段的每一个，分别确定下列参数的至少之一：

(a)该四级分选片段所包含所述捕获区域的BAF数值平均值以及该BAF数值平均值在 所有所述四级分选片段的捕获区域BAF数值中的排名；和

(b)该四级分选片段所包含所述捕获区域的拷贝数平均值以及所述拷贝数平均值在所 有所述四级分选片段的捕获区域拷贝数中的排名；

(8-2)针对所述四级分选片段的每一个，分别进行所述四级分选片段中的全部BAF数 值在待测样本与对照样本之间的同分布检验，并得到所述四级分选片段的检验指标数值， 如果所述检验指标数值低于平衡性阈值，则将所述四级分选片段归类为平衡片段，如果所 述检验指标数值高于不平衡阈值，则将所述四级分选片段归类为不平衡片段；和

(8-3)针对所述四级分选片段的每一个，如果所述四级分选片段中所包含的BAF离群 点去除前后的密度比值小于预定的纯合阈值，则将所述四级分选片段归类为纯合片段，否 则，将所述二级分选片段归类为非纯合片段。

在本文中多处所使用的术语“阈值”如无特别说明，本领域技术人员可以通过已知相 互关系的样本进行平行试验，进行确定，例如可以通过已知为纯合片段的平行试验，可以 确定纯合阈值。

S90：确定四级分选片段的BAF均值、拷贝数均值及其排序

在该步骤中，基于所述四级分选片段中所包含捕获区域的BAF数值和LRR数值，针对 所述多个四级分选片段的每一个，分别确定BAF均值，拷贝数均值，以及所述四级分选片段的BAF均值在全部四级分选片段中的排序，所述四级分选片段的拷贝数均值在全部四级分选片段中的排序。

在步骤(9)中，进一步包括：

(9-1)在所述四级分选片段中，以所述片段拷贝数平均值百分比序数作为第一轴，以 所述片段BAF平均值作为第二轴，进行平面坐标轴聚类作图，以便确定待测样本中的肿瘤 纯度；和

(9-2)在所述四级分选片段中，选择所有平衡片段，以所述平衡片段的拷贝数平均值 百分比序数作为第一轴，以所述平衡片段的拷贝数平均值作为第二轴，进行平面坐标轴聚 类做图；和

(9-3)基于所述四级分选片段的肿瘤纯度计算结果以及步骤(9-2)得到的平面坐标 轴聚类图，确定待测样本中的倍性校正因子。

S100：确定肿瘤纯度

在该步骤中，针对所述四级分选片段的BAF均值、拷贝数均值百分比序数进行平面坐 标轴聚类作图，以便确定所述待测样本的肿瘤纯度。

S110：确定倍性校正因子

在该步骤中，基于所述肿瘤纯度，将所述平衡片段的拷贝数均值及拷贝数百分比序数 进行聚类分析，以便确定所述待测样本的倍性校正因子。

关于拷贝数和B等位基因含量理论值(CN、B)与拷贝数及BAF实验值(CN*、BAF*)满足如下关系：

CN

(1)

其中scale_factor为倍性校正因子，purity为肿瘤纯度。对于已知拷贝数和B等位基因频率的区域，有：

对于相同基因型的片段，纯度计算值可以近似看作符合以纯度理论值为均值的正态分 布，故有：

通过对BAF-CN百分比序数二维数据进行聚类的方式，找到特定基因型对应的点簇，计 算BAF*均值并带入公式计算纯度。

已知纯度后，可以对所有被标记为“平衡”基因型的区段的拷贝数和拷贝数百分比序 数进行聚类，计算倍性校正因子(在文本中有时称为“倍性校正系数”，这两个术语可以互换使用)。如果标记为平衡的位点个数不足，则梯度增大阈值，比较点数变化。若无明 显变化，输出异常。

常规样本通常存在三种核型：二倍型为主导，二倍型和四倍型均衡，四倍型为主导， 其余情况判定为无法识别。无法识别情况多为纯度较低情形，因此程序默认按照二倍型为 主导处理，并报异常。

根据平衡基因型的密度分布，划分基因型区域。当满足平衡基因型密度在较大范围内 X方向只有一个分布峰，且x取值更小区域基因型密度在y方向仅有一个分布峰时，指认 为二倍型为主导。此时平衡基因型密度大于阈值的区域被指认为二倍型区域，拷贝数更低 的区域指认为单倍型区域，以单倍型区域计算纯度，以二倍型区域计算倍性校正因子。

当平衡基因型密度具有两个或以上峰时，指认为二倍型和四倍型均衡。根据两个分布 峰的位置确定单倍型和三倍型的区域位置，并以B基因型、BB基因型或ABB基因型中点数 符合要求的点簇为基准计算纯度。以AB基因型或AABB基因型点数多的区域计算倍性校正 因子。

当满足平衡基因型密度在较大范围内X方向只有一个分布峰，但x较小区域内y方向 有多个分布峰时，指认为四倍型为主导。程序尝试根据聚类点簇位置关系和密度分布的山 谷位置，划定单倍型至三倍型区域，并以B基因型、BB基因型或ABB基因型中点数符合要求的点簇为基准计算纯度。以AABB基因型区域计算倍性校正因子。

S120：对三级初选片段进行BAF数值校正和拷贝数校正

在该步骤中，基于所述突变样本的纯度和所述倍性校正因子，对所述多个三级初选片 段中的每一个分别进行BAF数值校正和拷贝数校正。

S130：确定A基因型参数值和B基因型参数值

在该步骤中，基于校正后的BAF数值和拷贝数数值，确定A基因型参数值和B基因型参数值。

S140：确定三级分选片段的变异类型

在该步骤中，基于所述A基因型参数值和B基因型参数值以及所述三级分选片段的长 度，确定所述三级分选片段的变异类型，所述变异类型包括：

在该步骤中，杂合性缺失、长片段拷贝数状态转变和跨端粒等位基因不平衡。

根据本发明的实施例，对段落分割程序确定的每个分段，可以计算出纯度校正后的BAF、 拷贝数，并标记其基因型是否平衡、是否为纯合。

如果标记为纯合，则A基因型为0，B基因型为校正后拷贝数取整。

如果标记为平衡，则判断校正后的拷贝数是否大于2，如否，A指定为1，B指定为1，如是，寻找与该值最接近的偶数，A和B指定为该偶数的一半。

如果不属于如上两种情况，则首先计算校正后的拷贝数和校正后的B值，枚举B＝max (B-10，1)～B+min(B+1，8)，A＝max(CN-B-8，0)～min(CN-B+8，B)范围内所有情况， 并计算各情况下BAF与该段落校正后BAF的差值。

对于上一步列表，去除被标记为明显不平衡，但A＝B的选项，去除|A+B-CN|大于a·CN+b 的选项，a和b根据历史数据检测情况推断，默认为0.15和0.5。

剩余选项中，筛选掉计算CN与该段实际CN校正值差值大于阈值、计算BAF与该段实际BAF校正值差值大于阈值，或两差值的乘积大于阈值的选项。

剩余选项中，去计算BAF与实际BAF校正值差值最小的一组选项，作为A、B的值

如上筛选步骤中，如果经过筛选后余下的选项数为0，则选择筛选前BAF差值和CN差 值相乘最小的一组，作为A、B的值。

所有片段A、B值计算完成后，合并相邻A、B值均相同的段落，输出最终检测结果。

根据突变检测结果，可以统计各种衡量基因组不稳定的指标。这里以业界常用的具有 一定临床科学意义的指标，杂合性缺失(LOH)、长片段拷贝数状态转变(LST)跨端粒等位基因不平衡(TAI)为例进行演示。这里杂合性缺失LOH定义为：A＝0且B！＝0，且区间 长度大于15M。根据部分文献报道，当整个染色体只有一个片段且满足上述条件时，不记 作杂合性缺失，发明人提供选项开启或关闭针对这一条件的过滤。

长片段拷贝数状态转变LST定义为A+B！＝2，且长度大于10M的区域。计算前需要首先合并所有长度小于3M的片段。

跨端粒等位基因不平衡TAI定义为统计每个染色体跨端粒且不过着丝粒的一段，如果 A！＝B且长度大于11M阈值，则计入一次。

此外，样本整体基因组倍性值为染色体各个片段区域拷贝数的加权平均。

在本发明的第二方面，本发明还提出了一种测序数据分析设备，其可以用于实施上述 测序分析方法，根据本发明的实施例，该测序数据分析设备包括：测序数据获取模块10、 BAF-LRR数值确定模块20、离群点去除模块30、分段处理模块40、断点集合合并模块50、二级分选片段合并模块60、四级分选片段获取模块70、四级分选片段归类模块80、四级 分选片段参数确定模块90、肿瘤纯度确定模块100、倍性校正因子确定模块110、三级分 选片段校正模块120、基因型参数值确定模块130、变异类型确定模块140。

具体的，根据本发明的实施例，该测序数据分析设备包括：

测序数据获取模块10，该模块用于针对待测样本，获取来自多个捕获区域的测序数据， 所述多个捕获区域的每一个均含有至少一个SNP位点；

BAF-LRR数值确定模块20，该模块用于针对所述多个捕获区域的每一个，确定所述捕 获区域的BAF数值和LRR数值，所述BAF数值表征所述捕获区域中中位SNP位点突变基因型频率，所述LRR数值是通过公式Log2(待测样本的测序深度/对照样本的测序深度)确 定的；

离群点去除模块30，该模块用于基于所述捕获区域的BAF数值和所述LRR数值，对所 述多个捕获区域进行离群点去除，以便获得经过过滤处理的所述多个捕获区域；

分段处理模块40，该模块用于基于所述捕获区域的所述BAF数值和所述LRR数值，对 至少一条染色体进行至少一轮分段处理，以便分别获得BAF一级分选片段和LRR一级分选 片段；

断点集合合并模块50，该模块用于针对所述BAF一级分选片段和LRR一级分选片段， 进行断点集合合并，并基于所述断点的组合，确定二级分选片段；

二级分选片段合并模块60，该模块用于基于预定的合并阈值，对所述二级分选片段进 行迭代式合并处理，以便获得三级分选片段；

四级分选片段获取模块70，该模块用于基于预定片段长度阈值对所述三级分选片段进 行过滤和切割处理，以便获得多个四级分选片段，所述四级分选片段的长度为3～11Mb，可 选的，所述四级分选片段的长度为4.5～5.5Mb；

四级分选片段归类模块80，该模块用于基于所述四级分选片段中所包含捕获区域的 BAF数值，将所述四级分选片段与对照样本中相应区域的BAF数值进行同分布检验，并基于预定的检验阈值和各所述四级分选片段中BAF离群点去除前后的密度比值，将各所述四级分选片段归类为平衡片段、非平衡片段、纯合片段或者非纯合片段，其中，所述四级分 选片段中BAF离群点去除前后的密度比值小于预定的纯合阈值，则所述四级分选片段为纯 合片段，所述四级分选片段中BAF离群点去除前后的密度比值不小于所述预定的纯合阈值，则所述四级分选片段为非纯合片段；

四级分选片段参数确定模块90，该模块用于基于所述四级分选片段中所包含捕获区域 的BAF数值和LRR数值，针对所述多个四级分选片段的每一个，分别确定BAF均值、拷贝数均值以及所述四级分选片段的BAF均值在全部四级分选片段中的排序，所述四级分选片段的拷贝数均值在全部四级分选片段中的排序；

肿瘤纯度确定模块100，该模块用于针对所述四级分选片段的BAF均值、拷贝数均值 百分比序数进行平面坐标轴聚类作图，以便确定所述待测样本的肿瘤纯度；

倍性校正因子确定模块110，该模块用于基于所述肿瘤纯度，将所述平衡片段的拷贝 数均值及拷贝数百分比序数进行聚类分析，以便确定所述待测样本的倍性校正因子；

三级分选片段校正模块120，该模块用于基于所述待测样本的肿瘤纯度和所述倍性校 正因子，对所述多个三级分选片段中的每一个分别进行BAF数值校正和拷贝数校正；

基因型参数值确定模块130，该模块用于基于校正后的BAF数值和拷贝数数值，确定A 基因型参数值和B基因型参数值；

变异类型确定模块140，该模块用于基于所述A基因型参数值和B基因型参数值以及 所述三级分选片段的长度，确定所述三级分选片段的变异类型，所述变异类型包括：

杂合性缺失、长片段拷贝数状态转变和跨端粒等位基因不平衡。

根据本发明的实施例，所述BAF-LRR数值确定模块被设置为针对多个所述捕获区域的 每一个，通过下列步骤确定所述BAF数值：

(2-1)确定所述捕获区域中的每一个所述SNP位点的次等位基因频率，以便获得捕获 区域内所有位点的次等位基因频率；

(2-2)确定所述每个捕获区域中所有位点的次等位基因频率的中位数值，作为所述捕 获区域的BAF数值，

和/或

所述BAF-LRR数值确定模块被设置为针对每一个所述捕获区域，通过下列步骤确定所 述LRR数值：

(2-i)针对所述捕获区域，确定所述待测样本的归一化后的测序深度，以及对照样本 归一化后的测序深度；和

(2-ii)按照下列公式确定所述捕获区域的LRR数值：

Log2(待测样本的测序深度/对照样本的测序深度)。

根据本发明的实施例，所述离群点去除模块被设置为通过下列步骤进行所述过滤处理：

根据本发明的实施例，在进行离群点去除之前，预先对异常的BAF数值和异常的LRR 数值进行校正处理，所述校正处理包括：

确定BAF数值或者LRR数值异常的捕获区域，在所述捕获区域的前后各延伸不小于5 个捕获区域，并将延伸后所得到区域的平均BAF数值或LRR数值作为异常数值的校正数值。

根据本发明的实施例，所述BAF数值离群点是通过下列步骤确定的：

选择BAF数值小于0.05或者大于0.95的数值点作为BAF数值离群点去除，并以BAF＝0.5 为中心进行镜面折叠，即|BAF数值-0.5|+0.5；；

和/或

所述LRR数值离群点是通过下列步骤确定的：

针对每个捕获区域，分别构建以该捕获区域为中心的LRR数值-位置索引二元正态分 布；

针对预定捕获区域，以所述预定捕获区域为中心，前后延伸预定数目的相邻捕获区域， 优选延伸100个所述相邻捕获区域，以获得LRR数值分析区域；

基于所述LRR数值分析区域中所包含所述捕获区域的对应LRR数值-位置索引二元正态 分布，确定每个所述LRR数值-位置索引二元正态分布在所述预定捕获区域位置的概率密 度；

叠加计算全部所述概率密度，以便确定所述预定捕获区域的概率密度总数值；

如果所述概率密度总数值低于预定离群点阈值，则确定所述预定捕获区域为LRR数值 离群点。

根据本发明的实施例，所述分段处理模块被设置为适于：

选择完整染色体序列的作为起始分析区域，并基于所述起始分析区域作为根字节，分 别基于LRR数值和BAF数值进行所述分段处理，

其中，所述分段处理进一步包括：

(4-1)，并计算待分段处理区域的参数平均值

(4-2)针对所述待分段处理区域中的多个所述捕获区域，确定每个所述捕获区域的参 数与所述参数平均值的差异

(4-3)针对所述待分段处理区域，沿着所述待分段处理区域的延伸方向，由第一个所 述捕获区域开始，针对所述捕获区域的每一个进行所述差异数值的一次累积加和

(4-4)在所述待分段处理区域中，确定所述累积加和的极值点，选择与所述待分段处 理区域端点之一距离均超过预定预定距离阈值的极值点作为分割点，利用所述分割点将所 述起始分析区域进行分割为多个子区域；

(4-5)针对所述多个子区域重复步骤(4-1)～(4-4)至少一次，直到新生成段落不再超过预定长度阈值，以便获得多个含有至少50所述捕获区域的段落，以便获得所述分选片段，

其中，所述预定长度阈值足以包括不小于50个所述捕获区域。

根据本发明的实施例，在步骤(4-5)之后，进一步包括：

(4-6)对长度小于预定长度阈值或SNP数量小于阈值的临近段落进行合并，并将合并 得到的片段作为所述分选片段，

所述合并进一步包括：

相邻的分选片段分别基于BAF，LRR数值进行同分布的循环统计检验，并将分布偏差小 于同分布检验阈值的相邻段落进行至少这一轮合并。

根据本发明的实施例，所述四级分选片段归类模块适于：

(8-1)针对多个所述四级分选片段的每一个，分别确定下列参数的至少之一：

(a)该四级分选片段所包含所述捕获区域的BAF数值平均值以及所述BAF数值平均值 在所有所述四级分选片段的捕获区域BAF数值中的排名；和

(b)该四级分选片段所包含所述捕获区域的拷贝数平均值以及所述拷贝数平均值在所 有所述四级分选片段的捕获区域拷贝数中的排名；

(8-2)针对所述四级分选片段的每一个，分别进行所述四级分选片段中的全部BAF数 值行待测样本与对照样本之间的同分布检验，并得到所述四级分选片段的检验指标数值， 如果所述检验指标数值低于平衡性阈值，则将所述四级分选片段归类为平衡片段，如果所 述检验指标数值高于不平衡阈值，则将所述四级分选片段归类为不平衡片段；和

(8-3)针对所述四级分选片段的每一个，如果所述四级分选片段中所包含的BAF离群 点去除前后的密度比值小于预定的纯合阈值，则将所述四级分选片段归类为纯合片段。

根据本发明的实施例，所述四级分选片段参数确定模块适于：

(9-1)在所述四级分选片段中，以所述片段拷贝数平均值百分比序数作为第一轴，以 所述片段BAF平均值作为第二轴，进行平面坐标轴聚类作图，以便确定待测样本中的肿瘤 纯度；和

(9-2)在所述四级分选片段中，选择所有平衡片段，以所述平衡片段的拷贝数平均值 百分比序数作为第一轴，以所述平衡片段的拷贝数平均值作为第二轴，进行平面坐标轴聚 类做图；和

(9-3)基于所述四级分选片段的肿瘤纯度计算结果以及步骤(9-2)得到的平面坐标 轴聚类图，确定待测样本中的倍性校正因子。

利用该测序数据分析设备，能够有效地实施上述测序数据分析方法，从而能够有效地 对基因组进行测序所得到的测序数据进行分析，得到相应的基因组不稳定性例如拷贝数异 常等。

在本发明的第三方面，本发明还提出了一种终端设备，根据本发明的实施例，该终端 设备包括：存储器和处理器；所述存储器用于存储程序代码；所述处理器，调用所述程序代码，当程序代码被执行时，用于执行如前面所述的方法。

在本发明的第四方面，本发明还提出了一种计算机可读存储介质，根据本发明的实施 例，该计算机可存储介质包括指令，当其在计算机上运行时，使得计算机执行如前面所述 的方法。

本发明的附加方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明 显，或通过本发明的实践了解到。

探针设计实施例

预期探针密度评估

考虑到芯片须包含足够数量的区域，维持单位长度基因组有足够的数据点数来保证 变异检测准确性。更高的探针密度意味着更高的分辨率，降低密度必然带来分辨率的损失，但也带来成本方面的优势。因此，本案例基于WGS实测样本进行了探针密度评 估，具体如下：

以2500bp为单位，将WGS划分为约113万个区间。通过均匀间隔抽取算法抽取区 域到目标数量，具体操作为：输入起始点(默认为第一个点)，总区域数(113万)和 预期区域数，那么分段间隔stepSize＝总区域数/预期区域数。

通过设置不同梯度预期区域数，如：60W，50W，40W，30W，20W，10W，5W，然后基 于抽取的不同梯度区域检测到的SNP变异位点评估各染色体区域实际分辨能力，通过 多例WES样本的实测分析，我们认为当预期区域低于10万时，染色体基因型分辨能力 会有较大幅度降低。因此，我们得到一个预期探针密度为10W。

候选SNP位点筛选

SNP位点选择主要依据来自人群频率数据库，本案例中参考的是采用的千人数据库 和dbSNP数据库。筛选条件主要为人群频率，默认使用东亚人群频率以提高位点选取 成功率。筛选阈值可以设定多档，默认为高频：≥10％，中高频：≥5％。且只筛选常染 色体及X染色体上的SNP位点。

基于人群频率可以初步获得一个候选SNP文件，主要包含：染色体，染色体坐标，人群频率，频率等级这些信息。

基于预期探针密度，对候选SNP位点进行二次筛选

基于第一步中得到的基因组上均匀分布的区域(10W)，可以计算出预期间隔。然后过滤筛选区域内的高频SNP位点(第二步中得到的数据库中候选高频人群突变位点)， 对于没有找到高频位点的区域，则扩大搜索范围(左右增加若干bp，默认1000)再次 进行搜索。随后，删除没有高频SNP位点的区间，并对存在高频SNP位点的区间，以 人群频率最高的位点为中心，左右各延长一定区域(默认左延60，右延59)。

填补gap区域

基于二次筛选的SNP位点，计算临近的SNP位点的距离，若二者距离大于预期间隔的1.5倍，则定义为gap区域，从gap区域的中心位置逐步扩展搜索范围(默认每次 扩展500BP)搜索高频SNP位点，如果扩展到达预期间隔阈值，则停止搜索，保持该区 域为空。对所有非空区域，对于该区域选定的SNP位点，左右各延长一定区域(默认 左延60，右延59)。

提交探针设计区域

合并第3,4步骤中的所有SNP位点，左右各延长一定区域(默认左延60，右延59) 作为最终候选探针区域，并提交探针供应商进行探针合成。

实施例1—本发明方法(GSA方法)中拷贝数检测及纯度校正模块

为了体现发明人针对大片段基因型检测所做的优化，发明人选择行业评价较优但未做 针对性优化的公开软件进行对比。PureCN为ABSOLUTE算法的一个比较完善的实现，该算 法是目前根据拷贝数数据进行纯度校正的主流方法。PureCN软件更新较为活跃，目前发明 人使用稳定的1.8.1版本进行测评。除了拷贝数计算本身，本发明的方法与PureCN纯度计 算的优劣亦在本对比中一个指标。

这里发明人将展示每个捕获区域所计算得到的BAF和LRR的散点图。软件给出的段落 划分依次使用红绿蓝三种背景色体现，颜色仅用于区分临近段落，本身并无特殊含义。本 实施例采用两对配对的样本，以及样本与配对按比例混合的情形。发明人将从分段、拷贝 数检测和纯度校正结果等若干方面比较。

分段倾向：

对比染色体分段前，发明人首先作如下规定：PureCN由于没有专门为染色体级别的长 片段CNV检测进行优化，因而并不区分短CNV和长CNV，对于局部少数点的数据波动，PureCN 会进行划分，导致结果分段细碎，不能很好地体现染色体区域的整体趋势。以图5为例， 左图示意了该染色体所有分段，右图则隐去了长度小于3M的段落。在实际对比中，发明人 只对比PureCN长度大于3M的段落。结果显示在图5中，左图为全部分段，右图为去除3M以下段落后的分段，背景的红绿蓝色并没有实际意义，仅用来便于区分和计数。

参考图6～16和图17～22接下来，发明人首先以样本S1为例，比较两个方法拷贝数检 测的差异。如无特殊说明，每个染色体均展示左右两个分段图，左为PureCN结果，右为GSA结果。首先发明人观察一号染色体，可以看到0-50M区域拷贝数方面存在一些数据波动，可能是真实的拷贝数变化，或存在亚克隆，也可能来自NGS数据的测试偏差。两个软 件都能一定程度上去除这些离群点的影响，但在具体表现上还是存在一些差异。首先， 0-3.7M区域、16-22M区域，以及74-78M区域的拷贝数分布上与周边并不相同(如图3所 示)，容易观察到，74-78M区域的BAF分布也不同于周边(如图4所示)，PureCN的分段 程序并没有区分出这些片段，而GSA方法可以进行区分。从后续表格中可以看到，PureCN 对于0-109M整个区域给出了3拷贝的结论，而观察其BAF，集中分布在0.5附近，说明其 AB基因比是接近1：1的，显然不可能是3拷贝，这应该是其对整个片段不加划分，并错 误地进行平均的结果。而GSA给出的结果看，0-120M大部分区域均为AABB型，4拷贝，0-3.7M 区域给出了10拷贝的结果，而16-22M区域给出了6拷贝的结论。74-78M区域的BAF不支 持AB基因型比为1：1，所以给出了5拷贝的结论。可以认为GSA方法给出的结论是更合 理的。

表1PureCN和GSA方法1号染色体各段落位置和拷贝数情况对比

类似的，很容易从散点图上看到，164-178M区域的拷贝数不同于周边，周边区域AB比为1：1，而此处的拷贝数要更低，GSA计算出150-164M和178-201M两区域的拷贝数为 4，164-178M区域的拷贝数为3，而PureCN给出的从143-228M的拷贝数都为3，这是明显 不合理的。综合来看，GSA可以正确地划分染色体区域，并给出更加合理的拷贝数计算。

接下来分析2号染色体，从散点图可以看到，PureCN最大的分段问题在于76M位置没 有划分，本发明的算法可以做出正确的区分，另外PureCN在21M、90-92M、和147M处进 行了断开，从数据分段看前后并不存在显出差距，发明人认为不应该断开更加合理。需要 指出的是，GSA方法有对染色体着丝粒gap区域做特别处理，如果gap前后区域分布相同 是会做合并处理的。但是，虽然PureCN对于染色体进行了分段，但四个片段都给出了拷贝 数3的结论，发明人则第一个分段判定为3拷贝，后续较大的片段均指认为4拷贝。基于 BAF分布接近0.5，以及拷贝数整体分布更偏大的理由，可以认为GSA方法的拷贝数计算更 加合理。

表2PureCN和GSA方法2号染色体各段落位置和拷贝数情况对比

接下来分析3号染色体的情况，PureCN将0-100MB当作一整个段落，但实际上可以看 到在染色体的一些区域，如22-25MB和31-35MB，64-67MB等区域，拷贝数的分布是明显不同于其他区域的，这并不是局部波动，而是整段具有一致性趋势的行为。GSA算法可以正 确识别这些段落，并给出更高拷贝数(ABBB基因型，4拷贝)的结果，而PureCN将整段计 算为3拷贝。另外PureCN没有正确处理着丝粒附近的GAP，GSA方法则做了正确处理。当 然GSA在后半也存在一些误判，比如没有正确区分175M及以后区域中的断点。下面展示了 BAF散点图在该区域的放大，实际上从数据上，这一区段落之间的过渡趋势较为缓和，并 非那么容易找到明确的断点。不过实际上仍然可以调整分段阈值，以获得更精细的划分。 目前所采用的默认参数实际上是对于常规数据做出一定权衡的结果，如果数据异常，实际 上可以根据实际情况做一些调整。

表3PureCN和GSA方法3号染色体各段落位置和拷贝数情况对比

其余段落大体上符合类似的趋势。总体上，GSA算法比PureCN给出的结果更加可靠。这里 仅列出各染色体对比结果，不再做细致讨论。

表4PureCN和GSA方法4-22号染色体各段落位置和拷贝数情况对比

样本S2分段性质表现类似，各染色体分段情况如图17～22，左为prueCN(去掉了长度 小于3M的分段)，右为GSA。

GSA程序在分段方面的表现较PureCN的改善是一目了然的。总体来说，GSA分段更能 体现染色体级别的整体趋势，数据波动(可能由于局部拷贝数不同，也可能由于测序质量 原因造成的假信号)会被离群点去除程序消除。PureCN则分段较为零散，如1号染色体107-243MB的长片段(跨越着丝粒)，虽然有一些数据波动，但整体上观察，会发现其表 现出统一的基因型特征，但在PureCN的分段下会被断成多段。同时一些断点在PureCN结 果中没有被正确识别，如1号染色体约76M位置，3号染色体的多个位置，明显存在基因 型的转变，GSA程序正确识别了这个断点，但PureCN没有办法识别。类似的，PureCN对2 号染色体整体性分段有诸多断开，但是48-54M区域一个基因型明显不同于前后的段落，却 没有正确识别。而GSA软件的分段可靠性明显好于PureCN。对于染色体gap处理方面，可 以观察到4号染色体着丝粒gap前后基因型相同，但PureCN并没有对其进行合并。5号染 色体gap前后基因型不同，PureCN则没有进行正确的切分，GSA软件则可以正确处理这种 情况。其他染色体表现类似，不再复述。

表5PureCN和GSA方法S2样本各个染色体各段落位置和拷贝数情况对比

2.2.2纯度计算：

这里通过将高纯度样本中混入对照样本的方式实现纯度的梯度。对于稀释的样本，正 确的校正可以获得一致的结果。但是需要注意的是，对于NGS方法，混合并非简单的比例 相除，而是需要考虑dup率的变化。

假设样本中总reads数为N，那么可以将其分为两部分，有效reads数Ne1和重复reads 数Nd1。前者是N中不重复的reads的种类数，后者则为N与前者的差值，表示至少有一个重复的reads数。总的dup率r1＝Nd/N。Nd1中实际上也可以获取到一个唯一的集合Ne2，相应的存在一个对应的dup集合Nd2，二阶的dup率r2＝Nd2/Nd1。Nd2中又可以存在独立 的集合Ne3，对应重复集Nd3和三阶dup率r3。以此类推。由此可以构成一个有穷级数(Nei,Ndi,ri)。理论上，有一些reads比其他reads更容易产生dup，因而越高阶，ri将越大， 并且可以遇见的是，在一阶dup率不高的情况下，Nei占N的比例会快速下降，有效的阶 数不会太多。

当发明人从N中抽取Nc的reads时，需要考虑其组成成分。Nc一部分来自Ne1，这部分是没有重复的，但另一部分来自Nd1，这部分与Ne1整个集合存在重复，但是其中一部 分可能是未被选择的一部分(Ne1-Ne1c)的重复，这部分当中的非重复集合Ne2c实际上不 与Ne1c重复。类似的，Ne3c也不会与Ne1c和Ne2c重复，实际上的dup率rc应等于 (Nc-Ne1c-Ne2c-Ne3c-...)/Nc。但是实际上不难发现，越高阶的Neic占比实际上是越低的， 不做特别严格探讨时，高阶项实际上可以忽略。

为了简化计算，这里发明人忽略三阶以上的项，当从总reads数N中抽取Nc时，有效reads数可以约略使用如下公式计算：

实际纯度可以如下计算：

据此发明人可以根据样本本身的病理纯度和计算的到的混合比，计算混样后的理论纯 度。继而发明人可以将软件计算的纯度与理论值进行对比。PureCN算法尝试计算不同纯度 和倍性下，理论数据与实测数据的拟合优度，高于阈值的结果都会被保留。其问题在于， 不一定是拟合最好的结果就是最合理的，而算法本身并不提供判断，需要人工根据各种情 况综合判断。GSA算法尝试通过聚类找到最优的纯度倍性，直接可以找到最优结果。

表6PureCN和GSA算法对于混样样本纯度和倍性的计算

从表6的结果可以看到，两个样本和各自三个混样样本共八个中，只有一半的情况， PureCN计算得到的第一选择是比较接近理论值的，而非第一选择的情况，有可能回是比较 靠后的选择较为接近，这样实际操作中，是比较难以执行的。而GSA算法总的来说给出的 选择都是比较接近理论频率的。当然实际上计算的到的纯度与病理纯度可能并不完全相同， 更加准确的表述应该是主克隆对应的细胞含量。但是从图23结果中也可以看到，发明人的 方法计算的到的纯度与稀释得到的理论纯度维持良好的线性关系，这说明其计算能够正确 地反应稀释带来的变化，而PureCN的结果并不能具有这样的自洽性。

对于倍性的形况，如图24所示，GSA算法给出的各个样本的倍性的波动是较小的—— 因为倍性是样本的本质属性，不应该随混入正常样本数量变化而变化，因而波动越小越说 明校正算法能够正确发挥作用。须知，纯度越低，样本的校正越困难。可以看到对于S1样 本，后两个混合样本PureCN的倍性计算出现了很大偏差，而GSA第三个虽然也产生了偏差， 但仍然比PureCN表现要好。S1总的来说是高倍性的样本，平均倍性接近3倍体，而PureCN 在第三个混样的计算结果已经接近2倍体，说明此时校正是失败的。而对于S2，由于样本 平均倍性就是接近2的，校正难度要小得多，两个算法都能得到波动较小的结果。但总的 来说，GSA给出的波动仍然更小。综合来说，发明人的GSA算法在纯度计算和校正方面明显由于PureCN。

具体实施案例2——基因组不稳定性及同源同组基因修复基因突变检测综合评估的试剂盒

(1)样本选取

随机选取3例商业化细胞系样本，利用本发明所设计的方法和试剂盒，进行实施例验 证。

(2)样本片段化

本试剂盒兼容100～400ng DNA样本，使用不含核酸酶超纯水(NF水)将其体积补齐至75μl。将液体转移至打断管内，放入核酸打断仪(Bioruptor pico)中。打断仪预冷 至4℃，打断参数为：

(3)末端修复

1)提前取出末修缓冲液室温解冻，涡旋并短暂离心。酶试剂在使用之前取出；

2)按下表加入各组分，吹吸混匀并置于冰盒上；

3)PCR仪上设置末端修复程序如下，将末端修复反应体系置于PCR仪器中反应；

4)末端修复程序结束后，1.5倍磁珠纯化(磁力板纯化，可人工纯化或自动化纯化)。

加入150μl磁珠A，最后用28.1μl DNA溶解液回溶于0.2ml管中。

(4)末端加’A’

1)取出加’A’缓冲液室温解冻、涡旋并短暂离心。酶试剂用之前取出；

2)按下表加入各组分，吹吸混匀并置于冰盒上；

3)在PCR仪上设置加‘A’程序如下，将加”A”反应体系置于PCR仪器中；

4)加’A’程序结束后，不纯化，直接进行下一步反应；

(5)接头连接

1)取出接头连接缓冲液、接头室温解冻，涡旋并短暂离心。酶试剂用之前取出；

2)按下表加入各组分，吹吸混匀，置于冰盒上；

样本与接头编号对应关系如下：

接头序列如下：

3)PCR仪上设置接头连接程序如下，将接头连接反应体系置于PCR仪器中；

4)接头连接程序结束后，磁珠A纯化(板纯，可人工纯化或自动化纯化)。加入30μl无核酸酶水和50μl磁珠A纯化，最后用41μl DNA溶解液回溶并转移至新的0.2ml PCR管中。

(6)文库扩增

1)取出标签引物、PCR反应液4℃解冻，涡旋并短暂离心后置于冰盒中；

2)按下表配制反应体系，置于冰盒暂存；

3)在PCR仪上设置程序如下，将上一步反应体系置于PCR仪器中，开始反应；

样本、接头引物对应关系如下表：

引物序列如下：

4)PCR程序结束后，磁珠A纯化(板纯，可人工纯化或自动化纯化)。加入100μl磁 珠A纯化，最后用32μl无核酸酶水回溶于0.2ml PCR管中保存。

(7)文库混样(Pooling)

根据Qubit结果进行文库混样，按下表所示：

(8)杂交文库准备

1)提前取出所需用量的磁珠B，室温平衡至少30min。

2)取20μL杂交试剂1到新1.5mL离心管中，加入上一步准备好的杂交底物。

3)加入130μL磁珠B，涡旋10s充分混匀。

4)室温静置10min，瞬时离心后，置于磁力架上静置1分钟至液体澄清。

5)小心吸取并弃去上清液，加入190μL 80％乙醇，静置30s至液体澄清。

6)小心吸取并弃去上清液。保持离心管在磁力架上，打开管盖，室温静置5min风干。

7)每管分别加入2μL杂交试剂2、11.4μL无核酸酶水，从磁力架上取下后，充 分涡旋混匀。

(9)芯片杂交

1)将杂交试剂3、杂交试剂4提前取出，常温解冻；将HRD探针从-20℃冰箱中拿出，放在冰上解冻。按下表配置杂交mix；

注：以上组分可以预先混合到一起保存，并命名为“杂交液”。

2)分别取43μL准备好的杂交反应液加入到上一步13.4μL的杂交底物中(带磁珠)。充分混匀，室温静置2min后，置于磁力架上至液体澄清。

3)取出杂交探针，冰上融化。取出对应杂交反应数的新PCR管，标记为本次杂交文库 号，分别加入4μL探针。

4)从步骤(2)转移56.4μL上清液至步骤(3)含有4μL HRD探针的PCR管中。涡 旋10s充分混匀。注意核对标记。

5)瞬时离心后置于PCR上，运行下表的杂交反应程序。

注意：PCR仪热盖需设置为105℃。

(10)洗脱前准备

1)杂交洗脱液准备

提前把恒温金属浴仪调至55℃，每个杂交反应准备100μl 1*Wash Buffer I和400μl 1*Stringent Wash Buffer进行预热。

2)准备杂交磁珠

a)按下表8准备1X磁珠C洗液：

表8磁珠C洗液准备

b)取出磁珠C，室温平衡30min。使用前涡旋15s，充分混匀。

c)取新的1.5mL离心管，每管加入50μL磁珠C。置于磁力架上1min至液体澄清， 小心吸取并弃去上清液。

d)每管加入100μL的1X磁珠C洗液。从磁力架上取下，充分混匀。把离心管放回磁力架上静置至少1min至液体澄清，小心吸取并弃去上清液。重复此步骤一次，共 洗2次。

e)每管加入50μL的1X磁珠C洗液。从磁力架上取下，涡旋10s充分混匀。转移到 新的1.5mL离心管中。

f)将离心管置于磁力架上静置至少1min，至液体澄清，小心吸取并弃去上清液。

g)注意：应立即进入下一步。避免时间过久，磁珠过于干透。

3)杂交文库与磁珠C结合

经杂交后，继续保持杂交混合液于55℃，待磁珠C清洗完成后，迅速将步骤二的 杂交反应混合液转移至上一步准备好的磁珠C中，盖好管盖，涡旋10s充分混匀，置 于PCR仪上55℃孵育15min(热盖105℃)。

(11)洗脱过程

a)经过15min孵育后，从PCR仪中取出杂交混合液、1X洗脱液1。

b)向每个含有60.4μL杂交混合液-磁珠C的离心管中，分别加入100μL预热到55℃的1X洗脱液1。涡旋10s充分混匀。置于磁力架上至液体澄清，小心吸取并弃去 上清液。

c)每管分别加入200μL预热到55℃的1X洗脱液S。从磁力架上取下后，涡旋10s 充分混匀，并放回PCR仪中55℃孵育5min(热盖105℃)。取出后置于磁力架上 至液体澄清，小心吸取并弃去上清液。重复此步骤一次。

d)每管分别加入200μL的1X洗脱液1(室温)，涡旋10s充分混匀。室温静置1min， 瞬时离心后置于磁力架上至液体澄清。小心吸取并弃去上清液。

e)每管分别加入200μL的1X洗脱液2，涡旋10s充分混匀。瞬时离心后转移到新 的1.5mL离心管中。室温静置1min后，置于磁力架上至液体澄清。小心吸取并弃 去上清液。

f)每管分别加入200μL的1X洗脱液3，涡旋10s充分混匀。室温静置1min，瞬时 离心后置于磁力架上至液体澄清。小心吸取并弃去上清液。

g)每管分别加入20μL的PCR-grade water/无核酸酶水，用移液器充分吹打混匀。瞬时离心后，将液体全部转移至新PCR管中(含磁珠)。

(12)捕获文库PCR反应及纯化

1)按照下列的配比准备捕获后PCR反应混合物并进行反应：

2)PCR产物磁珠纯化

准备工作：将磁珠B提前室温平衡30min，vortex 15s充分震荡混匀备用；

用1.0x磁珠(50ul)纯化，最终溶解在32μl ddH2O中。

引物序列如下表：

(13)测序及数据分析

将上步构建好的文库经电泳检测合格后进行华大MGISEQ-2000测序仪测序。取构建的 单链环状DNA文库进行DNA纳米球制备、MGISEQ-2000上机测序。测序过程严格按照MGISEQ-2000的标准操作流程进行上机操作。

下机数据经过过滤、比对、去重、质控后，完成基因组snv变异分析，并通过对杂合性(LOH)，端粒等位基因失衡(TAI)和大规模状态转变(LST)进行测算，获得同源重组 缺陷评分(HRDScore)以评估基因组不稳定性。

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

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