CN107077533B - 测序数据处理装置和方法 - Google Patents

Abstract

一种测序数据处理装置(100)，装置包括：数据接收单元(10)，用于接收所述测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每个读段都包含缺口；处理器(20)，用于执行数据处理程序，执行所述数据处理程序包括实现将所述测序数据与参考序列作比对，获得比对结果，以及消除所述比对结果中的每个读段的缺口，获得通用比对结果；以及，至少一个存储单元(30)，用于存储数据，其中包括所述数据处理程序。还提供了一种测序数据处理系统、方法，一种计算机可读存储介质、一种检测CNV的方法及设备。

Description

测序数据处理装置和方法

技术领域

本发明涉及生物信息领域，具体的，本发明涉及测序数据处理装置和方法，更具体地，本发明涉及一种测序数据处理装置、一种测序数据处理系统、一种测序数据的处理方法、一种计算机可读存储介质、一种检测CNV的方法以及一种CNV检测设备。

背景技术

存在于血清、血浆或其他体液中的cfDNA(cell-free DNA)，是一种有效的生物标记物，它可以应用于多种突变检测中，比如癌症、胎儿染色体变异等基因突变导致的疾病。由于缺乏高敏感度和准确度的定量分析技术，此前的大量研究关注的都是一些已知的疾病相关基因，比如色素瘤-GNAQ基因(Metz,Claudia HD,et al.Ultradeep sequencingdetects GNAQ and GNA11 mutations in cell-free DNA from plasma of patientswith uveal melanoma.Cancer medicine 2.2(2013):208-215.)，21三体-21号染色体(Liao,Gary JW,et al."Noninvasive prenatal diagnosis of fetal trisomy 21 byallelic ratio analysis using targeted massively parallel sequencing ofmaternal plasma DNA."PLoS One 7.5(2012):e38154.)等。

新一代测序技术454(Roche)、Solexa(Illumina)和SOLiD(ABI)等的诞生，使得测序通量迅速提升的而测序成本急剧下降，这为cfDNA检测提供了新的思路。目前大规模并行测序(Massively Parallel Sequencing,MPS)是最主流的cfDNA检测技术，它被广泛应用于血浆DNA分子诊断、胎儿染色体异倍体、全基因组核型分析，甚至胎儿全基因组测序中。拷贝数变异(Copy-Number Variations,CNV)是指在人类基因组中广泛存在的，从1000bp到数百万bp范围内的缺失、插入、重复和复杂多位点的变异。拷贝数变异是许多人类疾病(如癌症、遗传性疾病、心血管疾病)的重要生物标志，已成为许多疾病研究的热点。尤其是对肿瘤的拷贝数变异检测可以发现肿瘤DNA在整个染色体组的缺失或倍增。目前已有的CNV检测平台有基于大插入片段的比较基因组杂交(CGH)、代表性寡核苷酸微阵列分析(ROMA)等。这些平台对于小CNV(20kb以下)的检测能力不足，而且存在操作繁琐，成本高等问题。

发明内容

本发明旨在至少在一定程度上解决上述技术问题之一。

依据本发明的第一方面，本发明提出了一种测序数据处理装置，该装置包括：数据接收单元，用于接收所述测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读长对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读段对中的两个读段都来自所述染色体片段的正链或所述染色体的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂；处理器，用于执行数据处理程序，执行所述数据处理程序包括实现将所述测序数据与参考序列作比对，获得比对结果，以及消除所述比对结果中的每个读段的缺口，获得通用比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果；以及，至少一个存储单元，用于存储数据，其中包括所述数据处理程序。这里所说的分别来自一条染色体片段的两个位置的读段对，可以通过构建末端文库(pair-end library)或者配对末端文库(mate-pair library)，对所构建的文库进行测序来获得，在本发明的一个实施例中，利用Complete Genomics公司(CG)的文库构建方法及其测序平台，获得多对读段对，一对读段对之间的距离是由读段的长度以及酶的识别位点和切割位点的距离来控制的。CG平台通过酶切构建多接头配对末端文库，利用特有的组合探针连接测序(cPAL)技术对所构建的环状文库进行测序，测读出接头两旁的碱基，因为其是利用酶切连接一个接头的两段来进行配对末端文库构建的，由于每一种酶都有一个首选的切割距离，而在实际酶切时经常比首选距离多一个位置或少一个位置，这样使得读段中经常带有缺口(gap)，缺口常为+1或者-1，和/或，建库时倘若使用同一种酶多次酶切，每次的酶切位置易发生变化，酶切位置的变化也会使获得的读段带有缺口，例如在构建多接头环状文库时，利用Alu酶两次酶切来连接多个接头的不同部分，读测这些接头旁的碱基时，会产生带+3/-3的缺口的读段。在本发明中缺口的大小还可以是0。以CG平台目前的双接头(two adaptors,2-AD)测序文库为例，2-AD测序输出总长为60bp，可分为两对读段对(mate-paired reads)，每对读段对中的读段在10bp的位置都有小的gap，在20bp位置有一个无效测序位点N，一对读段对的两个读段之间的在基因组上的距离一般小于2000bp。来自多接头文库中的多个读段，一个读段可以和任一其它读段组成一对读段对。这里所说的“正链”和“负链”是组成染色体片段的互补的两条链，是相对的，称一条链为正链就可以称其互补链为负链，在本发明的一个实施例中，将与参考序列匹配的链称为正链，将另一条链称为负链。在本发明中，比对可以利用已知比对软件进行，比如SOAP、BWA等，也可以利用CG平台的比对软件TeraMap进行。在本发明的一个实施例中，比对是利用TeraMap进行的，所得的比对结果的格式为TeraMap。在本发明的一个实施例中，所说的消除比对结果中每个读段的缺口是指，对带负缺口的读段去除掉其负缺口即去除掉重叠的碱基，对带正缺口的读段以N替代正缺口的大小，N为A、T、C或G，例如，对带负缺口比如为-2nt的读段，基于缺口处该读段可分成两部分，两部分的末端有2nt重叠，比如该读段的两部分分别为ATCGCTTAAG和AGTACGATTC，消除其负缺口即重叠的AG，获得对应的读段为ATCGCTTAAGTACGATTC。

在本发明的一个实施例中，本发明的一方面的方法中的作比对包括：将每对读段对的左臂和右臂分别与所述参考序列比对，获得一级左比对结果和一级右比对结果；分别以所述一级左比对结果和所述一级右比对结果的其中一个为参考，对另一个作比对，获得二级左比对结果和二级右比对结果；基于所述二级左比对结果和所述二级右比对结果，获得多个所述读段对的比对结果，或者获得多个所述左臂的比对结果和多个所述右臂的比对结果。这样经过两次比对，可以获得读段对比对结果，在本发明的一个实施例中，第一次比对是与参考序列作全局比对，以该次左臂/右臂比对结果为基准对右臂/左臂比对结果进行的第二次比对为局部比对，这样，能够将分别来自二级左比对结果和二级右比对结果中的比对到同一染色体上且之间的距离符合预期的两个读段配对成一对读段对，获得读段对比对结果。

在本发明的一个实施例中，所说的作比对包括，设置所述缺口的大小以使每个左臂或者每个右臂与所述参考序列进行多次比对，以获得最佳比对结果。例如，将所述每个左臂或者所述每个右臂的缺口分别设置为-3nt、-2nt、-1nt、0nt、1nt、2nt、3nt、4nt、5nt、6nt和7nt，获得对应的多个读段，分别将所述对应的多个读段与所述参考序列比对，将最优比对的序列作为该左臂/右臂，这里对于比对结果的好坏可以基于所利用的比对软件对比对结果的默认评判。

在本发明的一个实施例中，执行所述数据处理程序还包括实现，在所述消除比对结果中的每个读段的缺口之前，提取所述比对结果中的唯一比对结果以替换所述比对结果，所述唯一比对结果包括唯一比对上所述参考序列的多个读段对，并且每一所述读段对比对到所述参考序列的相同染色体，每一所述读段对中的两个读段的距离符合预期的其来自的所述染色体片段的两个位置之间的距离。

在本发明的一个实施例中，执行所述数据处理程序还包括实现，修正使所述唯一比对结果中的每一对读段对比对到所述参考序列的相同染色体的正链。例如，对于分别比对上一染色体的正负链的一对读段，将比对上负链的读段变成其反向互补链，这样以其反向互补链来替代该读段得以实现所说的修正。

