模拟目标疾病仿真测序文库的方法及其应用
技术领域
本发明涉及生物技术领域,具体地,涉及模拟目标疾病仿真测序文库的方法及其应用。
背景技术
在基因测序领域,从一代、二代到三代甚至四代,测序技术在短时间内有了飞跃的发展。随着高通量测序的普及以及相应技术的快速更新迭代,诞生了许多全新的检测方法和流程。它们大多还不够完善,性能和适用范围等也缺少准确而清晰的定位,这些都需要通过测试去确定。一般来说,测试(“benchmarking”)会检测真实的测序下机数据,在检出变异之后还需要后续相应的实验验证,这导致该方法费时费力且费钱。数据模拟的方法可以很好地避开这些缺点,既不需要真实的下机数据,也不需要后续繁琐的实验验证步骤,迅速、高效且效费比高。目前可用于变异模拟的软件有很多,但多针对全基因组,无法模拟捕获测序下机数据,并且能模拟的变异类型有限,同时,模拟的算法服从,处理速度慢。
由此,同时能用于模拟全基因组和捕获测序下机数据的方法有待进一步研究。
发明内容
本发明旨在至少解决现有技术中存在的技术问题之一。为此,本发明的一个目的在于提出一种模拟目标疾病仿真测序文库的方法,既能模拟全基因组下机数据,也能模拟捕获测序下机数据,并能模拟多种变异类型以及目标疾病样本的特征,同时对测序数据进行拟合,使其仿真程度更高。
需要说明的是,本发明是基于发明人的下列工作而完成的:
基于现有的变异模拟方法,为了实现对捕获测序下机数据的模拟,首先,发明人将变异模拟限制在了捕获测序的目标区域,并且根据捕获测序的特点对变异的种类及输出等做了特别处理;其次,针对变异后的测序序列进行单独校检;此外,发明人进一步对下机数据进行深度波动拟合,更加真实地模拟出了捕获测序条件下的下机数据。
因而,根据本发明的一个方面,本发明提供了一种模拟目标疾病仿真测序文库的方法。根据本发明的实施例,该方法包括:获取碱基序列;基于胚系变异数据库,在所述碱基序列上添加胚系变异模拟生成胚系变异序列;基于目标疾病的体系变异数据库,在所述胚系变异序列的基础上迭代添加体系变异生成至少一条含有胚系和体系变异的序列;从所述胚系变异序列和所述含有胚系和体系变异的序列中选取部分变异作为纯合变异,分别以所述胚系变异序列和所述含有胚系和体系变异的序列为模板,得到胚系变异二倍体和含有胚系和体系变异的二倍体;以及将所述胚系变异二倍体和所述含有胚系和体系变异的二倍体转换为下机序列,以便得到胚系下机序列集和含有胚系和体系变异的下机序列集;分别对所述胚系下机序列集和所述含有胚系和体系变异的下机序列集中的所述下机序列进行深度拟合,以便得到拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集;以及按预定比例从所述拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集选项部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列,所述部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列构成所述目标疾病仿真测序文库。
根据本发明的实施例的模拟目标疾病仿真测序文库的方法,能根据需要得到不同体系和胚系变异特征、杂合/纯合比例和不同患病序列纯度,并能准确地模拟出接近真实的捕获测序条件下目标区域的深度波动,从而更加真实地模拟出了捕获测序条件下的下机数据。并且,该方法既能模拟全基因组下机数据,也能模拟捕获测序下机数据,适用范围广。同时,该方法运行速度快,能够在较短的时间内生成所需的模拟序列,并且模拟得到的序列的仿真程度高。
在此基础上,本发明提供了一种模拟目标疾病仿真测序文库的系统。根据本发明的实施例,该系统包括:碱基序列获取装置,用于获取碱基序列;胚系变异序列生成装置,所述胚系变异序列生成装置与所述碱基序列获取装置相连,基于胚系变异数据库,用于在所述碱基序列上添加胚系变异模拟生成胚系变异序列;体系变异序列生成装置,所述体系变异序列生成装置与所述胚系变异序列生成装置相连,基于目标疾病的体系变异数据库,用于在所述胚系变异序列的基础上迭代添加体系变异生成至少一条含有胚系和体系变异的序列;二倍体生成装置,所述二倍体生成装置分别与所述胚系变异序列生成装置和所述体系变异序列生成装置相连,用于从所述胚系变异序列和所述含有胚系和体系变异的序列中选取部分变异作为纯合变异,分别以所述胚系变异序列和所述含有胚系和体系变异的序列为模板,得到胚系变异二倍体和含有胚系和体系变异的二倍体;以及下机序列转换装置,所述下机序列转换装置与所述二倍体生成装置相连,用于将所述胚系变异二倍体和所述含有胚系和体系变异的二倍体转换为下机序列,以便得到胚系下机序列子集和含有胚系和体系变异的下机序列子集;深度拟合装置,所述深度拟合装置与所述下机序列转换装置相连,用于分别对所述胚系下机序列集和所述含有胚系和体系变异的下机序列集中的所述下机序列进行深度拟合,以便得到拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集;以及纯度模拟装置,所述纯度模拟装置与所述深度拟合装置相连,用于按预定比例从所述拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集选取部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列,所述部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列构成所述目标疾病仿真测序文库。
根据本发明的实施例的模拟目标疾病仿真测序文库的系统,能根据需要得到不同体系和胚系变异特征、杂合/纯合比例和不同患病序列纯度,并能准确地模拟出接近真实的捕获测序条件下目标区域的深度波动,从而更加真实地模拟出了捕获测序条件下的下机数据。并且,该系统既能模拟全基因组下机数据,也能模拟捕获测序下机数据,适用范围广。同时,该系统运行速度快,能够在较短的时间内生成所需的模拟序列,并且模拟得到的序列的仿真程度高。
本发明的附加方面和优点将在下面的描述中部分给出,部分将从下面的描述中变得明显,或通过本发明的实践了解到。
附图说明
本发明的上述和/或附加的方面和优点从结合下面附图对实施例的描述中将变得明显和容易理解,其中:
图1显示了根据本发明一个实施例的模拟目标疾病仿真测序文库的方法的流程示意图;
图2显示了根据本发明一个实施例的模拟目标疾病仿真测序文库的系统的结构示意图;
图3显示了根据本发明一个实施例的胚系变异序列生成装置的结构示意图;
图4显示了根据本发明一个实施例的体系变异序列生成装置的结构示意图;
图5显示了根据本发明一个实施例的深度拟合装置的结构示意图;
图6显示了根据本发明一个实施例的变异校验装置的结构示意图。
具体实施方式
下面详细描述本发明的实施例,所述实施例的示例在附图中示出,其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的,仅用于解释本发明,而不能理解为对本发明的限制。
在本发明的描述中,术语“纵向”、“横向”、“上”、“下”、“前”、“后”、“左”、“右”、“竖直”、“水平”、“顶”、“底”等指示的方位或位置关系为基于附图所示的方位或位置关系,仅是为了便于描述本发明而不是要求本发明必须以特定的方位构造和操作,因此不能理解为对本发明的限制。
需要说明的是,术语“第一”、“第二”仅用于描述目的,而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此,限定有“第一”、“第二”的特征可以明示或者隐含地包括一个或者更多个该特征。进一步地,在本发明的描述中,除非另有说明,“多个”的含义是两个或两个以上。
根据本发明的一个方面,本发明提供了一种模拟目标疾病仿真测序文库的方法。根据本发明的实施例的模拟目标疾病仿真测序文库的方法,能根据需要得到不同体系和胚系变异特征、杂合/纯合比例和不同患病序列纯度,并能准确地模拟出接近真实的捕获测序条件下目标区域的深度波动,从而更加真实地模拟出了捕获测序条件下的下机数据。并且,该方法既能模拟全基因组下机数据,也能模拟捕获测序下机数据,适用范围广。同时,该方法运行速度快,能够在较短的时间内生成所需的模拟序列,并且模拟得到的序列的仿真程度高。
根据本发明实施例的模拟目标疾病仿真测序文库的方法,能够在较短的时间内生成所需的含有变异的序列,比如在915个总长约200kb的目标区域内随机模拟200个单核苷酸位点变异(SNV)和300个插入缺失标记(indel),并同时生成对应的变异序列用时<30s。此外,迭代的层次约多,每次迭代的平均用时会更少,比如一次序列生成需要约21s,而4次迭代也不会超过22s。
根据本发明实施例的模拟目标疾病仿真测序文库的方法,通过迭代实现多层次变异多背景的叠加,并允许在任意层次添加指定突变。
根据本发明实施例的模拟目标疾病仿真测序文库的方法,能够获得含有一些特殊的变异组合的序列用于测试。
参考图1,根据本发明的实施例,对该模拟目标疾病仿真测序文库的方法进行解释说明,具体如下:
S100获取碱基序列
根据本发明的实施例,获取碱基序列。
其中,本文中的术语“碱基序列”是指使用一串字母表示的真实的或者假设的携带基因信息的DNA或者RNA分子的一级结构。
根据本发明的实施例,该碱基序列可以是捕获测序得到的。捕获测序不同于全基因组测序,它只检测感兴趣的部分区域,因此成本低、检测速度快、测序深度高。同时,由于探针的亲和性差异,捕获测序的reads数据也具有自己独特的深度波动特征。为了能够模拟出接近真实的捕获测序的下机数据,它一般需要具有以下几个特点:(1)能够指定变异,参入部分真实存在的正常的或肿瘤的变异后模拟的数据会更真实;(2)能够限制模拟的区域,实现捕获的初步特征,如果目标区域一个变异都没有那也就失去了实际意义;(3)能够实现包括snv、indel、CNV、SV在内的所有变异;(4)能够进一步地模拟出reads和肿瘤样本的一些特征,包括GC倾向性、体系/胚系差异、染色体倍数、异质性和肿瘤纯度等。相对于全基因的变异模拟,在捕获测序条件下,有些变异所导致的结果是重叠的。比如目标区域外到目标区域内的复制和移位最终表现出的都只是一段已知序列的插入。
进一步地,针对捕获测序的序列特征,即相对于全基因,捕获测序指是针对特定区域的一段序列,部分变异可能跨越目标区域,同时位于目标区域内和目标区域外。进而,为了判定变异是否跨越目标区域,也就是越界,发明人在标记后的特征串上设置了越界判定部,也就是说,标记后的特征串包括目标片段部和越界判定部,其中,越界判定部位于该目标片段部的两端,用于判定目标片段的两端是否存在越界的变异。
S200生成胚系变异序列
根据本发明的实施例,基于胚系变异数据库,在所述碱基序列上添加胚系变异模拟生成胚系变异序列。其中,“胚系变异”也就是胚细胞变异。
根据本发明的实施例,所述模拟生成胚系变异序列的方法包括:所述碱基序列进行胚系变异状态标记,以便得到胚系变异标记后的特征串;选取待添加的胚系变异;将所述待添加的胚系变异整合至所述胚系变异标记后的特征串上,以便得到添加胚系变异后的特征串和胚系变异信息;以及将所述添加胚系变异后的特征串进行碱基还原,以便得到胚系变异序列。由此,通过对碱基序列的变异状态进行标记,设定碱基的变异类型,从而对各种变异进行模拟,变异模拟的方法简单,生成速度快,并能根据需要设计特殊的变异组合用于测试,变异模拟后的序列的仿真程度高。
根据本发明的实施例,该变异状态标记为数值标记,也就是利用不同的数值标记不同的变异类型,不仅标记简单,而且,可以通过数值的加和进行多层次变异的叠加。
由于二进制只有0和1两个数字,与有和无的判定相符,根据本发明的实施例,该变异类型的设置以二进制数值表示有和无,得到二进制表示的该设置后的标记位;该数值转换包括将该二进制表示的设置后的标记位转换成十进制数字,该十进制数字为该碱基的变异类型标记。由此,通过二进制进行变异类型有和无的判断,便于通过数值进制的转换进行变异类型的标记,易于对不同变异类型的叠加的表示,变异状态标记易于存储和解读。
为了便于理解,举例说明变异状态标记的方法,对每一个长度为N的目标片段,生成一个相应的长度为N+2,初始值为2的变异状态标记位数组。两侧的额外两个标记位为越界判定部,用于记录变异的越界状态。标志位以十进制整数形式存储,但以二进制格式解读。如图2所示,采用了二进制的前6位来标记各种状态。例如,根据指定的变异信息(90个snp和100个inDel)对相应的标志位进行操作。比如3号染色体的第1000个碱基由G突变成了T,根据记录差值的公式得到并记录它们的差值2。同时将碱基G对应的标志位加32(SNV对应的标记位为二进制的第六位,转化为十进制后为32)。对于短删除,记录删除起始和结束位点,由于短删除对应的标记位为二进制的第一位,转化为十进制为1,而变异状态标记位的初始值为2,所以同时对所有待删除碱基的标志位减1。而对于短插入,记录插入位置左侧第一个碱基的坐标,同时对该碱基的标记位加8(短插入对应的标记位为二进制的第四位,转化为十进制后为8)。
根据本发明的实施例,所述胚系变异包括指定胚系变异和随机胚系变异,所述体系变异包括指定体系变异和随机体系变异。由此,既包括指定胚系变异,又包括随机胚系变异,变异具有多样性,更符合人类基因特点。其中,需要说明的是,指定变异可以来源于目标疾病的数据库,从而使模拟得到的序列具有目标疾病的变异特征,而随机变异由随机产生,使变异具有个体特征,并且,无论指定变异还是随机变异,都分别具有体系变异和胚系变异,使变异即具有目标疾病的一般变异特征,又具有个体的变异特征,变异模拟序列的仿真程度更高。
S300生成含有胚系和体系变异的序列
根据本发明的实施例,基于目标疾病的体系变异数据库,在所述胚系变异序列的基础上迭代添加体系变异生成至少一条含有胚系和体系变异的序列。其中,“体系变异”也就是体细胞变异。其中,需要说明的是,添加体系变异的方法与前述添加胚系变异的方法相似,在此不再赘述。
根据本发明的实施例,生成所述含有胚系和体系变异的序列包括:将所述添加胚系变异后的特征串进行体系变异状态标记,以便得到体系标记后的特征串;选取待添加的体系变异;将所述待添加的体系变异整合至所述体系标记后的特征串上,以便得到添加体变异后的特征串和胚系变异信息;以及将所述添加体系变异后的特征串进行碱基还原,以便得到所述含有胚系和体系变异的序列。由此,在添加胚系变异后的特征串添加体系变异,形成含有胚系和体系变异的序列。
根据本发明的实施例,所述胚系变异状态标记和所述体系变异状态标记均包括:将碱基序列上的碱基赋予若干标记位,每个变异类型特异对应至少一个标记位,可以通过在各变异类型对应的标记位是进行有和无的标记来确认各碱基的变异类型,并且,各标记位是独立的,也可以是相互关联的,也就是说,可以同时对一个碱基的不同标记位进行标记,从而进行多层次变异的叠加,同时,针对变异模拟的特点,对部分变异类型进行关联,使其不能同时发生,例如同一碱基上SNV(单核苷酸突变)和删除变异不同时发生;通过有和无对字符的各标记位进行变异类型的设置,得到设置后的标记位,同一碱基可以存在多种变异类型,从而实现碱基任意层次的指定变异;将设置后的标记位进行数值转换,得到数值化的碱基的变异类型标记。
根据本发明的实施例,所述胚系变异和所述体系变异的类型为选自单核苷酸位点变异、删除、插入、连续重复、倒置、易位和拷贝数变异中的至少一种,其中,删除还包括长删除和短删除,插入还包括长插入和短插入,而短的掺入和删除,也就是插入缺失标记(InDel)。
根据本发明的实施例,所述体系变异含有目标疾病的变异。由此,该添加的变异含有目标疾病的变异。也就是,通过添加目标疾病的变异,使变异模拟后的序列含有目标疾病的变异特征,与目标疾病的序列的相似度更高。
进一步需要说明的是,在指定变异的选取过程中,可以检查各变异之间是否存在互斥关系,也就是在同一碱基位点是否可以发生多层的变异,该位点的变异类型是否会影响后续该位点的其他变异的设置,例如在删除位点再进行其他变异就没有意义。
为了便于理解变异整合到特征串上的方法,根据本发明的实施例,对该整合进行解释说明,体系变异和胚系变异均包括指定变异整合和随机变异整合,具体地,整合包括:将指定变异与变异标记特征串上的对应的标记位进行第一整合,得到添加指定变异的特征串;将随机变异与添加指定变异的特征串上随机选取的字符进行第二整合,并相应修改随机选取的字符的变异类型标记,得到变异特征串。具体地,第一整合是按照指定变异的原发位点的位置在变异标记特征串上的对应的标记位添加变异,而随机变异是在添加指定变异的特征串上随机生成的,位置和变异的类型均是随机产生的,产生随机变异后,记录该随机变异的变异位置、变异类型和具体变异序列信息等。
其中,需要说明的是,在进行第一整合时,需要将第二整合中的指定变异的标记位标记为不可以变异位点,例如,可以将标记位的数值设置为0,以防止第一整合中的随机变异在第二整合中的指定位点产生变异,影响第二整合指定变异的添加。
一般情况下,变异整合的顺序根据变异片段的长度进行设置,先进行变异片段较长的相应的变异类型的模拟,再进行变异片段较短的变异类型的模拟,例如变异模拟的顺序可以为长删除、连续重复、倒置、插入缺失标记和单核苷酸位点变异。
根据本发明的实施例,将添加变异后的特征串进行碱基还原,得到变异模拟后的序列。其中,碱基还原指的是将特征串上的非碱基特征符号转化为对应的碱基,从而得到模拟后的碱基序列。
在当前的序列模拟工具中,均是从参考基因组序列模拟变异的fasta序列。并且,在模拟变异的过程中,需要在fasta序列模拟步骤便进行校验,以保障此步骤模拟的正确性。进而,根据本发明的实施例,该方法进一步包括:变异校验装置,所述变异校验装置与所述胚系变异序列生成装置、所述体系变异序列生成装置和所述二倍体生成装置相连,所述变异校验装置用于对所述含有胚系和体系变异的序列和所述胚系变异序列进行变异校验,以便得到校验后的含有胚系和体系变异的序列和校验后的胚系变异序列。
根据本发明的实施例,所述变异校验装置包括:序列读取件,所述序列读取件用于顺序读取变异序列,所述变异序列为所述含有胚系和体系变异的序列和所述胚系变异序列;变异区确认件,所述变异区确认件与所述序列读取件相连,所述变异区确认件基于所述体系变异信息和所述胚系变异信息确定所述变异序列上的每一变异区域;区域划分件,所述区域划分件与所述变异区确认件相连,所述区域划分件基于所述每一变异区域确定所述变异序列中交替排列的对应于非变异区域的非变异字符串和对应于变异区域的变异字符串;序列校验件,所述序列校验件与所述区域划分件相连,按照所述变异序列的读取顺序,对于所述非变异字符串和所述变异字符串,交替地基于所述原始序列和所述变异信息校验所述待校验序列,包括:(1)对于所述非变异字符串,对所述原始序列与所述待校验序列进行第一比对,以及(2)对于所述变异字符串,基于所述变异区域的所述变异信息和所述原始序列生成模拟变异段,并对所述模拟变异段和所述待校验序列进行第二比对。由此,基于所述原始序列和所述变异信息校验所述待校验序列以确定所述待校验序列的变异是否正确,增加了变异模拟过程的可靠性。
具体地,因为变异信息包括原始序列按照位置存储的所有变异区域的变异信息,在顺序读取所述原始序列时,可以基于所述变异信息确定所述原始序列上的每一变异区域。并且,由于所述变异区域在所述原始序列上分段排列,针对如上所述的以字符串形式的原始序列,可以基于所述每一变异区域确定所述原始序列中交替排列的对应于非变异区域的非变异字符串和对应于变异区域的变异字符串。也就是说,字符串的长度可以是根据突变位点进行分隔得到的,例如,两段突变序列之间的为非突变序列为一段字符串,两段突变序列分别为两段独立的字符串,从而,通过突变位点将一段序列分隔为非变异区域的非变异字符串和对应于变异区域的变异字符串。
接下来,可以按照所述原始序列的读取顺序,对于所述非变异字符串和所述变异字符串,交替地基于所述原始序列和所述变异信息校验所述待校验序列,包括:对于所述非变异字符串,对所述原始序列与所述待校验序列进行第一比对,以及,对于所述变异字符串,基于所述变异区域的所述变异信息和所述原始序列生成模拟变异段,并对所述模拟变异段和所述待校验序列进行第二比对。
其中,所述第一比对和所述第二比对均是以字符串为单位进行的。由此,相当于碱基的逐一比对,本发明通过对非变异字符串和变异字符串分别按照字符串进行比较,其中,非变异字符串直接将原始序列上的对应位置的序列与待校验序列上的相应非变异字符串进行比较,而变异字符串需要先将原始序列上的对应位置的序列添加变异信息,生成变异序列,再将生成的变异序列与变异字符串进行比较。由此,基于每个变异起始结束位点切割并处理每条序列为多字符串,再基于字符串比较的方式,将最大效率的快速检查每条序列上每个位置的碱基的正确性。
S400生成二倍体
根据本发明的实施例,从所述胚系变异序列和所述含有胚系和体系变异的序列中选取部分变异作为纯合变异,分别以所述胚系变异序列和所述含有胚系和体系变异的序列为模板,得到胚系变异二倍体和含有胚系和体系变异的二倍体。换句话说,就是在含有胚系变异的序列的基础上添加体系变异,使序列同时具有胚系变异和体系变异,然后从胚系变异和体系变异选取部分变异作为纯合变异,合成二倍体序列,该二倍体序列为配对的二倍体序列,各二倍体序列构成二倍体集合。由此,通过在胚系变异数据库和体系变异数据库中选取变异,在测序序列上利用前述的对测序序列进行变异模拟的方法生成胚系变异和体系变异,使测序序列获得体系和胚系变异特征,并特异选取部分变异实现杂合和纯合二倍体,得到不同的纯度目标疾病患者测序文库。由此,该方法可以根据需要进行具有特定体系和胚系变异特征和预定纯合比例的测序文库的模拟,并且得到的测序文库的仿真程度高,有利于对相关测序装置的变异检测准确度的检测。对应系列中整合添加的变异,可以根据具体的需求进行调整,根据本发明的一些实施例,按照变异类型的预定比例,从胚系变异序列集合和含有胚系和体系变异的序列集合中选取变异。
其中,纯合变异的比例可以根据需要进行设定,例如,可以按照临床目标疾病的通常纯合变异的比例进行设定,以使变异序列集合的仿真度更高。
S500生成下机数据
根据本发明的实施例,将所述胚系变异二倍体和所述含有胚系和体系变异的二倍体转换为下机序列,得到胚系下机序列集和含有胚系和体系变异的下机序列集。例如,可以利用ART软件进行下机序列转换。
根据本发明的一些具体实施例,可以ART软件将二倍体的fasta数据生成对应的fastq测序文件,也就是下机序列。由此,利用ART软件生成对应的fastq测序文件,并模拟好接近真实的深度分布。
S600深度拟合
根据本发明的实施例,分别对所述胚系下机序列集和所述含有胚系和体系变异的下机序列集中的所述下机序列进行深度拟合,以便得到拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集。
根据本申请实施例的深度拟合的方法包括:获取参比序列和待拟合序列;基于所述目标序列的碱基对应的所述参比序列的基准参数和所述待拟合序列的基准参数,获得所述碱基的比较参数,其中,所述待拟合序列的基准参数基于多个碱基获得;基于所述比较参数确定所述碱基是否合格以获得合格碱基和不合格碱基;对所述不合格碱基的所述比较参数进行迭代消减以获得目标参数,所述迭代消减的迭代次数基于所述比较参数,且每次迭代消减与至少一个碱基的至少一个比较参数相关联;基于所述合格碱基和所述不合格碱基的所述目标参数确定所述目标序列的拟合曲线。
具体地,在本申请实施例中,参比序列是指目标区域真实测序fastq数据,并且,可以用samtools depth模块处理所述目标区域真实测序fastq数据得到目标区域测序深度分布,目标区域内碱基的真实深度即所述参比序列的基准参数。
另外,可以从参考序列上获得目标区域fasta序列,例如,参考序列可以是hg19序列。此外,针对高通量测序的特点,可以在目标区域的至少一侧,例如,在目标区域的两侧延长若干单位,例如200单位。接下来,目标区域fasta序列输入ART生成fastq数据,得到目标区域的深度分布,即待拟合序列的基准参数。
这里,在所述目标区域的至少一侧添加延长段的目的是:程序初始化阶段需要从ART得到的基于目标fasta序列的sam文件,且希望在目标区域内的深度波动尽可能水平,但是由于算法和序列模拟的限值,得到的数据深度波动呈现出两端较低的,而中间部分较为水平,类似一个梯形,不符合对初始数据的要求;因此,加入延长段能使目标区域居于中间,从而使深度波动趋于水平,满足初始数据的要求。
在根据本申请实施例的对目标序列进行拟合的方法中,所述基准参数为碱基深度。
并且,在根据本申请实施例的对目标序列进行拟合的方法中,所述待拟合序列通过在所述目标序列的至少一侧延长预定单位而获得。
进一步地,基于所述目标序列的碱基对应的所述参比序列的基准参数和所述待拟合序列的基准参数,获得所述碱基的比较参数,其中,所述待拟合序列的基准参数基于多个碱基获得。
因此,在基于所述目标序列的碱基对应的所述参比序列的基准参数和所述待拟合序列的基准参数,获得所述碱基的比较参数中,可以直接获得所述目标序列的碱基对应的所述参比序列的基准参数和所述待拟合序列的基准参数之间的差值作为所述碱基的比较参数,以下称为D-value。当然,本领域技术人员可以理解,所述比较参数也可以是其它能够体现所述参比序列的基准参数和所述待拟合序列的基准参数之间的差异的参数。具体的,在本发明的一些实施例中,该差值可以是碱基的深度,也就是待拟合序列碱基的读段数与参比序列碱基的读段数的差值。
在根据本申请实施例的对目标序列进行拟合的方法用于深度拟合的情况下,其可以用于模拟数据的深度分布,且模拟数据用于检测工具的性能,深度影响测序数据突变的检出,靶向测序数据和一般高通量测序数据的深度波动存在较大差异,因此更准确的模拟数据深度分布能更准确反应工具的对于这类数据检测的性能。
根据本发明的实施例,基于所述比较参数确定所述碱基是否合格以获得合格碱基和不合格碱基的方法如下:可以设定一偏差值,例如记为std,然后将上述步骤获得D-value与std比较,如果-std≤D-value≤std,则认为所述碱基为合格碱基,否则,所述碱基为不合格碱基。这里,一般来说,std越小,拟合程度越高,同时,计算时间消耗也越大,经验上来说,0.1%-0.5%能在时间效率和拟合程度上取得平衡。
对所述不合格碱基的所述比较参数进行迭代消减以获得目标参数,所述迭代消减的迭代次数基于所述比较参数,且每次迭代消减与至少一个碱基的至少一个比较参数相关联。
如上所述,由于模拟的碱基的参数分布,例如深度分布与真实的碱基参数分布,例如深度分布之间存在差异,且这种差异是由于每个读段包含多个碱基,因此针对每个碱基的计算实际上是基于包含该碱基的全部多个读段,例如,如果覆盖该碱基的读段数为75,则对该碱基的计算是根据这75个读段所获得的。因此,通过消减具有差异的碱基所对应的读段,可以避免由于读段数目过多而造成的平均化,从而使得碱基的深度分布更加接近真实分布。因此,从这个角度来看,上述D-value与每个碱基应该删除的读段数目相关联。
在根据本申请实施例的深度拟合的方法中,对所述不合格碱基的所述比较参数进行迭代消减以获得目标参数包括:对所述目标序列的所述碱基进行评分以获得碱基分数;基于所述碱基分数对所述目标序列的所述碱基中的各个碱基进行排序以获得排序后的不合格碱基;以及,对所述排序后的不合格碱基的所述比较参数进行迭代消减以获得目标参数。
在根据本申请实施例的对目标序列进行拟合的方法中,基于所述碱基分数对所述目标序列的所述不合格碱基中的各个碱基进行排序以获得排序后的不合格碱基包括:按照所述碱基分数的由低到高的顺序对所述各个碱基进行排序以获得排序后的不合格碱基。
在这种情况下,所述不合格碱基按照与所述参比序列中的对应碱基的基准参数,即深度的差异进行排序,差异小的不合格碱基排在前面,而差异大的不合格碱基排在后面。因此,可以首先处理差异小的不合格碱基,并且结合下面描述,也将使得对于不合格碱基的整体的处理加快收敛速度,从而促进对目标序列的拟合,使拟合的曲线的拟合度更佳。
如上所述,在本申请实施例中,是要删除覆盖所述碱基A的多条读段中的一条或多条读段,从而使得覆盖碱基A的读段数目减小,以使得所计算的碱基A的参数,例如深度的波动加大。因此,类似于上述对于碱基的操作,首先对碱基A的各读段进行评分以获得读段分数,并按照所述读段分数排序,例如按照由高到低的顺序排序以进行迭代删除。
具体来说,迭代删除在满足下列条件之一的情况下终止:1)消减后的所述比较参数小于参数阈值,即,消减后的D-value在参数阈值,如上述std范围内,也就是说消减后的碱基的深度已经足够接近于参比序列中对应碱基的深度;以及2),所述读段覆盖的合格碱基数大于碱基阈值,即,如果所述读段已经覆盖足够多的合格碱基,则说明此读段应尽量保留,这里,在本申请实施例中,所述碱基阈值可设置为在读段范围内的自然数,例如,所述碱基阈值是大于等于三且小于等于五的自然数。
经过每次迭代消减之后,所述比较参数,即D-value的值减一,直到满足上述终止条件之后,所述D-value减去删除的读段的数目就可以获得所述目标参数。
因此,在根据本申请实施例的对目标序列进行拟合的方法中,所述排序后的所述不合格碱基的所述比较参数进行迭代消减以获得目标参数包括:将所述不合格碱基的各读段进行评分以获得读段分数;按照所述读段分数的由高到低的顺序对所述不合格碱基的各读段依次进行迭代删除,所述迭代删除在满足下列条件之一的情况下终止:消减后的所述比较参数小于参数阈值;以及,所述读段覆盖的合格碱基数大于碱基阈值;基于所述不合格碱基的所述比较参数减去删除的读段的数目以获得所述目标参数。
根据本申请实施例的对目标序列进行拟合的方法中的迭代删除的示例性过程的流程图。如图3所示,所述迭代删除的示例性过程包括:对所述目标序列的所述碱基进行评分以获得碱基分数;按照所述碱基分数的由低到高的顺序对所述各个碱基进行排序以获得排序后的不合格碱基;将所述不合格碱基的各读段进行评分以获得读段分数;按照所述读段分数的由高到低的顺序对所述不合格碱基的各读段依次进行迭代删除,其中,所述迭代删除在满足下列条件之一的情况下终止:消减后的所述比较参数小于参数阈值;以及,所述读段覆盖的合格碱基数大于碱基阈值;以及基于所述不合格碱基的所述比较参数减去删除的读段的数目以获得所述目标参数。
这里,需要注意的是,在对读段进行评分的过程中,由于每个读段覆盖多个碱基,因此需要首先设置碱基权重,再基于碱基权重获得读段权重。
也就是说,所述碱基的权重与所述碱基的比较参数,例如深度与所述参比序列中相应碱基的深度之间的差异成正比,且所述读段的权重与所述读段所覆盖的多个碱基的差异成正比。换句话说,所述读段覆盖的所有碱基与所述参比序列中相应碱基之间的差值越大,所述读段的权重越大。
因此,在根据本申请实施例的对目标序列进行拟合的方法中,所述不合格碱基的评分为所述不合格碱基的所述比较参数与所述基准参数的商与一之和。
并且,在根据本申请实施例的对目标序列进行拟合的方法中,所述读段分数为所述读段覆盖的所有碱基的碱基分数的平均值;以及,所述不合格碱基的各所述读段按照所述读段分数的由高到低的顺序进行所述迭代删除。
因此,为了促进对目标序列的排序,在本申请实施例中,在每次对所述排序后的不合格碱基的所述比较参数进行迭代消减之后,都基于所述消减后的所述比较参数对所述不合格碱基重新计算评分并排序,并且,所述重新计算的评分是计算前的评分与所述迭代次数的差值。
可以看到,通过调整不合格碱基的权重,相应的覆盖所述不合格碱基的权重的读段的读段分数也会改变。另外,可以理解的是,由于迭代之后碱基的碱基分数降低,使得在计算读段分数时,涉及碱基分数降低的碱基数目更多的读段,也就是与当前删减的读段关联性大的读段的读段分数降低更多,从而排序向前,优先计算。
这样,如果基于消减之后的比较参数,所述碱基仍然是不合格碱基,则进行迭代消除,直到所有碱基均满足迭代终止的条件为止,其中,迭代终止的条件也就是前述所述的:(1)消减后的所述比较参数小于参数阈值;(2)所述读段覆盖的合格碱基数大于碱基阈值,满足这两个条件之一即终止。当碱基因满足条件(2)而终止消减时,该碱基可以称为本发明实施例的拟合曲线中的“坏点”,在消减过程中必然会存在一部分无法通过的坏点,在此,设置一个参数,当总的合格的碱基比例超过某个阈值(如70%)的时候,可以直接终止迭代消减,也可以待全部不合格碱基都完成消减后,再拟合曲线,如果全部碱基消减结束后,总的合格的碱基比例仍然低于该阈值时,则该参比序列拟合失败,可以将该参比序列重新利用ART软件生成新的具有不同深度波动的序列进行重新拟合。最后,基于所述合格碱基和所述不合格碱基的所述目标参数确定所述目标序列的拟合曲线。
S700生成目标疾病仿真测序文库
根据本发明的实施例,按预定比例从所述拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集选项部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列,所述部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列构成所述目标疾病仿真测序文库。
其中,本文中的“预定比例”可以根据不同的疾病、人群和测试需要等多方面的因素由本领域技术人员自行设置,通常以接近实际目标疾病人群的患者与正常人的比例为宜。例如,以构建肿瘤患者的仿真测序文库,可以按肿瘤/正常的比例从normal及tumor的fasta文件中随机抽取reads,抽取完成之后将reads合并成新的fastq文件,以实现不同的肿瘤纯度,其中,新的fastq文件即构成了仿真测序文库。
模拟目标疾病仿真测序文库的系统
根据本发明的另一方面,本发明提供了一种模拟目标疾病仿真测序文库的系统。根据本发明的实施例的模拟目标疾病仿真测序文库的系统,能根据需要得到不同体系和胚系变异特征、杂合/纯合比例和不同患病序列纯度,并能准确地模拟出接近真实的捕获测序条件下目标区域的深度波动,从而更加真实地模拟出了捕获测序条件下的下机数据。并且,该系统既能模拟全基因组下机数据,也能模拟捕获测序下机数据,适用范围广。同时,该系统运行速度快,能够在较短的时间内生成所需的模拟序列,并且模拟得到的序列的仿真程度高。
根据本发明实施例的模拟目标疾病仿真测序文库的系统,能够在较短的时间内生成所需的含有变异的序列,比如在915个总长约200kb的目标区域内随机模拟200个单核苷酸位点变异(SNV)和300个插入缺失标记(InDel),并同时生成对应的变异序列用时<30s。此外,迭代的层次约多,每次迭代的平均用时会更少,比如一次序列生成需要约21s,而4次迭代也不会超过22s。
参考图2,根据本发明的实施例,对该系统进行解释,该系统包括:
碱基序列获取装置100:根据本发明的实施例,该碱基序列获取装置100用于获取碱基序列。根据本发明的实施例,该碱基序列可以是捕获测序得到的。
根据本发明的实施例,该碱基序列获取装置100为序列片段截取装置。由此,序列截取的速度快。
胚系变异序列生成装置200:根据本发明的实施例,该胚系变异序列生成装置200与所述碱基序列获取装置100相连,基于胚系变异数据库,用于在所述碱基序列上添加胚系变异模拟生成胚系变异序列。
参考图3,根据本发明的实施例,该胚系变异序列生成装置200包括:胚系变异状态标记件210、胚系变异选取件220、胚系变异整合件230和第一碱基还原件240,其中,胚系变异状态标记件210用于将所述碱基序列进行胚系变异状态标记,得到胚系变异标记后的特征串;胚系变异选取件220用于选取待添加的胚系变异;所述胚系变异整合件230与所述胚系变异状态标记件210和所述胚系变异选取件220相连,用于将所述待添加的胚系变异整合至所述胚系变异标记后的特征串上,以便得到添加胚系变异后的特征串;所述第一碱基还原件240与所述胚系变异整合件230相连,用于将所述添加胚系变异后的特征串进行碱基还原,以便得到胚系变异序列。由此,通过对碱基序列的变异状态进行标记,设定碱基的变异类型,从而对各种变异进行模拟,变异模拟的方法简单,生成速度快,并能根据需要设计特殊的变异组合用于测试,变异模拟后的序列的仿真程度高。
体系变异序列生成装置300:根据本发明的实施例,该体系变异序列生成装置300与胚系变异序列生成装置200相连,基于目标疾病的体系变异数据库,该体系变异序列生成装置300用于在所述胚系变异序列的基础上迭代添加体系变异生成至少一条含有胚系和体系变异的序列。
参考图4,根据本发明的实施例,该体系变异序列生成装置300包括:体系变异状态标记件310、体系变异选取320、体系变异整合件330和第二碱基还原件340,其中,该体系变异状态标记件310与所述胚系变异整合件200相连,用于将所述添加胚系变异后的特征串进行体系变异状态标记,得到体系标记后的特征串;体系变异选取320用于选取待添加的体系变异;所述体系变异整合件330与所述体系变异状态标记件310和所述体系变异选取件320相连,用于将所述待添加的体系变异整合至所述体系标记后的特征串上,得到添加体变异后的特征串;所述第二碱基还原件340与所述体系变异整合件330相连,所述第二碱基还原件340用于将所述添加体系变异后的特征串进行碱基还原,得到所述含有胚系和体系变异的序列。由此,在添加胚系变异后的特征串添加体系变异,形成含有胚系和体系变异的序列。
根据本发明的实施例,所述胚系变异状态标记件210和所述体系变异状态标记件310均包括:碱基标记模块、变异类型设置模块和数值转换模块,其中,碱基标记模块用于将所述碱基序列上的碱基赋予若干标记位,每个变异类型特异对应至少一个标记位,可以通过在各变异类型对应的标记位是进行有和无的标记来确认各碱基的变异类型,并且,各标记位是独立的,也可以是相互关联的,也就是说,可以同时对一个碱基的不同标记位进行标记,从而进行多层次变异的叠加,同时,针对变异模拟的特点,对部分变异类型进行关联,使其不能同时发生,例如同一碱基上SNV(单核苷酸突变)和删除变异不同时发生;所述变异类型设置模块与所述碱基标记模块相连,用于通过有和无对所述字符的各所述标记位进行变异类型的设置,得到设置后的标记位;所述数值转换模块与所述变异类型设置模块相连,用于将所述设置后的标记位进行数值转换,得到数值化的所述碱基的变异类型标记。
根据本发明的实施例,所述变异类型设置模块以二进制数值表示所述有和无,得到二进制表示的所述设置后的标记位;所述数值转换模块用于将所述二进制表示的所述设置后的标记位转换成十进制数字,所述十进制数字为所述碱基的变异类型标记。
根据本发明的实施例,该系统进一步包括变异校验装置800,所述变异校验装置800与所述胚系变异序列生成装置200、所述体系变异序列生成装置300和所述二倍体生成装置500相连,该变异校验装置800用于对所述含有胚系和体系变异的序列和所述胚系变异序列进行变异校验,得到校验后的含有胚系和体系变异的序列和校验后的胚系变异序列。
参考图6,根据本发明的实施例,所述变异校验装置800包括:序列读取件810、变异区确认820、区域划分件830、序列校验件840,其中,所述序列读取件810用于顺序读取变异序列,所述变异序列为所述含有胚系和体系变异的序列和所述胚系变异序列;所述变异区确认件820与所述序列读取件810相连,所述变异区确认件820基于序列中交替排列的对应于非变异区域的非变异字符串和对应于变异区域的变异字符串;所述序列校验件840与所述区域划分件830相连,按照所述变异序列的读取顺序,对于所述非变异字符串和所述变异字符串,交替地基于所述原始序列和所述变异信息校验所述待校验序列,该校验包括:对于所述非变异字符串,对所述原始序列与所述待校验序列进行第一比对,以及对于所述变异字符串,基于所述变异区域的所述变异信息和所述原始序列生成模拟变异段,并对所述模拟变异段和所述待校验序列进行第二比对。由此,基于所述原始序列和所述变异信息校验所述待校验序列以确定所述待校验序列的变异是否正确,增加了变异模拟过程的可靠性。
二倍体生成装置400:根据本发明的实施例,该二倍体生成装置400分别与所述胚系变异序列生成装置200和所述体系变异序列生成装置300相连,用于从所述胚系变异序列和所述含有胚系和体系变异的序列中选取部分变异作为纯合变异,分别以所述胚系变异序列和所述含有胚系和体系变异的序列为模板,得到胚系变异二倍体和含有胚系和体系变异的二倍体。由此,通过在胚系变异数据库和体系变异数据库中选取变异,在测序序列上利用前述的对测序序列进行变异模拟的方法生成胚系变异和体系变异,使测序序列获得体系和胚系变异特征,并特异选取部分变异实现杂合和纯合二倍体,得到不同的纯度目标疾病患者测序文库。由此,该方法可以根据需要进行具有特定体系和胚系变异特征和预定纯合比例的测序文库的模拟,并且得到的测序文库的仿真程度高,有利于对相关测序装置的变异检测准确度的检测。
下机序列转换装置500:根据本发明的实施例,该下机序列转换装置与所述二倍体生成装置相连,用于将所述胚系变异二倍体和所述含有胚系和体系变异的二倍体转换为下机序列,以便得到胚系下机序列子集和含有胚系和体系变异的下机序列子集。例如,可以利用ART软件进行下机序列转换。
根据本发明的实施例,该下机序列转换装置500为ART转换装置,可以ART软件将二倍体的fasta数据生成对应的fastq测序文件,也就是下机序列。由此,利用ART软件生成对应的fastq测序文件,并模拟好接近真实的深度分布。
深度拟合装置600:根据本发明的实施例,该深度拟合装置与所述下机序列转换装置相连,用于分别对所述胚系下机序列集和所述含有胚系和体系变异的下机序列集中的所述下机序列进行深度拟合,以便得到拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集。
参考图5,根据本发明的实施例,所述深度拟合装置600包括:序列获取件610、比较参数获取件620、碱基合格判断件630、迭代消减件640和变异拟合件650,具体地,序列获取件610用于获取参比序列和待拟合序列,其中,所述待拟合序列为所述胚系下机序列集和所述含有胚系和体系变异的下机序列集中的任意序列;所述比较参数获取件610与所述序列获取件620相连,用于基于所述参比序列的基准参数和所述待拟合序列的基准参数,获得碱基的比较参数,其中,所述待拟合序列的基准参数基于至少一个碱基获得;所述碱基合格判断件630与所述比较参数获取件620相连,用于基于所述比较参数确定所述碱基是否合格以获得合格碱基和不合格碱基;所述迭代消减件640与所述碱基合格判断件630相连,用于对所述不合格碱基的所述比较参数进行迭代消减以获得目标参数,所述迭代消减的迭代次数基于所述比较参数,且每次迭代消减与至少一个碱基的至少一个比较参数相关联;所述变异拟合件650与所述迭代消减件640相连,用于基于所述合格碱基和所述不合格碱基得到所述拟合后的胚系下机序列集和所述拟合后的含有胚系和体系变异的下机序列集。
进一步地,根据本发明的实施例,所述迭代消减件640包括:碱基评分模块、碱基排序模块和迭代消减模块,其中,碱基评分模块用于对所述目标序列的所述碱基进行评分以获得碱基分数;所述碱基排序模块与所述碱基评分模块相连,用于按照所述碱基分数的由低到高的顺序对所述各个碱基进行排序以获得排序后的不合格碱基;所述迭代消减模块与所述碱基排序模块相连,用于对所述排序后的不合格碱基的所述比较参数进行迭代消减以获得目标参数。
根据本发明的实施例,所述迭代消减模块640包括:读段评分子模块,用于将所述不合格碱基的各读段进行评分以获得读段分数;迭代删除子模块,所述迭代删除子模块与所述读段评分子模块相连,用于按照所述读段分数的由高到低的顺序对所述不合格碱基的各读段依次进行迭代删除,所述迭代删除在满足下列条件之一的情况下终止:消减后的所述比较参数小于参数阈值;以及所述读段覆盖的合格碱基数大于碱基阈值;目标参数获取子模块,所述目标参数获取子模块与所述迭代删除子模块相连,用于基于所述不合格碱基的所述比较参数减去删除的读段的数目以获得所述目标参数。
纯度模拟装置700:根据本发明的实施例,该纯度模拟装置与所述深度拟合装置相连,用于按预定比例从所述拟合后的胚系下机序列集和拟合后的含有胚系和体系变异的下机序列集选取部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列,所述部分拟合后的胚系下机序列和部分拟合后的含有胚系和体系变异的下机序列构成所述目标疾病仿真测序文库。
在本说明书的描述中,参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中,对上述术语的示意性表述不一定指的是相同的实施例或示例。而且,描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。
尽管已经示出和描述了本发明的实施例,本领域的普通技术人员可以理解:在不脱离本发明的原理和宗旨的情况下可以对这些实施例进行多种变化、修改、替换和变型,本发明的范围由权利要求及其等同物限定。