具体实施方式
下面参照附图对本发明进行更全面的描述,其中说明本发明的示例性实施例。本发明的示例性实施例及其说明用于解释本发明,但并不构成对本发明的不当限定。
以下对至少一个示例性实施例的描述实际上仅仅是说明性的,决不作为对本发明及其应用或使用的任何限制。
随着大规模并行测序技术成本的降低,使用计算机模拟方法与微阵列杂交技术相比,能够简单地增加测序深度即可提高检测敏感度,以便于检测更小拷贝数的变异;另外,在某些基因组上,使用计算机模拟方法的深度测序技术可以弥补样本不纯的缺陷;此外,双端加标签的测序片段可以检测到结构重排信息。
本发明的目的是提供检测拷贝数变异的方法和装置,旨在高通量测序数据背景下,通过生物信息学方法检测目标样本与对照样本间发生拷贝数变异的区域,以备后续分析。
本发明基于全基因组标准,针对各染色体检测拷贝数变异区域。目前所有拷贝数变异检测都只能检测缺失和扩增,本发明也只针对拷贝数变异的这两种类型,检测的原理为:如果染色体某一区域发生了拷贝数变异,则高通量测序时该区域的序列片段分布将发生变化,即,拷贝数缺失-序列密度将变小,拷贝数扩增-序列密度将变大。在本发明中,检测拷贝数变异就是寻找这些发生变化的区域边界,可以将此边界定义为变点(即,某个碱基的位置),该变点左右的序列片段的拷贝数比率存在极大差异,再以拷贝数比率为标准从变点间区域(即,变点与变点之间的区域,此处操作是将基因组分段,找出侯选区域)筛选出最终的拷贝数变异区域。
图1是本发明检测拷贝数变异的方法的一个实施例的流程示意图。
如图1所示,该实施例可以包括以下步骤:
S102,在对每个样本进行高通量测序后,对测序所得的每个样本序列分别进行过滤,以去除掉不合格的序列和接头序列,其中,样本包括目标样本(即,变异组织)和对照样本(即,正常组织);
具体地,对高通量测序后的样本序列进行过滤,去除不合格的序列及接头序列,其中,不合格序列可以为下列情况中的至少一种:测序质量低于某一阈值的碱基个数超过整条序列碱基个数的一定比例(例如,50%)和序列中测序结果不确定的碱基(例如,Illumina GA测序结果中的N)个数超过整条序列碱基个数的一定比例(例如,10%)。其中,高通量测序技术可以为Illumina GA或者HiSeq测序技术,也可以为现有的其他高通量测序技术,低质量阈值可以由具体测序技术和测序环境确定。
S104,将过滤后的每个样本序列分别比对到参考基因组序列,对比对后的每个样本序列分别进行筛选以得到唯一比对的样本序列,确定每个唯一比对的样本序列相对于参考基因组序列的位置信息,并对位置信息进行排序;
具体地:(1)首先可以通过任何一种短序列映射程序(例如,短寡核苷酸分析包(Short Oligonucleotide Analysis Package,SOAP))将过滤得到的每个样本序列(即,由多个测序片段数据构成的序列)分别比对到参考基因组序列(例如,人类基因组参考序列)得到每个样本序列在参考基因组上的位置情况;(2)然后,对比对结果进行一系列的筛选,例如,去除比对到多个位置的序列(因为这个序列已无法准确唯一的提供比对位置信息)、去除重复出现的序列(因为这些序列可能是由于前期实验引入的误差,如由测序错误引起,为使检测结果更加精准,故去除),以得到唯一比对的序列结果;(3)最后,选取对照样本N、目标样本T以及对照样本的子集N1和N2(将N随机分为近似均等的两部分N1和N2(N=N1+N2),并确保分割后的N1和N2均包含全部的基因组信息)的序列相对于参考基因组序列的比对位置信息,并将位置信息从小到大排序以提高后续处理的效率,其中,比对位置信息可以按染色体、碱基位置大小进行排序。
S106,根据唯一比对的对照样本序列确定初始变点阈值和候选变点阈值,并将这两个阈值作为计算目标样本和对照样本间变异边界的标准,这样不仅可以降低噪音、减小测序误差对拷贝数变异检测的影响,而且还可以通过设定这两个阈值控制检测拷贝数变异的标准。
S108,根据初始变点阈值确定唯一比对的目标样本序列与唯一比对的对照样本序列间的目标样本拷贝数变异变点初始集。
S110,根据候选变点阈值通过迭代法合并目标样本拷贝数变异变点初始集,以得到目标样本拷贝数变异候选区域。
S112,根据设定的拷贝数比率阈值从目标样本拷贝数变异候选区域中过滤出目标样本的拷贝数变异区域;
举例说明,以目标样本拷贝数变异候选区域(即,变点间区域)的拷贝数比率为过滤标准,过滤出最终的目标样本拷贝数变异区域:严格的过滤标准,例如,拷贝数比率大于1.5为扩增,拷贝数比率小于0.5为缺失;宽松的过滤标准:例如,拷贝数比率大于1.3为扩增,拷贝数比率小于0.7为缺失,同时还将小于1000碱基对长度的序列片段过滤掉。
该实施例基于全基因组标准化,针对各条染色体寻找拷贝数变异,能够更加凸显各染色体的局部特征,为后续基因关联分析提供支持。
在本发明方法的另一实施例中,根据唯一比对的对照样本序列确定初始变点阈值和候选变点阈值的步骤可以包括:
步骤一,将唯一比对的对照样本序列随机分为两部分,并确保分割后的第一对照样本序列和第二对照样本序列均包含全部的基因组信息。
步骤二,根据设定的变点初始集参数确定第一对照样本序列和第二对照样本序列间的对照样本拷贝数变异变点初始集;
具体可以包括:
(1)将第一对照样本序列和第二对照样本序列比对到参考基因组序列,根据第一对照样本序列或第二对照样本序列确定与参考基因组比对上的起始位点集,如图2所示,分别示出了N1和N2均能比对上的第3条染色体的第200号碱基,以及N1比对上的第700号碱基,由于在N2上不一定能够找到这个700号碱基,所以经计算选取最接近700号碱基的一点,记录下每个样本序列的起始位点的位置信息,这个起始位点即为b点;
(2)根据设定的窗长滑动计算起始位点集中每个位点的拷贝数比率差异对应的P值,其中,拷贝数比率R的计算方法如下:
其中,xL,xR通过W(利用模拟数据确定W的最优值)确定,表示窗口的左右位置范围,T(xL,xR)表示目标样本T在该窗口内的序列数目,aT表示目标样本序列总数,N(xL,xR)表示对照样本N在该窗口内的序列数目,aN表示对照样本序列总数,所以拷贝数比率为在全基因组范围标准化之后的拷贝数比率,该拷贝数比率服从变化后的正态分布,可以用P值(即,概率)表示每个位点左右窗口的拷贝数比率差异,由于P值服从正态分布,其取值范围是0≤P≤1,所以P值越小表明拷贝数比率差异越大;
如图3所示,假设每个窗口有3000个(即,W=3000)序列,则有3000个b点,假设样本序列中一共有x个序列,则需要计算的P值个数为x-2W个,其中,由于N前端的W个序列上的起始位点和末端的W个序列由于没有前后窗口的比较,所以无法计算相应的P值。
(3)对起始位点集中每个位点对应的P值按从小到大的顺序排序;如图4所示,计算b0、b1、...bx-2W点对应的P值,进行排序并选出1000个点。
(4)根据设定的变点初始集参数按照P值从小到大(即,拷贝数比率差异从大到小)的顺序选取位点,每选取一个位点,将所选取位点的左右窗口中的所有位点对应的P值设置为1,并将选取的位点作为对照样本拷贝数变异变点初始集中的变点元素;如图5所示,在选取出的P值最小的1000个c点中,利用迭代方法删选出10个P值最小的变点,假设c3是P值最大的变点,首先去掉c3的信息,此时c2和c4间的窗口改变,再重新计算c2和c4点对应的P值。
步骤三,将对照样本拷贝数变异变点初始集中最小拷贝数比率差异对应的P值设置为初始变点阈值,并去除第一对照样本序列和第二对照样本序列中未被对照样本拷贝数变异变点初始集选中的位点信息。
步骤四,根据设定的变点集参数通过迭代法合并对照样本拷贝数变异变点初始集;
具体可以包括:
(a)在对照样本拷贝数变异变点初始集中,以各变点间区域为窗口计算第一对照样本与第二对照样本间的拷贝数比率差异对应的P值;
(b)将对照样本拷贝数变异变点初始集中每个变点对应的P值进行排序,去除P值最大(即,拷贝数比率差异最小)的变点信息,重新计算所去除的P值最大的变点的左右两个变点对应的P值,再去除对照样本拷贝数变异变点初始集中所剩变点中P值最大的变点信息,直至对照样本拷贝数变异变点初始集中所剩变点个数满足设定的变点集参数。
如图6所示,用上述方法继续比较c2点左右窗口以后得出新的P值,同理,得出c4点对应的新P值,这样,在剩下的999个点中继续寻找最大P值对应的变点,依此类推,最终将这1000个变点逐渐减少合并为10个。
步骤五,将合并后的对照样本拷贝数变异变点集中最小拷贝数比率差异对应的P值设置为候选变点阈值。
以下通过一个具体实例来进一步说明:
(1)将对照样本N随机分为近似均等的两部分N1和N2,并确保分割后的N1和N2均包含全部的基因组信息;
(2)以单个样本比对到参考基因组序列的起始位置b为起始点,其中,b为每个序列比对上的起始位点,每个测序片段分别进行比对,在整个基因组上存在多个比对上的b点,每个测序片段比对上一个b点;b点的确定方法为:对照样本N、目标样本T与参考基因组序列均比对上的位点,可以遵循SOAP比对原则,b点即为进行SOAP比对后得到的比对上的位点;
(3)确定b点左右两侧窗口大小,例如,窗长的条件为:以b点为一端包含相同的测序片段数目(W),可以用模拟数据确定W的最优值;
(4)比较对照样本N1与N2间的拷贝数比率在b点左右两窗口的差异,对照样本N1与N2进行比较时,选取其中任意一个确定起始位置和窗口位置,例如,选取对照样本N1确定起始位置和窗口位置,如果W的值为400,则在N1的起始位置就在第401个序列片段所在的b点,对应到N2,不一定能找到与N1相对应的点(这依据SOAP比对的结果),此时需要通过计算找到与N1的第401个序列片段所在b点最接近的点;
(5)每个b点(基于全基因组测序中各染色体第一个比对上的b点,每次向右滑动一个测序片段,找到下一个b点,直至整个基因组上的所有b点都找到)对应一个P值,按P值由小到大选b点,每选一个b点,将该b点左右窗口内所有b点对应的P值置1,即,该窗口的显著性用最显著的那一点来代替,在N1与N2进行比较时,用这样的方法取前i个点(即,P值从小到大排列对应的前i个b点,其中,i可以是例如1000)作为对照样本拷贝数变异变点初始集(在该变点初始集中,保留这i个b点的信息,去除其他b点信息,保留整个样本序列,此时,各个b点(即,变点)间距离可能变大,i个b点按样本序列中的相对物理位置排序),同时定义该变点初始集中最大的P值为p_bkp(即,p_bkp为初始变点阈值);
(6)在新的窗口(即,变点与变点间的窗口)中计算各个变点在更大区域里的显著性水平(整个样本序列的第一变点的左窗口从第一个测序片段开始直到第一个变点,整个样本序列的最后一个变点的右窗口为最后一个变点直至整个样本序列结束);
(7)将最大P值对应的变点信息删除,合并其左右两窗口,重新计算所删除的变点左右的两个变点对应的P值,以对照样本拷贝数变异变点初始集为起始点,初始集中各变点间区域为窗口,比较对照样本N1与N2间变点左右窗口内序列拷贝数比率的差异,同样用P值表示差异大小,P值越小差异越大,当比较对照样本N1与N2时,以模拟数据为参照,设定对照样本拷贝数变异变点集合并后变点的个数为f,利用迭代的方法按P值由大到小的顺序依次减少变点,直至达到设定个数f,此时得到对照样本拷贝数变异变点集,定义该变异变点集中最大的P值为p_merge(即,候选变点阈值)。
该实施例通过可变窗口(“可变窗口”即固定的序列个数,但碱基数目不一定相同,进而窗口大小不同,从而得到可变窗口)形式计算变点能够较准确地找到拷贝数变异的区域边界,通过简单地增加测序深度即可提高检测灵敏度,得到更加准确的拷贝数变异区域边界,同时能够得到更小的拷贝数变异。此外,该实施例通过比较两份对照样本(即,N1和N2)来设定阈值参数,不仅可以减小测序误差对寻找拷贝数变异的影响,而且还可以通过参数的设定控制检测拷贝数变异的标准。
在本发明方法的又一实施例中,根据初始变点阈值确定唯一比对的目标样本序列与唯一比对的对照样本序列间的目标样本拷贝数变异变点初始集的步骤可以包括:
步骤一,将唯一比对的目标样本序列和唯一比对的对照样本序列比对到参考基因组序列,根据唯一比对的对照样本序列确定与参考基因组比对上的起始位点集;
步骤二,根据设定的窗长滑动计算起始位点集中每个位点的拷贝数比率差异对应的P值,具体方法同上述,在此不再重复;
步骤三,起始位点集中,选取拷贝数比率差异对应的P值小于初始变点阈值的位点,将所选出的位点构成的集合作为目标样本拷贝数变异变点初始集,并去除唯一比对的对照样本序列和唯一比对的目标样本序列中未被选取的位点信息。
以下通过一个具体实例来说明:
(1)以单个样本比对到参考基因组序列的起始位置b为起始点,其中,b点的确定方法为:对照样本N、目标样本T与参考基因组序列均比对上的位点;
(2)确定b点左右两侧窗口的大小,例如,窗长的条件为:以b点为一端包含相同测序片段数目(W),可以用模拟数据确定W的最优值;
(3)比较对照样本N与目标样本T之间的拷贝数比率在b点左右两窗口的差异,对照样本N与目标样本T进行比较时,选取对照样本N确定起始位置集和窗口位置;
(4)用P值表示b点左右窗口的拷贝数比率差异,由于P值服从正态分布,所以P值越小表明拷贝数比率差异越大,每个b点对应一个P值,将P值(假设k个)按照从小到大的顺序排序,选取出其中最大P值对应的变点ki,将其信息删除,再计算ki-1和ki+1这两个变点所对应的P值,此时将仅剩下k-1个变点,在k-1个变点中,选取出其中最大P值对应的变点kj,将其信息删除,再计算kj-1和kj+1这两个变点的P值,此时将仅剩下k-2个变点,以此类推,直至剩下的所有变点的P值均小于p_bkp。如果没有变点的P值小于p_bkp,则认为在此情况下,该染色体没有发生拷贝数变异。
该实施例利用初始变点阈值确定目标样本拷贝数变异变点初始集可以减小测序误差对寻找拷贝数变异的影响。
在本发明方法的又一实施例中,根据候选变点阈值通过迭代法合并目标样本拷贝数变异变点初始集,以得到目标样本拷贝数变异候选区域的步骤可以包括:
步骤一,在目标样本拷贝数变异变点初始集中,以各变点间区域为窗口计算唯一比对的目标样本与唯一比对的对照样本间的拷贝数比率差异对应的P值;
步骤二,将目标样本拷贝数变异变点初始集中每个变点对应的P值进行排序,去除P值最大的变点信息,重新计算所去除的P值最大的变点的左右两个变点对应的P值,再去除目标样本拷贝数变异变点初始集中所剩变点中P值最大的变点信息,直至所剩变点中对应的最大P值小于候选变点阈值或不存在变点。
以下通过一个具体实例来说明:
(1)在目标样本拷贝数变异变点初始集的新窗口(即,变点与变点间的窗口)中计算各个变点在更大区域里的显著性水平,其中,整个样本序列的第一变点的左窗口从第一个测序片段开始直到第一个变点,整个样本序列的最后一个变点的右窗口为最后一个变点直至整个样本序列结束;
(2)将计算出的P值(假设m个)按照从小到大的顺序排序,选取出其中最大P值对应的变点mi,将其信息删除,再计算mi-1和mi+1这两个变点所对应的P值,此时将仅剩下m-1个变点,在m-1个变点中,选取出其中最大P值对应的变点mj,将其信息删除,再计算mj-1和mj+1这两个变点对应的P值,此时将仅剩下m-2个变点,以此类推,直至剩下的所有变点的P值均小于p_merge或不再有变点为止。
以下以肝癌癌症组织样本及其癌旁组织样本为例说明拷贝数变异的检测过程。
图7是本发明实施例提供的基于二元分割算法的拷贝数变异检测流程示意图,详述如下:
S202,采用Illumina hiseq高通量测序技术对待测序列进行测序,接收hiseq测序序列后,对测序序列进行过滤,去除不合格的序列,将样本接头序列从序列片段中去除,其中,不合格序列包括:测序质量值低于5的碱基个数超过整条序列碱基个数的50%或序列中测序结果中N的个数超过整条序列碱基个数的10%;
S204,采用短寡核苷酸分析包(SOAP)映射程序将高通量测序技术得到的癌症样本和癌旁样本测序片段比对到人类参考基因组序列上,筛选掉比对结果中多重比对的序列(即,比对结果中第四列不为1的序列),去除重复出现的测序序列(即,两端测序数据中不同序列编号而碱基序列一样的序列,(仅保留一对两端测序序列)),一端测序数据中碱基序列一样的序列(仅保留一份),以降低结果中的假阳性,最后,根据需求提取比对结果中所需要处理的染色体编号及位置信息,染色***置信息的获取方式:将癌旁样本序列数据分为两份(N1和N2),分别提取各染色体的位置信息,再将N1、N2中各染色体的位置信息合并为N的相应染色***置信息,同时提取患病样本中的各染色体的位置信息,将所有位置信息从小到大排序,信息格式为1.txt,2.txt....22.txt,其内容为该染色体的位置信息;
S206,确定初始变点阈值p_bkp和候选变点阈值p_merge,在用两份癌旁样本数据确定阈值之前,首先用模拟数据确定参数W、i、f的最优值(例如,W=3000,i=1000,f=10),即设定可变窗口包含的测序片段数目(3000个),选择i个初始变点(即,初始化1000个变点)并最终将它们合并为f(即,10)个变点;
以癌旁样本序列比对到参考基因组上序列的起始位置b为起始点,比较两份癌旁样本序列相同起始点左右两侧窗口长度的差异,其中,窗口满足的条件为:以起始点为一端包含W个序列,两份癌旁样本间拷贝数比率的差异用P值表示,这样每个b点对应一个P值,按P值由小到大选b点,每选一个b点,将该点左右窗口内所有点的P值设为1,即该窗口的显著性用最显著的那一点来代替,用这样的方法取前1000个b点作为变点的初始集,其中最大的P值即为后续分析需要的初始变点阈值p_bkp;
合并变点:以这个变点初始集为起始点,初始集各点间区域为窗口,比较两份癌旁样本间左右窗口内序列片段拷贝数比率的差异,且其服从正态分布,同样用P值表示差异大小,P值越小差异越大,用迭代的方法按P值由大到小的顺序依次减少变点,直至设定个数10个,此时得到的变点集为癌旁样本间的变点集,其对应的最大P值即为后续分析需要的候选变点阈值p_merge;
S208,检测初始变点:利用与上一步分析比较两份癌旁样本数据(N1和N2)相同的方法,比较癌症样本T和癌旁样本N,得到初始变点集和对应的P值,选取P值小于初始变点阈值p_bkp的变点合并,拷贝数变异的边界之处,其左右分布必存在差异,故在i=1000标准下(即,p_bkp门限下),这些候选变点即为可能的拷贝数变异边界;
S210,合并初始变点:同样用迭代的方法合并变点,选取P值小于候选阈值p_merge的变点作为癌症样本与癌旁样本间的拷贝数变异变点,通过合并变点,我们将局部的差异变点放在更大的范围内来检测这些变点在f=10标准下(即,p_merge门限下)是否为真正的CNV边界,从而筛选出准确的拷贝数变异区域及其边界;
S212,确定拷贝数变异区域:用拷贝数变异变点间区域的序列拷贝数比率判断拷贝数发生缺失或扩增,从而过滤得到拷贝数变异区域,该实施例使用宽松过滤标准,拷贝数比率大于1.3为拷贝数扩增,拷贝数比率小于0.7为拷贝数缺失,此外,还过滤掉小于1k的片段的区域,最终所得的区域即为癌症与癌旁相比的拷贝数变异区域。
此外,还可以进一步对各拷贝数变异区域进行基因注释,关联分析,基因富集等分析。
图8示出了图7实施例检测出的22号染色体局部拷贝数发生变异的区域,其中,横坐标“chromosome 22position(Mb)”为“染色体22的位置(Mb)”,纵坐标“Copy number ratio”为“拷贝数比率”,过滤条件为拷贝数比率为1.3和0.7,拷贝数比率大于1.3的部分为拷贝数发生变异区域,拷贝数比率在0.7到1.3之间的部分为没有发生拷贝数变异的区域。
本发明在应用时,推荐使用低深度:W=600,i=2000,f=60;高深度:W=3000,i=1000,f=10,得到全基因组上各个区域的拷贝数比率等信息,按照拷贝数缺失或增加在不同梯度上进行过滤(1.5、0.5及1.3、0.7两个过滤标准,也可按照需求自行设定阈值),得到预测的拷贝数变异区域,以备后续分析疾病易感基因。
本发明的上述实施例能够利用生物信息学方法的高灵敏度和特异性快速地检测拷贝数变异。
图9是本发明检测拷贝数变异的装置的一个实施例的结构示意图。
如图9所示,该实施例的装置10可以包括:
序列过滤模块11,用于对测序所得的每个样本序列分别进行过滤,以去除掉不合格的序列和接头序列,其中,样本包括目标样本(即,变异组织)和对照样本(即,正常组织),不合格序列可以为下列情况中的至少一种:测序质量低于某一阈值的碱基个数超过整条序列碱基个数的一定比例(例如,50%)和序列中测序结果不确定的碱基(例如,IlluminaGA测序结果中的N)个数超过整条序列碱基个数的一定比例(例如,10%);
序列比对模块12,与序列过滤模块11相连,用于将过滤后的每个样本序列分别比对到参考基因组序列,对比对后的每个样本序列分别进行筛选以得到唯一比对的样本序列,确定每个唯一比对的样本序列相对于参考基因组序列的位置信息,并对位置信息进行排序;
具体地,可以通过任何一种短序列映射程序将过滤得到的每个样本序列分别比对到参考基因组序列,然后,对比对结果进行一系列的筛选,例如,去除比对到多个位置的序列、去除重复出现的序列,以得到唯一比对的序列结果,最后,选取对照样本N、目标样本T以及对照样本的子集N1和N2的序列相对于参考基因组序列的比对位置信息,并将位置信息从小到大排序以提高后续处理的效率,其中,比对位置信息可以按染色体、碱基位置大小进行排序;
阈值确定模块13,与序列比对模块12相连,用于根据唯一比对的对照样本序列确定初始变点阈值和候选变点阈值,并将这两个阈值作为计算目标样本和对照样本间变异边界的标准,这样不仅可以降低噪音、减小测序误差对拷贝数变异检测的影响,而且还可以通过设定这两个阈值控制检测拷贝数变异的标准;
目标样本变点初始集确定模块14,与序列比对模块12和阈值确定模块13相连,用于根据初始变点阈值确定唯一比对的目标样本序列与唯一比对的对照样本序列间的目标样本拷贝数变异变点初始集;
目标样本变点集确定模块15,与目标样本变点初始集确定模块14和阈值确定模块13相连,用于根据候选变点阈值通过迭代法合并目标样本拷贝数变异变点初始集,以得到目标样本拷贝数变异候选区域;
拷贝数变异确定模块16,与目标样本变点集确定模块15相连,用于根据设定的拷贝数比率阈值从目标样本拷贝数变异候选区域中过滤出目标样本的拷贝数变异区域,其中,可以采用严格的过滤标准,例如,拷贝数比率大于1.5为扩增,拷贝数比率小于0.5为缺失;也可以采用宽松的过滤标准:例如,拷贝数比率大于1.3为扩增,拷贝数比率小于0.7为缺失,同时还将小于1000碱基对长度的序列片段过滤掉。
该实施例基于全基因组标准化,针对各条染色体寻找拷贝数变异,能够更加凸显各染色体的局部特征,为后续基因关联分析提供支持。
图10是本发明检测拷贝数变异的装置的另一实施例的结构示意图。
如图10所示,与图9中的实施例相比,该实施例的装置20中的阈值确定模块21可以包括:
序列分割单元211,用于将唯一比对的对照样本序列随机分为两部分,并确保分割后的第一对照样本序列和第二对照样本序列均包含全部的基因组信息;
初始变点阈值确定单元212,与序列分割单元211相连,用于根据设定的变点初始集参数确定第一对照样本序列和第二对照样本序列间的对照样本拷贝数变异变点初始集,将对照样本拷贝数变异变点初始集中最小拷贝数比率差异对应的P值设置为初始变点阈值,并去除第一对照样本序列和第二对照样本序列中未被对照样本拷贝数变异变点初始集选中的位点信息;
具体地,可以将第一对照样本序列和第二对照样本序列比对到参考基因组序列,根据第一对照样本序列或第二对照样本序列确定与参考基因组比对上的起始位点集;根据设定的窗长滑动计算起始位点集中每个位点的拷贝数比率差异对应的P值;对起始位点集中每个位点对应的P值按从小到大的顺序排序;根据设定的变点初始集参数按照P值从小到大的顺序选取位点,每选取一个位点,将所选取位点的左右窗口中的所有位点对应的P值设置为1,并将选取的位点作为对照样本拷贝数变异变点初始集中的变点元素;
候选变点阈值确定单元213,与初始变点阈值确定单元212相连,用于根据设定的变点集参数通过迭代法合并对照样本拷贝数变异变点初始集,并将合并后的对照样本拷贝数变异变点集中最小拷贝数比率差异对应的P值设置为候选变点阈值;
具体地,可以在对照样本拷贝数变异变点初始集中,以各变点间区域为窗口计算第一对照样本与第二对照样本间的拷贝数比率差异对应的P值;将对照样本拷贝数变异变点初始集中每个变点对应的P值进行排序,去除P值最大的变点信息,以所去除的P值最大的变点的左右两个变点间区域为新窗口重新计算左右两个变点对应的P值,再去除对照样本拷贝数变异变点初始集中所剩变点中P值最大的变点信息,循环迭代以上合并变点步骤,直至对照样本拷贝数变异变点初始集中所剩变点个数满足设定的变点集参数。
该实施例通过计算符合一定要求的对照样本来设定阈值参数,将设定的阈值参数作为计算目标样本和对照样本间变异边界的标准,不仅可以降低噪音、减小测序误差对检测拷贝数变异的影响,而且还可以通过设定参数控制检测拷贝数变异的标准。
图11是本发明检测拷贝数变异的装置的又一实施例的结构示意图。
如图11所示,与图9中的实施例相比,该实施例的装置30中的目标样本变点初始集确定模块31可以包括:
起始位点集确定单元311,用于将唯一比对的目标样本序列和唯一比对的对照样本序列比对到参考基因组序列,根据唯一比对的对照样本序列确定与参考基因组比对上的起始位点集;
第一拷贝数比率差异计算单元312,与起始位点集确定单元311相连,用于根据设定的窗长滑动计算起始位点集中每个位点的拷贝数比率差异对应的P值;
目标样本拷贝数变异变点初始集确定单元313,与第一拷贝数比率差异计算单元相连312,用于在起始位点集中选取拷贝数比率差异对应的P值小于初始变点阈值的位点,将所选出的位点构成的集合作为目标样本拷贝数变异变点初始集,并去除唯一比对的对照样本序列和唯一比对的目标样本序列中未被选取的位点信息。
图12是本发明检测拷贝数变异的装置的再一实施例的结构示意图。
如图12所示,与图9中的实施例相比,该实施例的装置40中的目标样本变点集确定模块41包括:
第二拷贝数比率差异计算单元411,用于在目标样本拷贝数变异变点初始集中,以各变点间区域为窗口计算唯一比对的目标样本与唯一比对的对照样本间的拷贝数比率差异对应的P值;
变点迭代单元412,与第二拷贝数比率差异计算单元411相连,用于将目标样本拷贝数变异变点初始集中每个变点对应的P值进行排序,去除P值最大的变点信息,重新计算所去除的P值最大的变点的左右两个变点对应的P值,再去除目标样本拷贝数变异变点初始集中所剩变点中P值最大的变点信息,直至所剩变点中对应的最大的P值小于候选变点阈值或不存在变点。
上述检测拷贝数变异的装置的具体实例可以参照前述检测拷贝数变异的方法的具体实例,在此不再重复。
虽然已经通过示例对本发明的一些特定实施例进行了详细说明,但是本领域的技术人员应该理解,以上示例仅是为了进行说明,而不是为了限制本发明的范围。本领域的技术人员应该理解,可在不脱离本发明的范围和精神的情况下,对以上实施例进行修改。本发明的范围由所附权利要求来限定。