在本发明的一个实施例中，执行所述数据处理程序还包括实现数据格式转换，所述数据格式转换包括转换所述比对结果或所述唯一比对结果的格式。在本发明得一个实施利中，要求通用比对结果的格式为SAM或BAM，以利于后续的基于比对结果或比对结果对数据进一步分析，SAM或BAM是常见的二进制格式，BAM是SAM的压缩格式。由于利用不同的比对软件，可能输出的比对结果或者唯一比对结果的格式不适用于现有的后续数据处理或者分析软件程序，例如前述的TeraMap格式的比对结果，其输出的数据格式不满足现有的大部分变异检测软件SOAPsnp、GATK或SOAPindel等对输入数据格式的要求，转换数据格式使获得具有通用数据格式的通用比对结果，便于对比对数据进一步分析处理。

依据本发明的第二方面，提供一种测序数据处理系统，其包括一主机和一显示器，该系统还包括本发明一方面或者任一具体实施方式中的测序数据处理装置。前述对测序数据处理装置的优点及技术特征的描述，同样适用本发明的这一系统，在此不再赘述。

依据本发明的第三方面，提供一种测序数据处理方法，该方法包括如下步骤：获取测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读段对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体片段的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂；将所述测序数据与参考序列比对，获得比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果；消除所述比对结果中的每个读段的缺口，获得通用比对结果。关于读段对的获取方式、读段包含的缺口、比对、消除缺口，比对结果和通用比对结果等的特点可以参照上述对本发明一方面或者任一具体实施方式中的装置中的相应技术特征的描述。例如，同样的，这里所说的分别来自一条染色体片段的两个位置的读段对，可以通过构建末端文库(pair-end library)或者配对末端文库(mate-pair library)，对所构建的文库进行测序来获得，在本发明的一个实施例中，利用Complete Genomics公司(CG)的文库构建方法及其测序平台，获得多对读段对，一对读段对之间的距离是由读段的长度以及酶的识别位点和切割位点的距离来控制的。CG平台通过酶切构建多接头配对末端文库，利用特有的组合探针连接测序(cPAL)技术对所构建的环状文库进行测序，测读出接头两旁的碱基，因为其是利用酶切连接一个接头的两段来进行配对末端文库构建的，由于每一种酶都有一个首选的切割距离，而在实际酶切时经常比首选距离多一个位置或少一个位置，这样使得读段中经常带有缺口(gap)，缺口常为+1或者-1，和/或，建库时倘若使用同一种酶多次酶切，每次的酶切位置易发生变化，酶切位置的变化也会使获得的读段带有缺口，例如在构建多接头环状文库时，利用Alu酶两次酶切来连接多个接头的不同部分，读测这些接头旁的碱基时，会产生带+3/-3的缺口的读段。在本发明中缺口的大小还可以是0。来自多接头文库中的多个读段，一个读段可以和任一其它读段组成一对读段对。这里所说的“正链”和“负链”是组成染色体片段的互补的两条链，是相对的，称一条链为正链就可以称其互补链为负链，在本发明的一个实施例中，将与参考序列匹配的链称为正链，将另一条链称为负链。在本发明中，比对可以利用已知比对软件进行，比如SOAP、BWA等，也可以利用CG平台的比对软件TeraMap进行。在本发明的一个实施例中，比对是利用TeraMap进行的，所得的比对结果的格式为TeraMap。在本发明的一个实施例中，所说的消除比对结果中每个读段的缺口是指，对带负缺口的读段去除掉其负缺口即去除掉重叠的碱基，对带正缺口的读段以N替代正缺口的大小，N为A、T、C或G，例如，对带负缺口比如为-2nt的读段，基于缺口处该读段可分成两部分，两部分的末端有2nt重叠，比如该读段的两部分分别为ATCGCTTAAG和AGTACGATTC，消除其负缺口即重叠的AG，获得对应的读段为ATCGCTTAAGTACGATTC。

在本发明的一个实施例中，获取所述测序数据包括构建测序文库，获得测序文库，所述测序文库为单链环状DNA文库，所述测序文库由所述染色体片段的一条链和至少一个预定DNA序列构成。所说的单链环状文库可以利用已知文库构建方法来构建，比如参考LifeTechnologies公司的SOLiD的配对末端文库的构建获得单接头环状双链文库，接着分离双链获得单链环状文库，在本发明的一个实施例中，单链环状文库是利用CG的建库技术来构建，文库构建可参考US7897344，获得多接头单链环状文库。

在本发明的一个实施例中，所述每对读段分别来自所述染色体片段的两端。通过参考改进CG的建库技术，利用一个接头的两部分分别连接于一段染色体片段的两端，单链分离，单链成环，获得1-接头单链环状文库，该1-接头单链环状文库由所述染色体片段的一条链和连接所述一条链的两端的一个预定DNA序列构成，滚环扩增形成DNA纳米球(DNB)，利用CG测序cPAL技术对DNB进行测序，关于DNB种植在芯片上以及cPAL技术可分别参考US8278039B2和US8518640B2。所说的预定DNA序列为已知序列，为前述的接头或者接头的一条链。所说的改进CG建库方法构建1-接头环状单链文库包括步骤：(1)提取待测核酸；(2)末端磷酸化所述核酸，获得末端磷酸化产物；(3)末端修复所述末端磷酸化产物，获得末端修复产物；(4)将第一序列和第二序列连接至所述末端修复产物的两端，获得第一连接产物；(5)利用第三序列对所述连接产物进行缺刻平移和扩增，获得扩增产物，所述第三序列为一对引物对，所述引物对的至少一条引物带有生物素标记；(6)利用所述生物素标记对所述扩增产物进行单链分离，获得单链产物；(7)利用第四序列环化所述单链产物，获得所述测序文库；其中，所述第四序列能够连接所述第一序列的一端和所述第二序列的一端，所述第一序列和/或所述第二序列的另一端为双脱氧核苷酸。所说的第四序列能够连接所说的第一序列和所说的第二序列形成一个所说的接头，缺刻平移是为消除连接在末端修复产物两端的第一序列和/或第二序列的另一端的双脱氧核苷酸造成的缺刻(nick)，利用至少一条引物带有生物素标记使扩增产物的至少一条链带有生物素标记，使后续易于基于该生物素标记分离获得单链产物。在本发明的一个实施例中，所说的改进CG建库方法构建1-接头环状单链文库包括步骤：(1)提取待测核酸；(2)末端修复所述核酸，获得末端修复产物；(3)末端磷酸化所述末端修复产物，获得末端磷酸化产物；(4)将第一序列和第二序列连接至所述末端磷酸化产物的两端，获得第一连接产物；(5)利用第三序列对所述连接产物进行缺刻平移和扩增，获得扩增产物，所述第三序列为一对引物对，所述引物对的至少一条引物带有生物素标记；(6)利用所述生物素标记对所述扩增产物进行单链分离，获得单链产物；(7)利用第四序列环化所述单链产物，获得所述测序文库；其中，所述第四序列能够连接所述第一序列的一端和所述第二序列的一端，所述第一序列和/或所述第二序列的另一端为双脱氧核苷酸。末端修复和末端磷酸化哪个步骤先进行不作限制。末端修复是为获得平末端核酸片段，使得能够连接其它核苷酸或序列。末端磷酸化是为了减少样本核酸片段的互相连接，使得核酸含量很低的样本也能够进行文库构建且满足文库上机量要求。单接头环状单链文库如图1所示，将所构建的单接头环状单链文库(1-AD)上机测序，1-AD测序输出读段对总长约为30bp，其中一读段12bp，一读段19bp，一读段对两个读段之间在基因组上的距离的中位数约为140bp。单接头建库量较小，适合cfDNA含量较少的情况，此外还有建库时间短及建库成本低的优点。

在本发明的一个实施例中，本发明的这一方法中的比对包括：将每对读段对的左臂和右臂分别与所述参考序列比对，获得一级左比对结果和一级右比对结果，分别以所述一级左比对结果和所述一级右比对结果的其中一个为参考，对另一个作比对，获得二级左比对结果和二级右比对结果，基于所述二级左比对结果和所述二级右比对结果，获得多个所述读段对的比对结果，或者获得多个所述左臂的比对结果和多个所述右臂的比对结果。这样经过两次比对，可以获得读段对比对结果，在本发明的一个实施例中，第一次比对是与参考序列作全局比对，以该次左臂/右臂比对结果为基准对右臂/左臂比对结果进行的第二次比对为局部比对，这样，能够将分别来自二级左比对结果和二级右比对结果中的比对到同一染色体上且之间的距离符合预期的两个读段配对成一对读段对，获得读段对比对结果。

在本发明的一个实施例中，所说的比对包括，设置所述缺口的大小以使每个左臂或者每个右臂与所述参考序列进行多次比对，以获得最佳比对结果。例如，将所述每个左臂或者所述每个右臂的缺口分别设置为-3nt、-2nt、-1nt、0nt、1nt、2nt、3nt、4nt、5nt、6nt和7nt，获得对应的多个读段，分别将所述对应的多个读段与所述参考序列比对，将最优比对的序列作为该左臂/右臂，这里对于比对结果的好坏可以基于所利用的比对软件对比对结果的默认评判。

