一种基于cnv结果判定样本降解的方法
阅读说明：本技术 一种基于cnv结果判定样本降解的方法 (Method for judging sample degradation based on CNV result ) 是由 贺洪鑫 梁萌萌 余伟师 栗海波 李珉 于 2021-09-24 设计创作，主要内容包括：本发明公开了一种基于CNV结果判定样本降解的方法,包括以下步骤：生成测试样本的目标覆盖CNN文件；构建一个参考CNN对照集；从测试样本的目标覆盖CNN文件中检测出拷贝数变异,生成包含拷贝数变异的具体染色体区间和变异类型信息的中间文件,获取候选参数指标；以历史CNV-seq样本的目标覆盖CNN文件作为实验数据,构建分类模型并进行评估,将能完全区分降解样本和正常样本的分类特征及相应阈值作为判断样本是否降解的最终参数指标；利用测试样本的CNV结果对照该最终参数指标,判断测试样本是否降解。该方法给出了能够准确甄别出降解样本的参数指标包括对应的阈值范围,能够自动、高效、准确地区分降解样本和正常样本。(The invention discloses a CNV result-based method for judging sample degradation, which comprises the following steps: generating a target coverage CNN file of the test sample; constructing a reference CNN control set; detecting copy number variation from a target coverage CNN file of a test sample, generating an intermediate file containing specific chromosome intervals of the copy number variation and variation type information, and acquiring candidate parameter indexes; taking a target coverage CNN file of a historical CNV-seq sample as experimental data, constructing a classification model and evaluating, and taking classification characteristics and corresponding thresholds which can completely distinguish a degraded sample from a normal sample as final parameter indexes for judging whether the sample is degraded; and comparing the final parameter index with the CNV result of the test sample to judge whether the test sample is degraded. The method provides a parameter index which can accurately discriminate the degraded sample and comprises a corresponding threshold range, and can automatically, efficiently and accurately discriminate the degraded sample from the normal sample.)
技术领域
本发明涉及生物学与精准医学高通量测序与变异检测技术领域，具体涉及一种基于CNV结果判定样本降解的方法。
背景技术
近年来，随着高通量测序技术(Next-Generation Sequencing，NGS)的不断发展，全基因组测序（Whole Genome Sequencing，WGS）、全外显子组测序（Whole ExomeSequencing，WES）以及拷贝数变异测序（Copy Number Variation Sequencing, CNVseq）等检测技术已被越来越多的人所熟知。
这些测序技术要成功实施，第一步便是样本DNA提取，可以说，DNA提取是打好基础的关键步骤。但是，由于受温度、湿度、PH值、氧化反应以及微生物侵染等外界因素的干扰，DNA在提取的过程中不可避免地会出现降解效应。众所周知，当测序样本出现降解时不仅预示着该样本可能存在污染，而且还会使得测序结果产生误差，并最终造成遗传解读时的错误结论。为了避免这些情况的发生，则需事先对出现降解的样本进行判定。
目前判定降解样本的常用方法是凝胶电泳法，又称跑胶法，它的判定依据主要是通过查看是否出现电泳图谱脱尾的现象，如果有脱尾，则预示着该样本为降解样本。采用该种方法往往会带来以下三点弊端：
a.采用跑胶法需设计实验和操作实验，特别是操作实验过程，耗时耗力；
b.在操作实验的过程中难免会出现人为失误，从而使得判定的结果不准确；
c.操作实验需要有相关专业技术背景的人才能实施，对判断人员要求较高。
另外，由于上述判定方法是在测序之前完成，而在判定之后到测序过程这段时间样本是否发生降解或者污染则无法判断。
发明内容
本发明的目的是针对现有技术中存在的上述问题，提供一种基于CNV结果判定样本降解的方法，直接根据质控指标参数的阈值即可自动甄别出发生降解的DNA样本，既简单又高效，省时省力，普适性高，不需操作人员具备专业知识。
本发明技术方案详述如下：
一种基于CNV结果判定样本降解的方法，包括以下步骤：
（1）根据参考基因组，生成测试样本的目标覆盖CNN文件；
（2）以历史CNV-seq样本的目标覆盖CNN文件作为对照样本，根据参考基因组将多个对照样本合并，构建一个参考CNN对照集；所述历史CNV-seq样本包括正常样本和降解样本；
（3）根据参考CNN对照集，从测试样本的目标覆盖CNN文件中检测出拷贝数变异，生成包含拷贝数变异的具体染色体区间和变异类型信息的中间文件；
（4）从中间文件中获取用于以CNV结果判断测试样本是否降解的候选参数指标；
（5）以步骤（2）的历史CNV-seq样本的目标覆盖CNN文件作为实验数据，将实验数据随机分成训练集和测试集，并将步骤（4）的候选参数指标进行随机组合，一种组合形式作为一个分类特征，构成多个分类特征；以训练集构建分类模型，用测试集验证分类模型的性能，从而对各个分类特征进行分类性能评估，将能够完全区分降解样本和正常样本的分类特征及相应阈值作为判断样本是否降解的最终参数指标；
（6）根据测试样本的CNV结果，采用步骤（5）的最终参数指标判断该测试样本是否降解。
CNV，拷贝数变异，是基因结构变异的重要组成部分，有基因组发生重排而导致，一般指长度为1kb以上的基因组大片段的拷贝数增加或者减少，主要表现为亚限位水平的缺失和重复。基于高通量测序技术的染色体异常检测，可以获得拷贝数变异结果，CNV结果会显示CNV对应的染色体区间信息、染色体的类型等信息，以及Filter_Count、stdev、segments、mad等参数。通过将测试样本的CNV结果与上述方法得出的适用于该测试样本的分类特征和相应阈值进行比较，即可获知该测试样本是否降解。
可选或优选的，上述基于CNV结果判定样本降解的方法，步骤（1）生成测试样本的目标覆盖CNN文件的具体方法如下：
a.根据参考基因组版本号获取对应的目标bed文件；
b.根据目标bed文件中给定的bed区域计算出比对后bam文件的覆盖度，获得目标覆盖CNN文件。
可选或优选的，上述基于CNV结果判定样本降解的方法，步骤（2）构建参考CNN对照集的具体方法如下：
获取参考基因组序列文件，为参考基因组序列文件创建索引，以加快对照样本的测序序列比对到参考基因组上的速度；
以历史样本CNV-seq的目标覆盖CNN文件作为对照样本，计算对照样本的测序深度和参考基因组序列每个区域的GC含量，将所有的对照样本合并成一个正常基因组的测序分布模型，即参考CNN对照集。
可选或优选的，上述基于CNV结果判定样本降解的方法，步骤（3）生成中间文件的具体方法如下：
a.以参考CNN对照集作为标准，将测试样本的目标覆盖CNN文件的区域覆盖度和GC含量的偏差进行矫正，获得拷贝数比率表文件；
b.从拷贝数比率表文件中推断出离散的拷贝数片段；
c.从离散的拷贝数片段中得到每个片段的绝对拷贝数；
d.以常染色体和女性性染色体的cn值不能等于2，男性性染色体的cn值不能等于1作为过滤条件，对每个片段的绝对拷贝数进行过滤，形成过滤后文件；
e.对过滤后文件进行拷贝数的变异类型判定，当cn值＜2，则变异类型为缺失，当cn值＞2，则变异类型为重复，最终生成包含拷贝数变异的具体染色体区间和变异类型信息的中间文件。
可选或优选的，上述基于CNV结果判定样本降解的方法，步骤（4）中所述候选参数指标包括stdev、segments、mad、检出CNV的总数目和检出CNV具体类型对应的数量中的至少两个。
可选或优选的，上述基于CNV结果判定样本降解的方法，步骤（5）的最终参数指标为stdev和检出CNV的总数目。
可选或优选的，上述基于CNV结果判定样本降解的方法，当检出CNV的总数目>50，和/或Stdev>0.5时，判断测试样本为降解样本。
与现有技术相比，本发明具有如下有益效果：
该方法直接根据最终参数指标和相应的阈值即可自动甄别出已获得CNV结果的测试样本是否降解，既简单又高效，省时省力，且所得结果并不提前于CNV检测，而是与CNV检测为相同的一份样本，相较于先跑胶检测样本是否降解，再将样本进行高通量测序获得CNV结果的方式，本发明的方法准确性更高。
该方法普适性较高，易上手，无任何相关专业知识背景的人员也可进行核查判定，对人员技术要求较低。
该方法流程部署简单，使用操作方便，只需部署相关计算节点即可完成全流程分析。对服务器计算资源要求较低，一台普通的8核心64G内存的服务器，能够允许同时运行几十个目标基因的处理任务。
附图说明
图1为实施例1基于CNV结果判定样本降解的方法整体流程示意图；
图2为实施例1步骤1的流程示意图；
图3为实施例1步骤2的流程示意图；
图4为实施例1步骤3的流程示意图；
图5为实施例1步骤4的流程示意图；
图6为实施例1步骤5的流程示意图；
图7a为实施例1步骤5中方案一Filter_Count+Stdev组合的SVM分类边界图；
图7b为实施例1步骤5中方案一Filter_Count+Mad组合的SVM分类边界图；
图7c为实施例1步骤5中方案一Stdev+Mad”组合的SVM分类边界图；
图8a为实施例1步骤5中方案二Filter_Count+Stdev组合的SVM分类边界图；
图8b为实施例1步骤5中方案二Filter_Count+Mad组合的SVM分类边界图；
图8c为实施例1步骤5中方案二Stdev+Mad”组合的SVM分类边界图；
图9a为实施例1步骤5中方案三Filter_Count+Stdev组合的SVM分类边界图；
图9b为实施例1步骤5中方案三Filter_Count+Mad组合的SVM分类边界图；
图9c为实施例1步骤5中方案三Stdev+Mad”组合的SVM分类边界图；
图10为实施例1验证步骤中正常样本1 CNV结果的染色体图；
图11为实施例1验证步骤中异常样本1 CNV结果的染色体图；
图12为实施例1验证步骤中异常样本2 CNV结果的染色体图。
具体实施方式
下面结合附图和较佳的具体实施例对本发明的技术方案进行详细解释和说明，以使本领域技术人员能够更好地理解本发明并予以实施。
实施例1
请参考图1，整体流程概述如下：
1、目标覆盖CNN文件的生成；
2、参考CNN对照集的构建；
3、拷贝数变异的检出；
4、用于自动质控的参数指标的获取；
5、用于区分降解和正常样本的特征和阈值的获取。
下面对每个部分的操作步骤进行详述。
1、目标覆盖CNN文件的生成
根据参考基因组，生成测试样本的目标覆盖CNN文件，主要是用于记录按照目标bed文件中给定的区域计算出的比对后bam文件的覆盖度，并联合后续构建的参考CNN对照集（下文以Reference.cnn命名）进行CNV的检出。
请参考图2，构建过程如下：
参考基因组可以从公共基因数据库如Ensembl,、NCBI等下载，目前参考基因组版本有两种，分别是hg38和hg19。
a. 首先获取参考基因组版本号，根据参考基因组版本号获取对应的目标bed文件；
b. 接着利用CNV分析软件CNVkit里面的coverage方法，该方法可以根据给定的bed区域计算bam文件的覆盖度，按照目标bed文件中给定的区域来计算出比对后bam文件的覆盖度，并最终获得目标覆盖cnn文件。
输入文件：
测试样本经比对后的bam文件、特定参考基因组版本的目标bed文件。
相关软件：
Cnvkit软件中的coverage方法。
输出文件：
测试样本目标覆盖CNN文件。
2、参考CNN对照集的构建
以历史CNV-seq样本（包括正常样本和降解样本）的目标覆盖CNN文件作为对照样本，根据参考基因组将多个对照样本合并，构建一个参考CNN对照集。该步骤用于以一定数量的历史CNV-seq样本的目标覆盖CNN文件为基础，构建参考CNN对照集，与测试样本的目标覆盖CNN文件联合进行CNV的检出。
请参考图3，构建过程如下：
a. 首先获取参考基因组版本号（hg19或hg38），根据参考基因组版本号从公共基因数据库（Ensembl,NCBI等中）下载对应的参考基因组序列文件，参考基因组序列文件为FASTA格式数据，下文以ref.fa简称表示参考基因组序列文件；
b. 利用序列比对软件中的索引构建模块对下载到的参考基因组序列文件进行参考基因组比对索引的构建，生成索引文件，主要包含ref.fa.amb、ref.fa.bwt、ref.fa.ann、ref.fa.fai、ref.fa.misa、ref.fa.pac、ref.fa.sa；
c.接着获取一定数量的历史cnvseq样本的目标覆盖CNN文件，并以此作为对照样本，利用cnvkit工具中的reference方法，通过计算对照样本的测序深度和参考基因组序列每个区域的GC含量将所有的对照样本合并生成一个正常基因组的测序分布模型，该测序分布模型即为参考CNN对照集，记为Reference.cnn文件。
输入文件：
一定数量的历史cnvseq样本的目标覆盖cnn文件、参考基因组序列文件；
相关软件：
序列比对索引构建软件，
cnvkit工具中的reference方法；
输出文件：
参考CNN对照集。
3、拷贝数变异的检出
根据参考CNN对照集，从测试样本的目标覆盖CNN文件中检测出拷贝数变异，生成包含拷贝数变异的具体染色体区间和变异类型信息的中间文件。该步骤的作用是以Reference.cnn文件作为参考CNN对照集，从测试样本的目标覆盖cnn文件中检测出拷贝数变异，并最终形成cnr、cns、call.cns、filter.cns、CNV.bed等中间文件，为下一步筛选候选参数指标做准备。
请参考图4，构建过程如下：
a.根据给定的对照集，利用cnvkit工具中的fix方法，将测试样本的目标覆盖CNN文件的区域覆盖度和GC含量的偏差进行校正，并输出拷贝数比率表文件（*.cnr）；
b.接着利用cnvkit工具中的segment方法，从步骤a输出的拷贝数比率表文件（*.cnr）中推断出离散的拷贝数片段（*.cns文件）；
c.利用cnvkit工具中的call方法从离散的拷贝数片段（*.cns文件）中得到每个片段的绝对拷贝数（*.call.cns文件）；
d.接着对*.call.cns文件进行过滤，过滤条件为常染色体和女性性染色体的cn值不能等于2，男性性染色体的cn值不能等于1，形成过滤后的filter.cns文件；
e.最后对过滤filter.cns文件中拷贝数的变异类型进行判定。变异类型有两种：重复gain和缺失loss，判定条件为：如果cn<2 则该变异类型为loss，即缺失，如果cn>2则为gain，即重复。最终生成包含拷贝数变异具体染色体区间和变异类型信息的CNV.bed文件，即中间文件。
输入文件：
测试样本目标覆盖cnn文件、Reference.cnn对照集；
相关软件：
cnvkit工具中的fix方法，
cnvkit工具中的segment方法，
cnvkit工具中的call方法；
输出文件：
cnr、cns、call.cns、filter.cns、CNV.bed文件。
4、用于自动质控的参数指标的获取
该步骤从中间文件中获取用于以CNV结果判断测试样本是否降解的候选参数指标。
请参考图5，构建过程如下：
a.利用cnvkit工具中的metrics方法将cnr和cns文件联合生成包含stdev、segments、mad等参数在内的质控指标；
b.对filter.cns文件进行计数，得到检出CNV的总数目；
c.对CNV.bed文件进行筛选计数，得到检出CNV具体类型（gain和loss）对应的数量；
d.将stdev、segments、mad等参数、检出CNV的总数目、检出CNV具体类型（gain和loss）对应的数量作为以CNV结果判断测试样本是否降解的候选参数指标。
输入文件：
cnr、cns、filter.cns、CNV.bed文件等中间文件；
相关软件：
cnvkit工具中的metrics方法；
输出文件：
包含各个能自动质控降解样本的候选指标的结果文件。
5、用于区分降解和正常样本的特征和阈值的获取
该步骤以历史CNV-seq样本的目标覆盖CNN文件作为实验数据，将实验数据随机分成训练集和测试集，并将候选参数指标进行随机组合，一种组合形式作为一个分类特征，构成多个分类特征；以训练集构建分类模型，用测试集验证分类模型的性能，从而对各个分类特征进行分类性能评估，将能够完全区分降解样本和正常样本的分类特征及相应阈值作为判断样本是否降解的最终参数指标。
请参考图6，详细过程如下：
（1）数据获取：
样本来源：历史cnvseq样本；
样本总数：510个（正常样本：489个，异常样本：21个）。
（2）将510个样本按照之前描述的详细步骤将bam文件分别生成cnr、cns、filter.cns和CNV.bed文件等中间文件，利用cnvkit工具中的metrics方法将cnr和cns文件联合生成包含Stdev、Segments、Mad这3个评价参数在内的文件，接着再将filter.cns文件进行计数获取每个样本检出的CNV总数（后文以Filter_Count命名），最终将这4个参数作为候选参数，参与后续的测试比较。每个样本都有对应这4个参数的数值。
（3）对4个候选参数进行随机组合，每一个组合作为一种分类特征，则共有15种情况，分别为: “Filter_Count”、“Segments”、“Stdev”、“Mad”、“Filter_Count+Segments”、“Filter_Count+Stdev”、“Filter_Count+Mad”、“Segments+Stdev”、“Segments+Mad”、“Stdev+Mad”、“Filter_Count+Segments+Mad”、“Filter_Count+Segments+Stdev”、“Filter_Count+Stdev+Mad”、“Segments+Stdev+Mad”以及“Segments+Stdev+Mad+Filter_Count”。
（4）经典的二分类模型SVM是一种线性分类器，属于机器学习方法中的一种，主要通过寻找一个最优的决策边界将两个类别进行区分。
利用上述模型SVM对这15个分类特征的区分性能进行评估。首先将510个样本的4个参数值按照15种的组合类型构建15个不同的矩阵，矩阵中的内容分别为各个样本在不同组合模式下对应的参数值以及样本的状态（以数字表示，正常样本用“0”表示，异常样本用“1”表示）。
接着将这15个矩阵依次进行SVM算法验证，验证方法为：将矩阵数据按照8比2的比例随机分成训练集和测试集，将训练集用于训练和构建二分类模型，将测试集对构建的二分类模型进行测试和验证，并最终引入Score、Intercept、Coefficients这3个参数对各个分类特征的分类性能进行评估，其中Score值越大表示该特征的准确率越高，表示其分类能力越好，15个特征对应的具体的数值结果如下表所示。表格中Score数值由大到小排列，选择排名前3的特征。
特征
Score
Intercept
Coefficients
Stdev+Mad
0.99
-1.69347366
1.85678699,1.99126666
Filter_Count+Stdev
0.97
-2.13857193
0.05151421,-0.27803567
Filter_Count+Mad
0.97
-1.95546197
0.04871313,-0.27713227
（5）接着设计三种不同的实验测试情况，继续对这3个分类特征的区分性能进行评估比较，并进行SVM分类边界可视化来直观地反映这3个特征的分类能力，具体实验方案如下：
方案一：将510个历史cnvseq样本随机分成训练集和测试集，其中80%为训练集，20%为测试集。
方案二：人工手动选取数据构建训练集和测试集，将408个（正常样本：392个，异常样本：16个）cnvseq样本作为训练集，将102个（正常样本：97个，异常样本：5个）cnvseq样本作为测试集。
方案三：基于方案二中的样本数据，首先将方案二中测试集中的5个异常样本拆开，将其中的2个放入训练集中，将剩余的3个继续保留在测试集中，此时的训练集中共有410个样本（正常样本：392个，异常样本：18个），测试集则共有100个样本（正常97个，异常3个）。
（6）依照上述设计的三种实验方案，依次利用SVM二分类模型对“Stdev+Mad”、“Filter_Count+Stdev”、“Filter_Count+Mad”这3个分类特征的区分性能进行评估，并画出SVM分类边界图。三种不同的实验方案下3种不同特征各自的SVM分类边界图如图7-9所示。
图7a~图7c为方案一的SVM分类边界图，从图中可以发现，“Filter_Count+Stdev”特征能将降解和正常CNV-seq样本明显区分开，“Filter_Count+Mad”特征效果次之，“Stdev+Mad”特征则无法将两者完全区分开，效果最差。
图8a~图8c为方案二的SVM分类边界图，从图中可以发现，“Filter_Count+Stdev”和“Filter_Count+Mad”特征能将出现降解和正常cnvseq样本明显区分开，“Stdev+Mad”特征则无法将两者完全区分开，效果最差。
图9a~图9c为方案三的SVM分类边界图，从图中可以发现，“Filter_Count+Stdev”特征依然能将降解和正常cnvseq样本明显区分开，“Filter_Count+Mad”特征区分效果次之，“Stdev+Mad”特征则无法将两者完全区分开，效果依然最差。
综上所述，三种实验方案下，“Filter_Count+Stdev”特征均显示出了较好的对于降解和正常cnvseq样本的区分能力，“Filter_Count+Mad”特征次之，“Stdev+Mad”特征效果不明显。因此，最终选择“Filter_Count+Stdev”特征组合作为自动质控CNV数过多的参数指标。
（7）最后，根据“Filter_Count+Stdev”特征的三个不同SVM分类边界图上临界样本的Filter_Count和Stdev的真实数值来最终确认出能够区分降解和正常cnvseq样本的中间阈值：Filter_Count>50或Stdev>0.5。即当Filter_Count>50或Stdev>0.5时，则表示该cnvseq样本为降解样本。
因此，将“Filter_Count+Stdev”特征，以及Filter_Count>50或Stdev>0.5，作为判断样本是否降解的最终参数指标。
根据测试样本的CNV结果，对应上述最终参数指标进行判断，即可获知该测试样本是正常样本还是已经发生了降解的样本。
为了进一步验证这一标准和阈值的准确性，挑选已知CNV-seq样本（从上述方法中用到的历史CNV-seq样本中选择的3个，其中1个正常样本，2个异常样本）。
验证方法为——通过查看这些已知CNV-seq样本实际检出的CNV状况是否与应用这一标准和阈值判断的结果相一致。
挑选的已知CNV-seq样本详细情况说明：
正常样本1：Filter_Count<50且Stdev<0.5；
异常样本1：Filter_Count>50；
异常样本2：Stdev>0.5。
3个样本实际检出CNV情况与对应的Filter_Count和Stdev参数值如下所示：
正常样本1（Filter_Count<50 且 Stdev<0.5）：
该样本Filter_count=36，Stdev=0.27，则可判定该CNV-seq样本为正常样本
实际检出CNV状况如图10所示，其中箭头所指的蓝色区域表示该CNV类型为重复，圈出来的红色区域则表示该CNV为缺失，由图可见该CNV-seq样本为正常样本，与特征和阈值判定的结果一致。
异常样本1（Filter_Count>50）：
该样本Filter_count=566，Stdev=0.497，则可判定该CNV-seq样本为降解样本
实际检出CNV状况如图11所示，其中箭头所指的蓝色区域表示该CNV类型为重复，圈出来的红色区域则表示该CNV为缺失，由图可见该CNV-seq样本异常CNV数过多，为降解样本，与特征和阈值判定的结果一致。
异常样本2（Stdev >0.5）：
该样本Filter_count=63，Stdev=1.59，则可判定该CNV-seq样本为降解样本
实际检出CNV状况如图12所示，其中箭头所指蓝色区域表示该CNV类型为重复，圈出来的红色区域则表示该CNV为缺失，由图可见该CNV-seq样本异常CNV数过多，为降解样本，与特征和阈值判定的结果一致。
综上所述，应用这一标准和阈值对这些已知CNV-seq样本判断的结果与这些已知样本实际检出的CNV状况相一致。
本文中应用了具体个例对发明构思进行了详细阐述，以上实施例的说明只是用于帮助理解本发明的核心思想。应当指出，对于本技术领域的普通技术人员来说，在不脱离该发明构思的前提下，所做的任何显而易见的修改、等同替换或其他改进，均应包含在本发明的保护范围之内。
- 上一篇：一种医用注射器针头装配设备
- 下一篇：一种基于图神经网络的乳酸菌抗菌肽预测方法