在本发明的一个实施例中，执行所述数据处理程序还包括实现，在所述消除比对结果中的每个读段的缺口之前，提取所述比对结果中的唯一比对结果以替换所述比对结果，所述唯一比对结果包括唯一比对上所述参考序列的多个读段对，并且每一所述读段对比对到所述参考序列的相同染色体，每一所述读段对中的两个读段的距离符合预期的其来自的所述染色体片段的两个位置之间的距离。

在本发明的一个实施例中，执行所述数据处理程序还包括实现，修正使所述唯一比对结果中的每一对读段对比对到所述参考序列的相同染色体的正链。例如，对于分别比对上一染色体的正负链的一对读段，将比对上负链的读段变成其反向互补链，这样以其反向互补链来替代该读段得以实现所说的修正。

在本发明的一个实施例中，执行所述数据处理程序还包括实现数据格式转换，所述数据格式转换包括转换所述比对结果或所述唯一比对结果的格式。在本发明得一个实施利中，要求通用比对结果的格式为SAM或BAM，以利于后续的基于比对结果或比对结果对数据进一步分析，SAM或BAM是常见的二进制格式，BAM是SAM的压缩格式。由于利用不同的比对软件，可能输出的比对结果或者唯一比对结果的格式不适用于现有的后续数据处理或者分析软件程序，例如前述的TeraMap格式的比对结果，其输出的数据格式不满足现有的大部分变异检测软件SOAPsnp、GATK或SOAPindel等对输入数据格式的要求，转换数据格式使获得具有通用数据格式的通用比对结果，便于对比对数据进一步分析处理。

依据本发明的第四方面，本发明提供一种计算机可读存储介质，其用于存储供计算机执行的程序，所述程序的执行包括完成前述本发明一方面的或者其任一具体实施方式中的测序数据处理方法。前述对本发明的测序数据处理方法的优点和技术特征的描述也适用于该计算机可读存储介质，在此不再赘述。所称存储介质可以包括：只读存储器、随机存储器、磁盘或光盘等。

依据本发明的第五方面，本发明提供一种检测拷贝数变异(CNV)的方法，该方法包括：a.获取待测样本的核酸；b.对所述核酸进行测序，获得测序数据；c.对所述测序数据进行处理，以获得通用比对结果；d.基于所述通用比对结果检测所述CNV；其中，c步骤本发明一方面的或者任一具体实施方式中的测序数据处理装置和/或方法进行的。上述对本发明的测序数据处理装置和/或方法的优点及技术特征的描述，也适用于本发明这一方面的CNV检测方法，在此不再赘述。

在本发明的一个实施例中，b步骤包括，对所述核酸进行测序文库构建，获得测序文库，所述测序文库为单链环状DNA文库，所述单链环状DNA文库的构建包括：末端磷酸化所述核酸，获得末端磷酸化产物；末端修复所述末端磷酸化产物，获得末端修复产物；将第一序列和第二序列连接至所述末端修复产物的两端，获得第一连接产物；利用第三序列对所述连接产物进行缺刻平移和扩增，获得扩增产物，所述第三序列为一对引物对，所述引物对的至少一条引物带有生物素标记；利用所述生物素标记对所述扩增产物进行单链分离，获得单链产物；环化所述单链产物，获得所述测序文库；其中，所述第四序列能够连接所述第一序列的一端和所述第二序列的一端，所述第一序列和/或所述第二序列的另一端为双脱氧核苷酸。在本发明的另一实施例中，先进行末端修复再进行末端磷酸化。末端修复是为获得平末端核酸片段，使得能够连接其它核苷酸或序列。末端磷酸化是为了减少样本核酸片段的互相连接，使得核酸含量很低的样本也能够进行文库构建且满足文库上机量要求。单接头环状单链文库如图1所示，单接头建库量较小，适合cfDNA含量较少的情况，此外还有建库时间短及建库成本低的优点。所说的第四序列能够连接第一序列和的第二序列形成一个所说的接头，缺刻平移是为消除连接在末端修复产物两端的第一序列和/或第二序列的另一端的双脱氧核苷酸造成的缺刻(nick)，利用至少一条引物带有生物素标记使扩增产物的至少一条链带有生物素标记，使后续易于基于该生物素标记分离获得单链产物。在本发明的一个实施例中，对构建的文库进行测序是利用组合探针锚定连接测序技术进行的，例如利用CG测序平台进行。

基于通用比对结果检测CNV可以利用目前已知的CNV检测方法，比如利用隐马可夫模型、环状二元分割、等级分割或核平滑算法等。在本发明的一个实施例中，d步骤包括：在所述参考序列上设置多个窗口，基于所述通用比对结果中匹配到所述窗口的读段的量与对照样本的通用比对结果中匹配到相同窗口的读段的量的差异具有显著性，判定所述待测样本核酸存在所述CNV，其中，所述窗口为所述参考序列的一部分。其中，窗口的大小可依据预检测的CNV的大小来调整设置，对照样本的通用比对结果可以通过本发明的一方面的方法或者其任一具体实施方式中的测序数据处理方法获得，差异是否为显著性的判断可以利用统计检验比如z-score(标准分数)计算z值来进行，当z值大于或小于某一预定阈值时则判定该窗口区域存在CNV，比如正常对照为二倍体(CNV＝2)，当z值为正数时说明待测样本的该窗口的CNV>2，为负数则说明待测样本的该窗口的CNV<2，在本发明的一个实施例中，设定预定阈值为3，即当z值得绝对值大于3时则确定该窗口发生CNV。所说的读段的量可以为一个数目，也可以为一个比值，例如，也可以基于待测样本的窗口的测序深度与对照样本相应窗口的测序深度的差异，使用z-score(标准分数)进行检验检测出拷贝数变异，所说的窗口的测序深度＝比对到该窗口的读段的量/该窗口的大小。在本发明的一个实施例中，考虑到实际测序过程中读段(reads)中GC含量会对测序深度有一定的影响[Alkan,Can,Jeffrey MKidd,Tomas Marques-Bonet,Gozde Aksay,Francesca Antonacci,FereydounHormozdiari,Jacob O Kitzman,et al.“Personalized Copy Number and SegmentalDuplication Maps Using next-Generation Sequencing.”Nature Genetics 41,no.10(October 2009):1061–67]，先进行GC含量校正，消除GC含量对测序深度的影响。所说的GC含量校正，可以利用多个对照样本的测序数据，取多个窗口算窗口的GC含量和平均测序深度，对GC-测序深度的数据进行二维回归分析，例如利用局部加权回归散点平滑法(lowess回归)建立二者的关系，根据回归得的关系对各个窗口的测序深度进行GC含量校正。所说的测序深度和GC含量的关系可以通过以下来建立：获得多个对照样本核酸的测序数据，所述测序数据由多个读段组成；在所述参考序列上设置多个窗口，将所述多个对照样本的测序数据分别与所述参考序列的窗口比对，计算各个对照样本的测序数据中比对上每个窗口的读段的数目，获得每个窗口的测序深度，所述窗口为所述参考序列的一部分，所述窗口的测序深度＝各个对照样本的比对上所述窗口的读段的总数目/(对照样本个数*所述窗口的大小)；基于每个窗口的测序深度和该窗口的GC含量，利用二维回归分析法建立所述测序深度和GC含量的关系。

在本发明的一个实施例中，d步骤包括：在所述参考序列上设置多个窗口，计算窗口的测序深度，窗口的测序深度＝所述通用比对结果中比对到所述窗口的读段的数量/所述窗口的大小；利用测序深度和GC含量的关系校正所述窗口的测序深度，获得窗口的校正测序深度；基于所述窗口的校正测序深度与对照样本的相同窗口的校正测序深度的差异具有显著性，判定所述待测样本核酸存在所述CNV，其中，所述窗口为所述参考序列的一部分。较佳地，前述的对照样本的个数不小于30个，样本数目达到30个使样本数据呈现满足特定分布符合适于利用多数统计检验方法来检验，例如，利用t检验、z检验等来统计检验一般要求多个样本数据符合正态分布。所说的对照样本的相同窗口的校正测序深度是利用所述测序深度和GC含量的关系校正对照样本的相同窗口的测序深度获得的，所述对照样本的相同窗口的测序深度＝所述对照样本的测序数据中比对到所述窗口的读段的数目/所述窗口的大小。前述对照样本的测序数据、比对结果等可以通过参照前述本发明一方面的或者任一具体实施方式中的测序数据处理方法获得，可以与待测样本的测序数据、比对结果同时获得，也可以预先获得保存备用。

依据本发明的第六方面，本发明提供一种CNV检测设备，该设备用以执行本发明一方面的CNV检测方法的全部或部分步骤，所述设备包括：核酸获取装置，用以获取待测样本的核酸；测序装置，用以对来自所述核酸获取单元的核酸进行测序，获得测序数据；数据处理装置，用于对来自所述测序装置的测序数据进行处理，以获得通用比对结果；检测装置，用于基于来自所述数据处理装置的通用比对结果检测所述CNV；其中，所述数据处理装置包括，数据接收单元，用于接收来自所述测序装置的测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读长对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂，处理器，用于执行数据处理程序，执行所述数据处理程序包括实现将所述测序数据与参考序列作比对，获得比对结果，以及消除所述比对结果中的每个读段的缺口，获得通用比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果，以及，至少一个存储单元，用于存储数据，其中包括所述数据处理程序。前述对本发明一方面的或者任一其具体实施方式中的CNV检测方法的优点和技术特征的描述，同样适用本发明这一方面的CNV检测设备，在此不再赘述，而且，本领域普通技术人员可以理解，本发明的这一装置中的全部或部分单元，可选择的、可拆卸的包含一个或多个子单元以执行或实现前述本发明CNV检测方法的各个具体实施方式。

通过CG平台单接头测序获得测序数据，成本更低速度也更快。利用本发明的数据处理装置、系统和/或方法，开发TeraMap2Sam转换软件，将CG平台TeraMap的比对结果转化为通用的SAM格式，使能够直接使用Samtools，GATK等众多优秀的开源软件进行变异检测，使后续分析的选择更加广泛。利用本发明的CNV检测方法和/或设备开发的CNV检测程序基于标准分数方法做CNV分析，速度快，分辨率高。

附图说明

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

图1是本发明的一个实施例中的单接头环状单链文库的结构示意图；

图2是本发明的一个实施例中的测序数据处理装置的结构示意图；

图3是本发明的一个实施例中的测序数据处理系统的结构示意图；

图4是本发明的一个实施例中的测序数据处理方法的流程图；

图5是本发明的一个实施例中的测序数据处理方法的流程图；

图6是本发明的一个实施例中的CNV检测方法的流程图；

图7是本发明的一个实施例中的CNV检测方法的流程图；

图8是本发明的一个实施例中的CNV检测设备的结构示意图；

图9是本发明的一个实施例中的单接头文库构建和测序得流程图；

图10是本发明的一个实施例中的Teramap2Sam软件的算法流程图。

具体实施方式

下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中，自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的，仅用于解释本发明，而不能理解为对本发明的限制。需要说明的是在本文中所使用的术语“第一”、“第二”、“第三”、“第四”或者“一级”、“二级”等仅为方便描述指代，而不能理解为指示或暗示相对重要性，也不能理解为之间有先后顺序关系。在本发明的描述中，除非另有说明，“多个”的含义是两个或两个以上。

图2显示本发明的测序数据处理装置的一个实施例的装置结构示意图，该测序数据处理装置100包括：数据接收单元10、处理器20和一个存储单元30，处理器20与数据接收单元10和存储单元30连接，存储单元30和数据接收单元10连接。数据接收单元10，用于接收测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读长对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读段对中的两个读段都来自所述染色体片段的正链或所述染色体的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂。这里所说的分别来自一条染色体片段的两个位置的读段对，可以通过构建末端文库(pair-end library)或者配对末端文库(mate-pair library)，对所构建的文库进行测序来获得，在本发明的一个实施例中，利用Complete Genomics公司(CG)的文库构建方法及其测序平台，获得多对读段对，一对读段对之间的距离是由读段的长度以及酶的识别位点和切割位点的距离来控制的。CG平台通过酶切构建多接头配对末端文库，利用特有的组合探针连接测序(cPAL)技术对所构建的环状文库进行测序，测读出接头两旁的碱基，因为其是利用酶切连接一个接头的两段来进行配对末端文库构建的，由于每一种酶都有一个首选的切割距离，而在实际酶切时经常比首选距离多一个位置或少一个位置，这样使得读段中经常带有缺口(gap)，缺口常为+1或者-1，和/或，建库时倘若使用同一种酶多次酶切，每次的酶切位置易发生变化，酶切位置的变化也会使获得的读段带有缺口，例如在构建多接头环状文库时，利用Alu酶两次酶切来连接多个接头的不同部分，读测这些接头旁的碱基时，会产生带+3/-3的缺口的读段。在本发明中缺口的大小还可以是0。以CG平台目前的双接头(two adaptors，2-AD)测序文库为例，2-AD测序输出总长为60bp，可分为两对读段对(mate-paired reads)，每对读段对中的读段在10bp的位置都有小的gap，在20bp位置有一个无效测序位点N，一对读段对的两个读段之间的在基因组上的距离一般小于2000bp。来自多接头文库中的多个读段，一个读段可以和任一其它读段组成一对读段对。这里所说的“正链”和“负链”是组成染色体片段的互补的两条链，是相对的，称一条链为正链就可以称其互补链为负链，在本发明的一个实施例中，将与参考序列匹配的链称为正链，将另一条链称为负链。

处理器20，用于执行数据处理程序，执行所述数据处理程序包括实现将所述测序数据与参考序列作比对，获得比对结果，以及消除所述比对结果中的每个读段的缺口，获得通用比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果。比对可以利用已知比对软件进行，比如SOAP、BWA等，也可以利用CG平台的比对软件TeraMap进行。在本发明的一个实施例中，比对是利用TeraMap进行的，所得的比对结果的格式为TeraMap。在本发明的一个实施例中，所说的消除比对结果中每个读段的缺口是指，对带负缺口的读段去除掉其负缺口即去除掉重叠的碱基，对带正缺口的读段以N替代正缺口的大小，N为A、T、C或G，对缺口为0的读段不作处理，例如，对带负缺口比如为-2nt的读段，基于缺口处该读段可分成两部分，两部分的末端有2nt重叠，比如该读段的两部分分别为ATCGCTTAAG和AGTACGATTC，消除其负缺口即重叠的AG，获得对应的读段为ATCGCTTAAGTACGATTC。

存储单元30，用于存储数据，存储单元30中存储有上述数据处理程序，也存储来自数据接收单元10的测序数据和处理器20的处理的中间数据或结果。

图3显示本发明的测序数据处理系统的一个实施例中的系统结构示意图。该测序数据处理系统1000包括：测序数据处理装置100，主机200和显示装置300。主机200可为一音频/视频/信号发送源装置，比如电脑主机、大型机等，其用以传送显示装置300所需的显示数据。该主机200包含至少一个与测序数据处理装置100电性连接的接口，测序数据处理装置100接收从主机200输出的测序数据，并对测序数据进行处理，然后将处理的数据或结果输出到显示装置300。

图4显示本发明的测序数据处理方法的一个实施例的测序数据处理流程图。该测序数据处理方法包括步骤：S1获取测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读段对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体片段的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂；S2将所述测序数据与参考序列比对，获得比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果；S3消除所述比对结果中的每个读段的缺口，获得通用比对结果。关于读段对的获取方式、读段包含的缺口、比对、消除缺口，比对结果和通用比对结果等的特点可以参照上述对本发明一方面或者任一具体实施方式中的测序数据处理装置中的相应技术特征的描述。例如，同样的，这里所说的分别来自一条染色体片段的两个位置的读段对，可以通过构建末端文库(pair-end library)或者配对末端文库(mate-pair library)，对所构建的文库进行测序来获得，在本发明的一个实施例中，利用Complete Genomics公司(CG)的文库构建方法及其测序平台，获得多对读段对，一对读段对之间的距离是由读段的长度以及酶的识别位点和切割位点的距离来控制的。CG平台通过酶切构建多接头配对末端文库，利用特有的组合探针连接测序(cPAL)技术对所构建的环状文库进行测序，测读出接头两旁的碱基，因为其是利用酶切连接一个接头的两段来进行配对末端文库构建的，由于每一种酶都有一个首选的切割距离，而在实际酶切时经常比首选距离多一个位置或少一个位置，这样使得读段中经常带有缺口(gap)，缺口常为+1或者-1，和/或，建库时倘若使用同一种酶多次酶切，每次的酶切位置易发生变化，酶切位置的变化也会使获得的读段带有缺口，例如在构建多接头环状文库时，利用Alu酶两次酶切来连接多个接头的不同部分，读测这些接头旁的碱基时，会产生带+3/-3的缺口的读段。在本发明中缺口的大小还可以是0。来自多接头文库中的多个读段，一个读段可以和任一其它读段组成一对读段对。这里所说的“正链”和“负链”是组成染色体片段的互补的两条链，是相对的，称一条链为正链就可以称其互补链为负链。这里，将与参考序列匹配的链称为正链，将另一条链称为负链。比对可以利用已知比对软件进行，比如SOAP、BWA等，也可以利用CG平台的比对软件TeraMap进行。在本发明的一个实施例中，比对是利用TeraMap进行的，所得的比对结果的格式为TeraMap。在本发明的一个实施例中，所说的消除比对结果中每个读段的缺口是指，对带负缺口的读段去除掉其负缺口即去除掉重叠的碱基，对带正缺口的读段以N替代正缺口的大小，N为A、T、C或G，对缺口为0的读段不作处理，例如，对带负缺口比如为-2nt的读段，基于缺口处该读段可分成两部分，两部分的末端有2nt重叠，比如该读段的两部分分别为ATCGCTTAAG和AGTACGATTC，消除其负缺口即重叠的AG，获得对应的读段为ATCGCTTAAGTACGATTC。

图5是本发明的测序数据处理方法的一个实施例的数据处理流程图。该测序数据处理方法包括：S10获取测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读段对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体片段的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂；S20将所述测序数据与参考序列比对，获得比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果；S30提取所述比对结果中的唯一比对结果以替换所述比对结果，所述唯一比对结果包括唯一比对上所述参考序列的多个读段对，并且每一所述读段对比对到所述参考序列的相同染色体，每一所述读段对中的两个读段的距离符合预期的其来自的所述染色体片段的两个位置之间的距离；S40修正使所述唯一比对结果中的每一对读段对比对到所述参考序列的相同染色体的正链。例如，对于分别比对上一染色体的正负链的一对读段，将比对上负链的读段变成其反向互补链，这样以其反向互补链来替代该读段得以实现所说的修正；S50消除所述唯一比对结果中的每个读段的缺口，获得通用比对结果。

图6是本发明的CNV检测方法的一个实施例的检测流程图。该CNV检测方法包括步骤：S11获取待测样本的核酸；S12对所述核酸进行测序，获得测序数据；S13对所述测序数据进行处理，以获得通用比对结果；S14基于所述通用比对结果检测所述CNV；其中，S13是利用本发明一方面的或者任一具体实施方式中的测序数据处理装置和/或测序数据处理方法进行的。基于通用比对结果检测CNV可以利用目前已知的CNV检测方法，比如利用隐马可夫模型、环状二元分割、等级分割或核平滑算法等。

图7是本发明的CNV检测方法的一个实施例的检测流程图。该CNV检测方法包括步骤：S110获取待测样本的核酸；S120对所述核酸进行测序，获得测序数据；S130对所述测序数据进行处理，以获得通用比对结果，S130是通过上述本发明一方面的或者任一具体实施方式中的测序数据处理装置和/或测序数据处理方法进行的；S140在所述参考序列上设置多个窗口，计算窗口的测序深度，窗口的测序深度＝所述通用比对结果中比对到所述窗口的读段的数量/所述窗口的大小；S150利用测序深度和GC含量的关系校正所述窗口的测序深度，获得窗口的校正测序深度；S160基于所述窗口的校正测序深度与对照样本的相同窗口的校正测序深度具有显著差异，判定所述待测样本核酸存在所述CNV，其中，所述窗口为所述参考序列的一部分。前述的对照样本的个数不小于30个，样本数目达到30个使样本数据呈现满足特定分布符合适于利用多数统计检验方法来检验，例如，利用t检验、z检验等来统计检验一般要求多个样本数据符合正态分布。所说的对照样本的相同窗口的校正测序深度是利用所述测序深度和GC含量的关系校正对照样本的相同窗口的测序深度获得的，所述对照样本的相同窗口的测序深度＝所述对照样本的测序数据中比对到所述窗口的读段的数目/所述窗口的大小。前述对照样本的测序数据、比对结果等可以通过参照前述本发明一方面的或者任一具体实施方式中的测序数据处理方法获得，可以与待测样本的测序数据、比对结果同时获得，也可以预先获得保存备用。所说的测序深度和GC含量的关系可以通过以下来建立：获得多个对照样本核酸的测序数据，所述测序数据由多个读段组成；在所述参考序列上设置多个窗口，将所述多个对照样本的测序数据分别与所述参考序列的窗口比对，计算各个对照样本的测序数据中比对上每个窗口的读段的数目，获得每个窗口的测序深度，所述窗口为所述参考序列的一部分，所述窗口的测序深度＝各个对照样本的比对上所述窗口的读段的总数目/(对照样本个数*所述窗口的大小)；基于每个窗口的测序深度和该窗口的GC含量，利用二维回归分析法，例如利用Lowess回归建立所述测序深度和GC含量的关系。

图8是本发明的CNV检测设备的一个实施例的设备结构示意图。该设备2000包括：核酸获取装置200，用以获取待测样本的核酸；测序装置400，用以对来自所述核酸获取单元的核酸进行测序，获得测序数据；数据处理装置600，用于对来自所述测序装置的测序数据进行处理，以获得通用比对结果；检测装置800，用于基于来自所述数据处理装置600的通用比对结果检测所述CNV；其中，所述数据处理装置600包括，数据接收单元610，用于接收来自所述测序装置的测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读长对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂，处理器630，用于执行数据处理程序，执行所述数据处理程序包括实现将所述测序数据与参考序列作比对，获得比对结果，以及消除所述比对结果中的每个读段的缺口，获得通用比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果，以及，至少一个存储单元650，用于存储数据，其中包括所述数据处理程序。前述对本发明一方面的或者任一其具体实施方式中的CNV检测方法的优点和技术特征的描述，同样适用本发明这一方面的CNV检测设备，在此不再赘述，而且，本领域普通技术人员可以理解，本发明的这一装置中的全部或部分单元，可选择的、可拆卸的包含一个或多个子单元以执行或实现前述本发明CNV检测方法的各个具体实施方式。

以下实施例仅用于说明本发明的优选实施方式，实施例中未注明具体操作手段或条件的，可以按照本领域内的文献所描述的技术或条件(例如参考J.萨姆布鲁克等著，黄培堂等译的《分子克隆实验指南》，第三版，科学出版社)或者按照产品说明书进行。所用试剂或仪器未注明生产厂商者，均为可以通过市购获得的常规产品或服务。

实施例一

以下以肺癌患者的外周血血浆作为检测对象，样本来自西南医院，进行如下检测：

(一)文库建立及测序

建库及测序流程如图9所示，以下涉及的具体序列，都是从左到右为5’端至3’端，序列中的“//”中为末端修饰基团，“phos”表示磷酸化，“dd”表示双脱氧，“bio”表示生物素。

1、cfDNA的提取(采用SnoMag Circulating DNA Kit)：

1)取200ul血浆于1.5mlEP管，加入600ul buffer LSB。

2)加入20μlNanoMag Circulating Beads混匀，室温放置10min，每2-3min混匀一次。

3)将EP管置于磁力架上吸附1min，弃上清。

4)取下EP管加入150uL Buffer WA，混匀。

5)将EP管置于磁力架上吸附1min，弃上清。

6)取下EP管加入150uL 75％乙醇，混匀。

7)将EP管置于磁力架上吸附1min，弃上清。

8)重复6-7一次。

9)室温干燥磁珠5min。

10)加入32ul elution buffer混匀磁珠，室温静置5min。

11)将EP管置于磁力架上吸附1min，转移上清至新的1.5mlEP管。

2、文库的构建：

1)rSAP去磷酸化反应

反应条件：

2)T4 DNA Polymerase末端补平

反应条件：

12℃	20min
		4℃	hold

60ulAmpureXP beads纯化以上反应产物，22ulElution buffer洗脱。

3)第一序列和第二序列分别连接到末端补平的DNA片段的两端

反应条件：

20℃	15min
		4℃	hold

40ulAmpureXP beads纯化以上反应产物,22ulElution buffer洗脱。

第一序列的两条链为：TTGGCCTCCGACT/3-ddT/(SEQ ID NO：1)，/5phos/AAGTCGGAGGCCAAGCGGTCGT/ddC/(SEQ ID NO：2)。

第二序列的两条链分别为：/5Phos/GTCTCCAGTCGAAGCCCGACG/3ddC/(SEQ ID NO：3)，GCTTCGACTGGAGA/3ddC/(SEQ ID NO：4)。

4)缺刻平移(Nick Translation)

第三序列中的上游引物/5-bio/TCCTAAGACCGCTTGGCCTCCGACT(SEQ ID NO：5)，

第三序列中的下游引物5Phos/AGACAAGCTCxxxxxxxxxxGATCGGGCTTCGACTGGAGAC(SEQ ID NO：6)，中间“x”处为可变的标签序列区域，可以以N替代，N为A、T、C或G，当没有其它样本文库一起混合上机，只有一个样本文库上机，不需要标签序列，即第三序列可为5Phos/AGACAAGCTCGATCGGGCTTCGACTGGAGAC(SEQ ID NO：7)，在该示例中，由于是肿瘤游离核酸样本，混合核酸中的目标核酸(ctDNA)含量低，若多个这样的样本文库混合上机获得混合数据，需要拆分混合数据对应到各自样本，会损失一部分数据，且构建的是单接头环状文库读段相对短，要准确检测变异需要深度测序获得相对大量的测定数据，较佳的，单个样本文库上机。

反应条件：

60℃	5min
		37℃	0.1℃/secs-hold

向上步反应物中加入如下8ul Nick Translation mix

反应条件：

37℃	20min
		4℃	hold

40ulAmpureXP beads纯化以上反应产物,37.4ulElution buffer洗脱。

5)PCR with Pfx

反应条件：

50ul AmpureXP beads纯化以上反应产物，22ulElution buffer洗脱。

6)Qubit定量

利用Qubit dsDNA HS assay kit对PCR产物进行浓度测定。

7)链分离(Strand Separation)

a)多个文库混合，使DNA共约160ng。样品补1xTE至总体积为60ul。

b)提前准备以下试剂：4X BBB，Streptavidin Beads，0.3M MOPS acid，0.5％Tween20，1X BBB/Tween Mix，1X BWB/Tween Mix，0.1M NaOH。其中1X BWB/Tween Mix、0.1MNaOH，Streptavidin Beads需现配现用。

c)提前15min配置以下四种试剂

0.5％Tween20，1X BBB/Tween Mix，1X BWB/Tween Mix，0.1M NaOH.

其中0.5％Tween20配置方法同前述，其他三种配置方法如下：

d)1X BBB/Tween Mix

1X BBB	30ul
		0.5％Tween20	0.3ul
Total	30.3ul

e)1X BWB/Tween Mix

1X BWB	2000ul
		0.5％Tween20	20ul
Total	2020ul

f)0.1M NaOH

0.5M NaOH	15.6ul
		Water	62.40ul
Total	78.0ul

g)Streptavidin Beads洗涤方法如下：

·每个样品取30ul Streptavidin Beads：加入3-5倍体积的1XBBB，混匀后置于磁力架上静止吸附，调整不粘管的方向，使得beads在1XBBB洗液中前后游动，弃上清液后，重复上述操作一次，

·取出不黏管加入1倍体积(30ul)1X BBB/Tween Mix悬浮，混匀后室温静置。

h)向60ulPCR产物样品中加入20ul 4XBBB混匀，然后转移到上步骤含有30ul1XBBB/Tween Mix溶解的beads的不粘管中混匀，此110ul混合物室温下结合15-20min，中间轻轻弹匀一次。

i)将上述不粘管磁力架放置3-5min，弃去上清液，用1ml的1X BWB/Tween Mix洗涤2次，方法同Streptavidin Beads的洗涤方法

j)向上述beads中加入26ul 0.1M NaOH，吹打混匀后放置10min，再置于磁力架上3-5min，取上清到新的1.5ml EP管中。

k)向上述1.5mlEP管中加入13ul 0.3M MOPS，混匀备用。

l)此步骤产物可以冻存于-20℃。

8)环化(Splint Circulation)

a)向上一步得到的39ul的样品中加入10ul的20uM第四序列，第四序列为TCGAGCTTGTCTTCCTAAGACCGC(SEQ ID NO：8)；

b)提前5分钟准备连接酶反应混合液，配制如下：

Water	4.2ul
		10x TA Buffer(LK1)	6ul
100mM ATP	0.6ul
		600U/ul Ligase	0.2ul
Total	11ul

c)将连接酶反应混合液震荡充分混匀，离心后，向已经加入引物反应混合液的EP管中加入连接酶反应混合液11ul，震荡10s混匀，瞬时离心。

d)置于PCR仪中37℃孵育1.5h。

e)反应完成后，取出5ul样品，待6％变性胶电泳检测，剩余的约55ul体积，进入下一步酶反应。

9)酶切消化(Exo I and III)

a)提前5分钟左右准备引物反应混合液，配制如下：

10x TA Buffer(LK1)	1ul
		20U/ul Exo I	3ul
200/ul Exo III	1ul
		Total	5ul

b)将上述混合液震荡充分混匀，离心后，向上一步得到的55ul的样品中分别加入5ul的反应混合液；

c)震荡10s混匀离心，置于PCR仪中37℃孵育30min。

d)酶切30min完成后，向样品中加入2.5ul 500mM EDTA终止酶反应。

e)上述样品用PEG32 beads/tween20纯化，方法如下：

将上步骤样品59ul转移到1.5ml不粘管中，加入78ul的PEG32 beads/tween20(PEG32 beads：tween20＝100:1)，室温结合15min，期间吹打混匀一次；

f)不粘管置于磁力架3-5min后弃去上清，用700ul 75％乙醇洗涤两次，洗涤时将不粘管前后方向反转，使得beads在乙醇中游动，每次洗涤游动2-3次；

g)室温下晾干后用27ul TE/tween20回溶(TE:tween20＝500:1)，溶解时间共计15min，中间混匀一次；

h)上清转移到新的1.5mlEP管中，将最终得到产物用Qubit^TMssDNA Assay Kit定量。Buffer与染料比例为199:1混匀后votex并离心混合备用，取两份190ul稀释后染料工作液分别加入10ul的两种标准品votex并离心混合备用，取198ul稀释后染料工作液加入2ul样品，votex后并离心进行Qubit仪器定量。

i)浓度标准化(Normalization)

按照单链分子定量测定的浓度调整DNB制备使用的样本起始量统一调整为35.3ng－53ng，将对应体积样本(<60ul)转移至Biorad PCR板中，使用1XTE补齐使总体积不超过120ul。

终浓度为5.625-7.5fmol/ul，体积为120ul，则总量为35.3ng-53ng，1adapter测序中的DNB需要120fmol，7.5foml/ul，16ul。故需要把文库稀释至7.5fmol/ul。

10)CG 1-Adapter测序

利用CG平台的标准化流程测序。DNA纳米芯片是CG独创的一种高通量测序技术。该示例的对改进的单接头测序文库进行测序，较其他测序方案成本更低、速度更快，并集成质控确保其测序质量。

实施例二

对实施例一的下机数据进行处理。利用本发明的测序数据处理方法和/或CNV检测方法，基于CG平台测序技术，可以对超微量的cfDNA进行富集，文库建立、测序和数据分析工作。在该示例中，由于CG测序原理的特殊性，其测序的reads较短，且在特定位置存在重测序以及小gap的现象，难以不处理直接使用普通的比对软件对测序结果进行比对。针对reads的特殊结构，我们使用CG平台专有的TeraMap进行比对，其工作原理是：首先，它将读长中的两端(LeftArm，RightArm)分别做比对，其间TeraMap会尝试多种gap值来处理读长，以获得更多的比对结果；然后，将每一端的比对结果拿出来作为参考，对另一端做局部比对(例如4-AD，局部比对的范围是0～700bp)；如果两端可以良好比对到同一染色体，且insert-size符合期望(例如4-AD，一读段对的两读段的距离为0～700bp)，则只输出最佳比对结果，否则两端的多个比对结果全部输出。TeraMap是CG测序平台的比对软件，它可将CG特有序列比对到参考基因组上，其输出格式由三部分组成，简要说明如下：第一行，是reads序列信息；第二行和第三行，是reads比对情况说明；第四行和第五行，是reads比对结果详细信息。

第一行：

第二行：

第四行：

因为TeraMap比对存在gap问题，使得无法进行下游分析，依据本发明的方法开发Teramap2Sam软件，将TeraMap比对结果中gap去除并转换为SAM(sequence alignment/mapformat)。Teramap2Sam软件的主要流程可分为三部分，算法流程图如图10所示。

第一步：提取唯一比对结果。根据TeraMap输出结果matchCount判定是否唯一比对，同时要求插入片段长度满足要求以及两端read比对在用一条参考序列上。

第二步：去除gap。根据gaps字段判定reads中的gap位置，并修正read序列。

第三步：计算FLAG。根据双端read的比对方向，计算SAM文件中的FLAG参数，获得比对情况。

SAM是存储比对信息的一种较通用的格式，每一行是一个reads的比对结果，主要由十一个字段组成，其后还可添加更多字段包含更多信息，比如XT:A:U就是表示此reads为unique比对。简要说明如下：

在实际使用中为了节约存储资源，主要使用其二进制压缩格式(BAM)。此外CG又针对其reads结构开发了Assembly Software将reads重新组装，组装完成后进行变异检测分析等后续工作。

由于CG单接头reads的特殊结构存在reads太短(最短12bp)的缺点，在一些特殊数据的处理中CG原有的突变检测工具不再适用或者检测结果不佳。针对这种情况，我们首先开发工具将TeraMap的比对结果转化为通用的SAM/BAM格式，其中SAM/BAM是高通量测序中普遍使用的比对结果格式，所以我们采用这种通用的格式，然后再使用BAM数据检测拷贝数变异。目前已有的拷贝数变异检测方法有隐马可夫模型、环状二元分割、等级分割、核平滑算法等。我们根据总长达1,000,000bp的多个窗口的reads深度分布，使用z-score(标准分数)得到拷贝数变异结果。

考虑到实际测序过程中reads中GC含量会对测序深度有一定的影响，我们对比对结果(BAM)进行GC含量校正，消除GC含量对深度的影响。具体的，取总长达1,000,000bp的多个窗口算窗口的GC含量和平均测序深度，对GC-测序深度的数据进行lowess回归，根据回归曲线对GC含量进行校正。

标准分数(standard score)也叫z分数(z-score),是一个分数与平均数的差再除以标准差的过程。用公式表示为：z＝(x-μ)/σ。其中x为某一具体分数，μ为平均数，σ为标准差。Z值的量代表着原始分数和母体平均值之间的距离，是以标准差为单位计算。在原始分数低于平均值时Z则为负数，反之则为正数。在该示例中，通过对2000bp窗口内reads计数(原始分数)和总体reads平均值(多个正常对照样本)之间的距离使用标准差进行度量，可以有效检测出拷贝数变异。Z值为正数时反应为拷贝数大于2(正常样本是2倍体)，比如重复，z值为负数时反应拷贝数小于2，比如缺失。将该实施例中的上述的CNV检测方法编写成程序，并将该程序命名为calcu_zscore_query，将z绝对值大于3的区域判断为发生CNV。

较之传统方法，我们使用的基于CG单接头测序的方法可以实现超微量建库测序，建库只需要1-10ng核酸，需要外周血量2-5ml，并且CG的标准化流程简单快速，TeraMap比对结果转换为SAM格式之后比闭源的TeraMap格式更加通用，可以使用Samtools等软件进行处理。此外，使用z-score(标准分数)可以快速检测出CNV，50乘全基因组数据的CNV分析只需4小时，作为对比，CONTRA软件[http://sourceforge.net/projects/contra-cnv/]需要1天以上时间。

该示例中利用TeraMap进行比对。测序完成后使用CG平台的集成工具makeADF得到原始reads，然后用TeraMap进行比对，将测序得到的reads比对的参考序列上。得到的比对结果使用TeraMap2Sam转换为通用的SAM格式。表1展示结果。

表1

Claims (35)

1.一种测序数据处理装置，其特征在于，包括，

数据接收单元，用于接收所述测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读段对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂；

处理器，用于执行数据处理程序，执行所述数据处理程序包括实现将所述测序数据与参考序列作比对，获得比对结果，以及消除所述比对结果中的每个读段的缺口，获得通用比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，

所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果；以及，

至少一个存储单元，用于存储数据，其中包括所述数据处理程序；

所述作比对包括，

将每对读段对的左臂和右臂分别与所述参考序列比对，获得一级左比对结果和一级右比对结果，

分别以所述一级左比对结果和所述一级右比对结果的其中一个为参考，对另一个作比对，获得二级左比对结果和二级右比对结果，

基于所述二级左比对结果和所述二级右比对结果，获得多个所述读段对的比对结果，或者获得多个所述左臂的比对结果和多个所述右臂的比对结果；

消除所述比对结果中的每个读段的缺口包括，

若所述读段包含正缺口，以N填补所述正缺口的大小，

若所述读段包含负缺口，去除所述负缺口，其中，

N为A、T、C或G。

2.如权利要求1所述的装置，其特征在于，所述作比对包括，设置所述缺口的大小以使每个左臂或者每个右臂与所述参考序列进行多次比对。

3.如权利要求2所述的装置，其特征在于，所述每个左臂或者每个右臂与参考序列进行多次比对为，将所述每个左臂或者所述每个右臂的缺口分别设置为-3nt、-2nt、-1nt、0nt、1nt、2nt、3nt、4nt、5nt、6nt和7nt，获得对应的多个读段，分别将所述对应的多个读段与所述参考序列比对。

4.如权利要求1-3任一所述的装置，其特征在于，所述比对结果的格式为TeraMap。

5.如权利要求1-3任一所述的装置，其特征在于，执行所述数据处理程序还包括实现，在所述消除比对结果中的每个读段的缺口之前，提取所述比对结果中的唯一比对结果以替换所述比对结果，所述唯一比对结果包括唯一比对上所述参考序列的多个读段对，并且每一所述读段对比对到所述参考序列的相同染色体，每一所述读段对中的两个读段的距离符合所述染色体片段的两个位置的距离；

执行所述数据处理程序还包括实现，修正使所述唯一比对结果中的每一对读段对比对到所述参考序列的相同染色体的正链。

6.如权利要求1-3任一所述的装置，其特征在于，执行所述数据处理程序还包括实现，在所述消除比对结果中的每个读段的缺口之前，提取所述比对结果中的唯一比对结果以替换所述比对结果，所述唯一比对结果包括唯一比对上所述参考序列的多个读段对，并且每一所述读段对比对到所述参考序列的相同染色体，每一所述读段对中的两个读段的距离符合所述染色体片段的两个位置的距离；

执行所述数据处理程序还包括实现数据格式转换，所述数据格式转换包括转换所述比对结果或所述唯一比对结果的格式。

7.如权利要求1-3任一所述的装置，其特征在于，所述通用比对结果的格式为SAM或BAM。

8.一种测序数据处理系统，其包括一主机和一显示装置，其特征在于，所述系统还包括权利要求1-7任一所述的测序数据处理装置。

9.一种测序数据处理方法，其特征在于，包括如下步骤，

获取测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读长对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体片段的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂；

将所述测序数据与参考序列比对，获得比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，

所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果；

消除所述比对结果中的每个读段的缺口，获得通用比对结果；

所述比对包括，

将每对读段对的左臂和右臂分别与所述参考序列比对，获得一级左比对结果和一级右比对结果，

分别以所述一级左比对结果和所述一级右比对结果的其中一个为参考，对另一个作比对，获得二级左比对结果和二级右比对结果，

基于所述二级左比对结果和所述二级右比对结果，获得多个所述读段对的比对结果，或者获得多个所述左臂的比对结果和多个所述右臂的比对结果；

消除所述比对结果的每个读段的缺口包括，

若所述读段包含正缺口，以N填补所述正缺口的大小，

若所述读段包含负缺口，去除所述负缺口，其中，

N为A、T、C或G。

10.如权利要求9所述的方法，其特征在于，获取所述测序数据包括构建测序文库，获得测序文库，所述测序文库为单链环状DNA文库，所述测序文库由所述染色体片段的一条链和至少一个预定DNA序列构成。

11.如权利要求9所述的方法，其特征在于，所述每对读段分别来自所述染色体片段的两端。

12.如权利要求11所述的方法，其特征在于，所述获取测序数据包括测序文库构建，获得测序文库，所述测序文库为单链环状DNA文库，所述测序文库由所述染色体片段的一条链和连接所述一条链的两端的一个预定DNA序列构成。

13.如权利要求12所述的方法，其特征在于，构建所述测序文库包括，

(1)提取待测核酸；

(2)末端磷酸化所述核酸，获得末端磷酸化产物；

(3)末端修复所述末端磷酸化产物，获得末端修复产物；

(4)将第一序列和第二序列连接至所述末端修复产物的两端，获得第一连接产物；

(5)利用第三序列对所述连接产物进行缺刻平移和扩增，获得扩增产物，所述第三序列为一对引物对，所述引物对的至少一条引物带有生物素标记；

(6)利用所述生物素标记对所述扩增产物进行单链分离，获得单链产物；

(7)利用第四序列环化所述单链产物，获得所述测序文库；其中，

所述第四序列能够连接所述第一序列的一端和所述第二序列的一端，所述第一序列和/或所述第二序列的另一端为双脱氧核苷酸。

14.如权利要求12所述的方法，其特征在于，构建所述测序文库包括，

(1)提取待测核酸；

(2)末端修复所述核酸，获得末端修复产物；

(3)末端磷酸化所述末端修复产物，获得末端磷酸化产物；

(4)将第一序列和第二序列连接至所述末端磷酸化产物的两端，获得第一连接产物；

(5)利用第三序列对所述连接产物进行缺刻平移和扩增，获得扩增产物，所述第三序列为一对引物对，所述引物对的至少一条引物带有生物素标记；

(6)利用所述生物素标记对所述扩增产物进行单链分离，获得单链产物；

(7)利用第四序列环化所述单链产物，获得所述测序文库；其中，

所述第四序列能够连接所述第一序列的一端和所述第二序列的一端，所述第一序列和/或所述第二序列的另一端为双脱氧核苷酸。

15.如权利要求9-14任一所述的方法，其特征在于，所述比对包括，设置所述缺口的大小以使每个左臂或者每个右臂与所述参考序列进行多次比对。

16.如权利要求15所述的方法，其特征在于，所述每个左臂或者每个右臂与参考序列进行多次比对为，将所述每个左臂或者所述每个右臂的缺口分别设置为-3nt、-2nt、-1nt、0nt、1nt、2nt、3nt、4nt、5nt、6nt和7nt，获得对应的多个读段，分别将所述对应的多个读段与所述参考序列比对。

17.如权利要求9-14任一所述的方法，其特征在于，所述比对结果的格式为TeraMap。

18.如权利要求9-14任一所述的方法，其特征在于，在所述消除比对结果中的每个读段的缺口之前，提取所述比对结果中的唯一比对结果以替换所述比对结果，所述唯一比对结果包括唯一比对上所述参考序列的多个读段对，并且每一所述读段对比对到所述参考序列的相同染色体，每一所述读段对中的两个读段的距离符合所述染色体片段的大小。

19.如权利要求18所述的方法，其特征在于，对所述唯一比对结果进行修正，以使所述唯一比对结果中的每一对读段对比对到所述参考序列的相同染色体的正链。

20.如权利要求18所述的方法，其特征在于，获得所述通用比对结果还包括，对所述比对结果或所述唯一比对结果进行数据格式转换。

21.如权利要求9-14任一方法，其特征在于，所述通用比对结果的格式为SAM或BAM。

22.一种计算机可读存储介质，其特征在于，用于存储供计算机执行的程序，所述程序的执行包括完成权利要求9-21任一方法。

23.一种检测CNV的方法，其特征在于，包括，

a.获取待测样本的核酸；

b.对所述核酸进行测序，获得测序数据；

c.对所述测序数据进行处理，以获得通用比对结果；

d.基于所述通用比对结果检测所述CNV；其中，c步骤利用权利要求1-8任一测序数据处理装置进行。

24.如权利要求23所述的方法，其特征在于，b步骤包括，对所述核酸进行测序文库构建，获得测序文库，所述测序文库为单链环状DNA文库。

25.如权利要求24所述的方法，其特征在于，所述测序文库构建包括，

末端磷酸化所述核酸，获得末端磷酸化产物；

末端修复所述末端磷酸化产物，获得末端修复产物；

将第一序列和第二序列连接至所述末端修复产物的两端，获得第一连接产物；

利用第三序列对所述连接产物进行缺刻平移和扩增，获得扩增产物，所述第三序列为一对引物对，所述引物对的至少一条引物带有生物素标记；

利用所述生物素标记对所述扩增产物进行单链分离，获得单链产物；

利用第四序列环化所述单链产物，获得所述测序文库，其中，

所述第四序列能够连接所述第一序列的一端和所述第二序列的一端，所述第一序列和/或所述第二序列的另一端为双脱氧核苷酸。

26.如权利要求24所述的方法，其特征在于，所述测序文库构建包括，

末端修复所述核酸，获得末端修复产物；

末端磷酸化所述末端修复产物，获得末端磷酸化产物；

将第一序列和第二序列连接至所述末端磷酸化产物的两端，获得第一连接产物；

利用第三序列对所述连接产物进行缺刻平移和扩增，获得扩增产物，所述第三序列为一对引物对，所述引物对的至少一条引物带有生物素标记；

利用所述生物素标记对所述扩增产物进行单链分离，获得单链产物；

利用第四序列环化所述单链产物，获得所述测序文库，其中，

所述第四序列能够连接所述第一序列的一端和所述第二序列的一端，所述第一序列和/或所述第二序列的另一端为双脱氧核苷酸。

27.如权利要求23-26任一所述的方法，其特征在于，所述测序是利用组合探针锚定连接测序技术进行的。

28.如权利要求23所述的方法，其特征在于，d步骤包括，

在所述参考序列上设置多个窗口，基于所述通用比对结果中匹配到所述窗口的读段的量与对照样本的通用比对结果中匹配到相同窗口的读段的量的差异具有显著性，判定所述待测样本核酸存在所述CNV，其中，

所述窗口为所述参考序列的一部分。

29.如权利要求28所述的方法，其特征在于，所述对照样本的通用比对结果是通过权利要求12-26任一测序数据处理方法获得的。

30.如权利要求23所述的方法，其特征在于，d步骤包括，

在所述参考序列上设置多个窗口，计算窗口的测序深度，窗口的测序深度＝所述通用比对结果中比对到所述窗口的读段的数量/所述窗口的大小；

利用测序深度和GC含量的关系校正所述窗口的测序深度，获得窗口的校正测序深度；

基于所述窗口的校正测序深度与对照样本的相同窗口的校正测序深度的差异具有显著性，判定所述待测样本核酸存在所述CNV，其中，

所述窗口为所述参考序列的一部分。

31.如权利要求28-30任一所述的方法，其特征在于，所述对照样本的个数不小于30个。

32.如权利要求30所述的方法，其特征在于，所述测序深度和GC含量的关系的建立包括，

获得多个对照样本核酸的测序数据，所述测序数据由多个读段组成；

在所述参考序列上设置多个窗口，将所述多个对照样本的测序数据分别与所述参考序列的窗口比对，计算各个对照样本的测序数据中比对上每个窗口的读段的数目，获得每个窗口的测序深度，所述窗口为所述参考序列的一部分，所述窗口的测序深度＝各个对照样本的比对上所述窗口的读段的总数目/(对照样本个数*所述窗口的大小)；

基于每个窗口的测序深度和该窗口的GC含量，利用二维回归分析法建立所述测序深度和GC含量的关系。

33.如权利要求32的方法，其特征在于，所述二维回归分析法为局部加权回归散点平滑法。

34.如权利要求30的方法，其特征在于，所述对照样本的相同窗口的校正测序深度是利用所述测序深度和GC含量的关系校正对照样本的相同窗口的测序深度获得的，所述对照样本的相同窗口的测序深度＝所述对照样本的测序数据中比对到所述窗口的读段的数目/所述窗口的大小。

35.一种CNV检测设备，其特征在于，包括，

核酸获取装置，用以获取待测样本的核酸；

测序装置，用以对来自所述核酸获取单元的核酸进行测序，获得测序数据；

数据处理装置，用于对来自所述测序装置的测序数据进行处理，以获得通用比对结果；

检测装置，用于基于来自所述数据处理装置的通用比对结果检测所述CNV；其中，

所述数据处理装置包括，

数据接收单元，用于接收来自所述测序装置的测序数据，所述测序数据包括多对读段对，每对读段对由两个读段组成，分别来自一条染色体片段的两个位置，每对读长对中的两个读段分别来自所述染色体片段的正链和负链，或者每对读长对中的两个读段都来自所述染色体片段的正链或所述染色体的负链，每个读段都包含缺口，将一对读段对的两个读段分别定义为左臂和右臂，

处理器，用于执行数据处理程序，执行所述数据处理程序包括实现将所述测序数据与参考序列作比对，获得比对结果，以及消除所述比对结果中的每个读段的缺口，获得通用比对结果，所述比对结果包括多个所述读段对的比对结果，和/或，所述比对结果包括多个所述左臂的比对结果和多个所述右臂的比对结果，以及，

至少一个存储单元，用于存储数据，其中包括所述数据处理程序；

所述作比对包括，

将每对读段对的左臂和右臂分别与所述参考序列比对，获得一级左比对结果和一级右比对结果，

分别以所述一级左比对结果和所述一级右比对结果的其中一个为参考，对另一个作比对，获得二级左比对结果和二级右比对结果，

基于所述二级左比对结果和所述二级右比对结果，获得多个所述读段对的比对结果，或者获得多个所述左臂的比对结果和多个所述右臂的比对结果；

消除所述比对结果中的每个读段的缺口包括，

若所述读段包含正缺口，以N填补所述正缺口的大小，

若所述读段包含负缺口，去除所述负缺口，其中，

N为A、T、C或G。

Images