CN110060734A - 一种高鲁棒性dna测序用条形码生成和读取方法 - Google Patents
一种高鲁棒性dna测序用条形码生成和读取方法 Download PDFInfo
- Publication number
- CN110060734A CN110060734A CN201910252268.7A CN201910252268A CN110060734A CN 110060734 A CN110060734 A CN 110060734A CN 201910252268 A CN201910252268 A CN 201910252268A CN 110060734 A CN110060734 A CN 110060734A
- Authority
- CN
- China
- Prior art keywords
- sequencing
- sequence
- code
- bar code
- barcode
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 38
- 238000001712 DNA sequencing Methods 0.000 title claims abstract description 30
- 238000012163 sequencing technique Methods 0.000 claims abstract description 70
- 238000003780 insertion Methods 0.000 claims abstract description 53
- 230000037431 insertion Effects 0.000 claims abstract description 53
- 238000012217 deletion Methods 0.000 claims abstract description 41
- 230000037430 deletion Effects 0.000 claims abstract description 41
- 125000004122 cyclic group Chemical group 0.000 claims abstract description 24
- 238000006467 substitution reaction Methods 0.000 claims abstract description 24
- 238000012937 correction Methods 0.000 claims abstract description 22
- 230000015572 biosynthetic process Effects 0.000 claims abstract description 9
- 238000003786 synthesis reaction Methods 0.000 claims abstract description 9
- 108020004414 DNA Proteins 0.000 claims description 15
- 238000012165 high-throughput sequencing Methods 0.000 claims description 9
- 230000008439 repair process Effects 0.000 claims description 8
- 230000002194 synthesizing effect Effects 0.000 claims description 8
- 230000005540 biological transmission Effects 0.000 claims description 6
- 238000001514 detection method Methods 0.000 claims description 5
- 238000013507 mapping Methods 0.000 claims description 5
- 102000012410 DNA Ligases Human genes 0.000 claims description 4
- 108010061982 DNA Ligases Proteins 0.000 claims description 4
- 102000004190 Enzymes Human genes 0.000 claims description 4
- 108090000790 Enzymes Proteins 0.000 claims description 4
- 230000009471 action Effects 0.000 claims description 4
- 238000000137 annealing Methods 0.000 claims description 4
- 238000005516 engineering process Methods 0.000 claims description 4
- 238000007476 Maximum Likelihood Methods 0.000 claims description 3
- 238000013138 pruning Methods 0.000 claims description 3
- 239000003550 marker Substances 0.000 claims 1
- 238000003672 processing method Methods 0.000 abstract description 4
- 238000002156 mixing Methods 0.000 abstract description 3
- 238000012216 screening Methods 0.000 abstract 1
- 238000011109 contamination Methods 0.000 description 3
- 238000004088 simulation Methods 0.000 description 3
- 230000008901 benefit Effects 0.000 description 2
- 238000010276 construction Methods 0.000 description 2
- 238000002515 oligonucleotide synthesis Methods 0.000 description 2
- 102000053602 DNA Human genes 0.000 description 1
- 230000006820 DNA synthesis Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 230000037433 frameshift Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000035772 mutation Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000011084 recovery Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000009897 systematic effect Effects 0.000 description 1
- 238000007671 third-generation sequencing Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B25/00—ICT specially adapted for hybridisation; ICT specially adapted for gene or protein expression
- G16B25/10—Gene or protein expression profiling; Expression-ratio estimation or normalisation
Landscapes
- Health & Medical Sciences (AREA)
- Genetics & Genomics (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Biophysics (AREA)
- Molecular Biology (AREA)
- Engineering & Computer Science (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Biotechnology (AREA)
- Evolutionary Biology (AREA)
- General Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Theoretical Computer Science (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
- Error Detection And Correction (AREA)
Abstract
本发明公开了一种高鲁棒性DNA测序用条形码生成和读取方法，属于基因组学的高通量测序领域，将已知伪随机序列与相同长度的短码长循环码逐比特组合成比特对，比特对映射为一个碱基，不同的碱基序列经过筛选后就构成了测序用条形码；进一步，提出用循环移位和动态规划判断接收测序用条形码发生插入、删节错误的位置，修正接收的纠错码字并进行译码的处理方法。混样测序中，该发明对测序条形码的各类错误均具有良好纠错能力，能够对抗合成、建库以及测序中的插入(insertion)、删节(deletion)与替代错误，具有很高的鲁棒性。
Description
技术领域
本发明涉及基因组学的高通量测序领域，尤其涉及一种高鲁棒性DNA测序用条形码生成和读取方法。
背景技术
下一代高通量测序(HTS)在单次上机测序中产生大量读段，并且测序速度不断提高，促进了生物医学领域的新发现。为充分利用测序设备的高通量测序能力，可以将多个包含数百万个不同DNA分子的DNA样品混合在一起，进行大规模并行测序，然后基于DNA分子上一段用于区分不同样品来源的独特的碱基序列—DNA测序用条形码(或者称为“索引”)，将测序后的所有读段重新分配到其各自的样品。这种多样本复用的测序方式实现了测序成本的进一步降低，但是也存在不同样本之间相互干扰的问题。因此，保证DNA测序条形码的鲁棒性，对保证高通量多样本并行测序的准确性具有重要意义。
研究发现，一些测序平台，例如Illumina公司的测序平台HiSeq 3000/4000、HiSeqX Ten以及NovaSeq，混样上机测序会出现测序条形码(或者称为索引)的“错误分配”问题，也称为“索引跳变”(index hopping)，样本错误分配率达到0.2％～6％；且对于PCR-free文库，即使采用双索引策略，平均污染率也只能降到0.08％。索引污染可能由于实验操作、PCR错误、测序错误、寡核苷酸合成错误引入，这些错误通过混淆样品来源来降低测序准确性。最常见的DNA合成错误是单碱基的删节(Deletion)，删节错误对译码尤其具有挑战性，因为它会在下游测序中引起移码；在大规模寡核苷酸合成期间，替代和插入错误也是存在的。而且，高通量测序本身也会存在10-3～10-4的错误率。三代测序的错误更是接近10％。因此，设计能对抗插入、删节与替代错误的高鲁棒性DNA测序用条形码对降低“索引跳变”率有重要意义。提高条形码的鲁棒性可通过纠错码来实现。通过纠正测序用条形码上的插入、删节和替代错误，降低索引污染率，增加正确识别的样本总数，从而提高整体测序效率。
目前，已经提出了不同结构的具有纠错功能的DNA测序用条形码。2008年，Hamady等人在Nature Method发表论文首次提出了基于二进制汉明码的测序用条形码；Krishnan等人对Hamady提出的条形码进行了改进，进一步提出了基于BCH码的条形码，具有更大的码字间最小距离和更高的纠错能力。但是，上述两种方法仅仅能纠正传统替代错误。Adey与Morrison首次提出了基于Levenshtein距离的测序用条形码，能够纠正插入、删节和替代错误，但是没有已知的系统的码构造方法可在编辑距离下给出最佳码字。Ashlock和Houghten在论文中提出利用贪心算法来找到具有高的编辑距离的适当大小的条形码集，可以满足生物约束，然而没有给出实际译码步骤。Buschmann和Bystrykh将传统的Levenshtein距离度量进行改进，称为“Sequence-Levenshtein”码，将其作为测序用条形码，该方法能够检测并恢复损坏的条形码的边界，并且平均校正比传统的Levenshtein或汉明码更多的随机突变。David和Steffenr利用Davey和Mackay提出的水印码作为测序用条形码，具有良好的同步能力，不需要任何标记就能够恢复测序读段中条形码的位置，而且具有良好的插入删节纠错能力，但是这种方法的译码复杂度较高，短码的性能并不具有优势，而测序条码一般只能使用非常短的序列。
进一步，实际测序中条形码紧邻DNA文本，此时条形码上单独的插入或删节错误会引起序列长度的变化，由于无“停止词”将条形码与DNA文本分开，表现为条形码末端发生等数目的删节或插入错误。
发明内容
本发明提供了一种高鲁棒性DNA测序用条形码生成和读取方法，本发明以较低的复杂度，获得了良好的纠正插入、删节错误的性能，详见下文描述：
一种高鲁棒性DNA测序用条形码生成和读取方法，所述方法包括以下步骤：
(1)用长度为k的信息序列mi表示不同的DNA测序样本，信息序列mi采用短码长分组循环码生成长度为n的编码码字c；然后产生长度为n的伪随机序列p；将伪随机序列p与编码码字c进行逐比特组合，按照预设规则将比特对映射为碱基得到长度为n的碱基序列b，即得到最多支持2k的不同样本混合测序的DNA测序条形码；
(2)不同样本的信息序列生成不同的测序条形码，对各个样本分别进行样本建库；利用引物合成技术合成测序条形码正链和反链，合成时在正链3’端添加一个T碱基；通过退火正链和反链，合成具有3’端T尾粘性末端的双链条形码DNA；利用DNA末端修复酶对已打断的待测样本进行尾端修复同时在3’端添加A尾粘性末端；在DNA连接酶的作用下将条形码添加至待测样本；混合样本上机进行高通量测序，得到包含不同错误的混合样品接收序列；
(3)根据解映射后得到的对应已知伪随机序列的部分，通过相关检测或动态规划法获取伪随机序列窗口，进而得到测序用条形码的窗口，接收窗口的测序用条形码解映射后，得到接收随机序列s和接收码字d；
(4)将s和d依次左循环移位i次得s(i)和d(i)，将p依次左循环移位k次，得到p(k)，比较p(k)和s(i)，利用动态规划标记插入/删节位置，对d(i)进行修正并纠错纠删译码，产生候选码字集选择译码成功的码字作为最终译码结果据此将测序读段分配至不同样本。
其中，所述预设规则为：00→A，01→T，10→G，11→C。
进一步，所述解映射规则为：A→00，T→01，G→10，C→11。
其中，所述利用动态规划标记插入/删节位置，对d(i)进行修正的操作具体为：
1)利用动态规划比较p(k)和s(i)，计算序列之间的Levenstein距离并同步记录其对应的插入删节错误的位置；
2)若某位置标记为插入，则删掉d(i)对应位置的比特；若标记为删节，则在d(i)对应位置补入比特0或比特1；得到修正后的码字dc (i)。
进一步地，所述纠错纠删译码具体为：
(1)对dc (i)右循环移位i次，输入译码器进行纠错纠删译码，产生译码码字
(2)计算译码码字的伴随式S(x)，若S(x)＝0，译码成功，执行步骤(5)；否则执行步骤(3)；
(3)判断k＝i+c是否成立，若成立，结束伪随机序列p的循环移位，执行步骤(4)；否则k增加1，执行伪随机序列左循环移位k次；
(4)判断i＝n-1是否成立，若成立，结束接收序列s和接收码字d循环移位，并执行步骤(5)；否则令i增加1，执行接收序列左循环移位i次；
(5)输出译码码字至此译码结束。
进一步地，所述利用动态规划比较p(k)和s(i)，计算序列之间的Levenstein距离并同步记录其对应的插入删节错误的位置具体为：
计算递归计算网格图中到达每一节点(i,j)的三种路径的最小距离值d(xi,yj)时，同步记录得到最小距离值所对应的路径，即上个节点的坐标值；
当i＝n,j＝n，终止递归，此时d(x,y)＝d(xn,yn)即为发送序列与接收序列的Levenshtein距离；
从节点(n,n)出发，利用pointer(i,j)中的信息，寻找一条到达节点(0,0)的最大似然路径并同步记录此路径上的插入、删节和替代错误；
若pointer(i,j)＝(i-1,j)表示发生删节；若pointer(i,j)＝(i-1,j-1)表示正确传输或替代，此时再判断xi和yj的值，若xi＝yj，则表示正确传输，若xi≠yj，则表示替代；若pointer(i,j)＝(i,j-1)，则表示发生插入。
其中，所述最小距离值d(xi,yj)具体为：
其中，表示删节的权重；若xi＝yj，w(xi,yj)＝0，若xi≠yj，w(xi,yj)表示替代错误的权重；表示插入的权重。
进一步地，
若则pointer(i,j)＝(i-1,j)；
若d(xi,yj)＝d(xi-1,yj-1)+w(xi,yj)，则pointer(i,j)＝(i-1,j-1)；
若则pointer(i,j)＝(i,j-1)。
本发明提供的技术方案的有益效果是：
1、本发明通过循环移位和动态规划相结合的方法，提高了判断插入、删节位置的准确性；
2、本发明提高了测序用条形码对混样测序中插入、删节错误的纠错能力，具有非常高的鲁棒性。
附图说明
图1为高鲁棒性DNA测序用条形码生成和读取方法的流程图；
图2为DNA测序用条形码处理方法的流程图；
图3为动态规划法标记插入、删节位置的流程图；
图4为本方法用于测序用条形码(BCH(31,11,5))时，在加入1个插入、1个删节和若干替代错误以及2个插入、2个删节和若干替代错误时，误帧率与错误数目关系曲线图；
图5为本方法用于测序用条形码(BCH(31,11,5))时，在两种错误模型下，误帧率与错误数目关系曲线图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面对本发明实施方式作进一步地详细描述。
针对上述问题，本发明实施例提出了一种低复杂度、能够对抗插入、删节和替代错误的高鲁棒性DNA测序用条形码的生成和读取处理方法，其对测序用条形码上随机发生的插入、删节错误具有良好的纠错能力。同时，本发明实施例中的测序条形码，对上述DNA文本中测序用条形码末端发生插入/删节错误的情形，也具有更加优越的性能，具有更好的实用性。
实施例1
参见图1，高鲁棒性DNA测序用条形码生成和读取方法的具体步骤为：
(1)用长度为k(k为正整数)的信息序列mi(i＝1,2,...,2k)表示不同的DNA测序样本，信息序列mi采用短码长分组纠错码生成长度为n的编码码字c；然后产生长度为n的伪随机序列p；将伪随机序列p与编码码字c进行逐比特组合成比特对，按照{00→A，01→T，10→G，11→C}规则将比特对映射为一个碱基，得到长度为n的碱基序列b，即得到最多可支持2k的不同样本混合测序的DNA测序条形码；
(2)不同样本的信息序列生成不同的测序条形码，对各个样本分别进行样本建库；利用引物合成技术合成测序条形码正链和反链，合成时在正链3’端添加一个T碱基；通过退火正链和反链，合成具有3’端T尾粘性末端的双链条形码DNA；利用DNA末端修复酶对待测样本进行尾端修复同时在3’端添加A尾粘性末端；在DNA连接酶的作用下，借助A/T末端粘性，将条形码添加至待测样本；在此之后完成添加测序接头等步骤完成样本建库；然后混合样本上机进行高通量测序，得到包含不同错误的混合样品接收序列；
(3)对接收序列进行解映射，解映射规则为A→00，T→01，G→10，C→11，然后根据解映射后得到的对应已知伪随机序列的部分，通过相关检测或动态规划法，得到伪随机序列窗口，从而得到测序用条形码的窗口，区分出样本以及样本的序列区域；接收窗口的测序用条形码解映射后，得到接收随机序列s和接收码字d；
(4)为识别插入/删节错误位置，将s和d依次左循环移位i次(i＝0,1,...,n-1)得s(i)和d(i)，将p依次左循环移位k次(k∈[i-c,i+c]，c为正常数)，得到p(k)，
然后比较p(k)和s(i)，利用动态规划标记可能的插入/删节位置，据此对d(i)进行修正，并进行纠错纠删译码，产生候选码字集最后选择译码成功的码字作为最终译码结果
参见图2，为识别插入/删节错误位置，将s和d依次左循环移位i次(i＝0,1,...,n-1)得s(i)和d(i)，将p依次左循环移位k次(k∈[i-c,i+c]，c为正常数)，得到p(k)，
然后比较p(k)和s(i)，利用动态规划标记可能的插入/删节位置，据此对d(i)进行修正，并进行纠错纠删译码，产生候选码字集最后选择译码成功的码字作为最终译码结果的具体步骤为：
(4.1)初始化接收序列和接收码字循环移位的次数i＝0；
(4.2)将接收序列s＝(s0,s1,...,sn-1)依次左循环移位i次，
得到s(i)＝(si,...,sn-1,s0,...,si-1)，相应地，将接收码字d＝(d0,d2,...,dn-1)依次左循环移位i次，得到与s(i)对应的d(i)＝(di,...,dn-1,d0,...,di-1)，其中n表示序列长度；
(4.3)初始化已知伪随机序列循环移位的次数k＝0；
(4.4)将伪随机序列p＝(p0,p1,...,pn-1)依次左循环移位k次(k∈[i-c,i+c]，c为常数)，得到p(k)＝(pk,...,pn-1,p0,...,pk-1)；
(4.5)利用动态规划法比较p(k)和s(i)，计算序列之间的Levenstein距离并同步记录其对应的插入删节错误的位置；
(4.6)根据标记的插入删节位置，对接收码字d(i)进行修正，即若某位置标记为插入，则删掉d(i)对应位置的比特；若标记为删节，则在d(i)对应位置补入比特0或比特1；最后得到修正后的码字dc (i)；
(4.7)对dc (i)右循环移位i次，然后输入译码器进行纠错纠删译码，产生译码码字
其中，二进制纠错码的纠错纠删译码方法是在纠错译码方法的基础上，加入判断语句，即先将d(i)在标记的删节位置全填上0后得dc (i)，将dc (i)右循环移位i次后进行纠错译码，若译码器译码成功，则输出译码码字；若译码器译码错误，这时把dc (i)在标记的对应删节位置全填上1后得dc (i)，将dc (i)右循环移位i次再送入译码器，译得的码字即为所求。
(4.8)计算译码码字的伴随式S(x)，若S(x)＝0，则说明译码成功，执行步骤(4.11)；若S(x)≠0，则译码失败，执行步骤(4.9)；
(4.9)判断k＝i+c是否成立，如果k＝i+c成立，结束伪随机序列p的循环移位，并执行步骤(4.10)；否则令k增加1，并返回步骤(4.4)，其中n为序列长度；
(4.10)判断i＝n-1是否成立，如果i＝n-1成立，则结束接收序列s和接收码字d循环移位，并执行步骤(4.11)；否则令i增加1，并返回步骤(4.2)，其中n为序列长度；
(4.11)输出译码码字至此译码结束。
参见图3，利用动态规划法比较p(k)和s(i)，计算序列之间的Levenstein距离并同步记录其对应的插入删节错误的位置的具体步骤为：
令x＝p(k)，y＝s(i)，则x＝(x1,x2,...,xl)和y＝(y1,y2,...,yl)分别表示长度为n的发送序列和接收序列；令xi和yj分别表示序列x与y的前i个和前j个元素，即xi＝(x1,x2,...,xi)，yj＝(y1,y2,...,yj)，0≤i,j≤n。
(4.5.1)初始化原序列序号i＝0和接收序列序号j＝0；
(4.5.2)初始化Levenshtein距离，令d(x0,y0)＝0；
(4.5.3)计算网格图中到达每一节点(i,j)的三种路径的最小距离值：
其中，表示删节的权重；若xi＝yj，则w(xi,yj)＝0，若xi≠yj，则w(xi,yj)表示替代错误的权重；表示插入的权重。这里取(其中xi≠yj)；
(4.5.4)在计算d(xi,yj)时，同步记录得到最小距离值所对应的路径，即上个节点的坐标值；
若则pointer(i,j)＝(i-1,j)；
若d(xi,yj)＝d(xi-1,yj-1)+w(xi,yj)，则pointer(i,j)＝(i-1,j-1)；
若则pointer(i,j)＝(i,j-1)；
(4.5.5)判断j＝n是否成立，如果j＝n，结束j的循环；否则，令接收序列序号j增加1，并返回步骤(3)；
(4.5.6)判断i＝n是否成立，如果i＝n，结束计算距离值，此时d(x,y)＝d(xn,yn)即为发送序列与接收序列的Levenshtein距离，并执行步骤(7)；否则，令原序列序号i增加1，并返回步骤(3)；
(4.5.7)回溯，从节点(n,n)出发，通过利用pointer(i,j)中的信息，寻找一条到达节点(0,0)的最大似然路径并同步记录此路径上的插入、删节和替代错误，即若pointer(i,j)＝(i-1,j)表示发生删节；若pointer(i,j)＝(i-1,j-1)表示正确传输或替代，此时再判断xi和yj的值，若xi＝yj，则表示正确传输，若xi≠yj，则表示替代；若pointer(i,j)＝(i,j-1)，则表示发生插入。
实施例2
在本发明的实施例中，通过对BCH(31,11,5)构造的DNA测序用条形码进行纠错性能仿真，验证了上述DNA测序条形码生成和处理方法在各类错误情况下均具有良好的纠错性能。
在步骤(1)中，用特定样本对应的长度为11的信息序列mi(i＝1,2,...,211)进行BCH(31,11,5)编码，生成长度为31的编码码字c；然后产生长度为31的伪随机序列p；将伪随机序列p与编码码字c进行逐比特组合成比特对，按照{00→A，01→T，10→G，11→C}规则将一组比特对映射为一个碱基，得到长度为31的碱基序列b，即得到最多可支持211的不同样本混合测序的DNA测序条形码；
在步骤(2)中，不同样本的信息序列生成不同的测序条形码，对各个样本分别进行样本建库；利用引物合成技术合成测序条形码正链和反链，合成时在正链3’端添加一个T碱基；通过退火正链和反链合成具有3’端T尾粘性末端的双链条形码DNA；利用DNA末端修复酶对待测样本进行尾端修复同时在3’端添加A尾粘性末端；在DNA连接酶的作用下，借助A/T末端粘性，将条形码添加至待测样本；在此之后完成添加测序接头等步骤完成样本建库；然后混合样本上机进行高通量测序，得到包含不同错误的混合样品接收序列；
在步骤(3)中，对接收序列进行解映射，解映射规则为A→00，T→01，G→10，C→11，然后根据解映射后得到的对应已知伪随机序列的部分，通过相关检测或动态规划法，得到伪随机序列窗口，从而得到测序用条形码的窗口，区分出样本以及样本的序列区域；接收窗口的测序用条形码解映射后，得到接收随机序列s和接收码字d；
在步骤(4)中，为识别插入/删节错误位置，将s和d依次左循环移位i次(i＝0,1,...,30)，得s(i)和d(i)，
将p依次左循环移位k次(k∈[i-2,i+2])，得到p(k)，然后比较p(k)和s(i)，利用动态规划标记可能的插入/删节位置，据此对d(i)进行修正，并进行纠错纠删译码，产生候选码字集最后选择译码成功的码字作为最终译码结果
图4中给出了两种错误情形下(i)随机加入1个插入，1个删节和若干替代错误(ii)随机加入2个插入，2个删节和若干替代错误，利用条形码处理方案得到的误帧率与错误数目关系图。图5给出了在两种错误模型下(i)测序用条形码上随机发生相同数目的插入和删节错误(ii)将条形码嵌入到DNA文本中时，条形码上仅发生插入或删节错误，末端补入等数目的删节或者插入错误，误帧率与错误数目关系图。
仿真结果表明，BCH(31,11,5)构造的测序用条形码在随机发生1个插入、1个删节和若干替代错误以及随机发生2个插入，2个删节和若干替代错误的情形下均具有良好的纠正插入、删节和替代错误的能力。而且在错误模型(ii)下，误码率均明显低于(i)，说明设计的测序用条形码在实际的场景中性能更为优越。
综上所述，本发明实施例提出了一种利用循环移位和动态规划判断接收测序用条形码上发生插入、删节位置，修正接收纠错码字并进行译码的处理方法。仿真结果表明，该发明对混样测序中测序条形码的各类错误均具有良好纠错能力。
本领域技术人员可以理解附图只是一个优选实施例的示意图，上述本发明实施例序号仅仅为了描述，不代表实施例的优劣。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
Claims (8)
1.一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，所述方法包括以下步骤：
(1)用长度为k的信息序列mi表示不同的DNA测序样本，信息序列mi采用短码长分组循环码生成长度为n的编码码字c；然后产生长度为n的伪随机序列p；将伪随机序列p与编码码字c进行逐比特组合，按照预设规则将比特对映射为碱基得到长度为n的碱基序列b，即得到最多支持2k的不同样本混合测序的DNA测序条形码；
(2)不同样本的信息序列生成不同的测序条形码，对各个样本分别进行样本建库；利用引物合成技术合成测序条形码正链和反链，合成时在正链3’端添加一个T碱基；通过退火正链和反链，合成具有3’端T尾粘性末端的双链条形码DNA；利用DNA末端修复酶对已打断的待测样本进行尾端修复同时在3’端添加A尾粘性末端；在DNA连接酶的作用下将条形码添加至待测样本；混合样本上机进行高通量测序，得到包含不同错误的混合样品接收序列；
(3)根据解映射后得到的对应已知伪随机序列的部分，通过相关检测或动态规划法获取伪随机序列窗口，进而得到测序用条形码的窗口，接收窗口的测序用条形码解映射后，得到接收随机序列s和接收码字d；
(4)将s和d依次左循环移位i次得s(i)和d(i)，将p依次左循环移位k次，得到p(k)，比较p(k)和s(i)，利用动态规划标记插入/删节位置，对d(i)进行修正并纠错纠删译码，产生候选码字集选择译码成功的码字作为最终译码结果据此将测序读段分配至不同样本。
2.根据权利要求1所述的一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，所述预设规则为：00→A，01→T，10→G，11→C。
3.根据权利要求1所述的一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，所述解映射规则为：A→00，T→01，G→10，C→11。
4.根据权利要求1所述的一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，所述利用动态规划标记插入/删节位置，对d(i)进行修正的操作具体为：
1)利用动态规划比较p(k)和s(i)，计算序列之间的Levenstein距离并同步记录其对应的插入删节错误的位置；
2)若某位置标记为插入，则删掉d(i)对应位置的比特；若标记为删节，则在d(i)对应位置补入比特0或比特1；得到修正后的码字dc (i)。
5.根据权利要求1所述的一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，所述纠错纠删译码具体为：
(1)对dc (i)右循环移位i次，输入译码器进行纠错纠删译码，产生译码码字
(2)计算译码码字的伴随式S(x)，若S(x)＝0，译码成功，执行步骤(5)；否则执行步骤(3)；
(3)判断k＝i+c是否成立，若成立，结束伪随机序列p的循环移位，执行步骤(4)；否则k增加1，执行伪随机序列左循环移位k次；
(4)判断i＝n-1是否成立，若成立，结束接收序列s和接收码字d循环移位，并执行步骤(5)；否则令i增加1，执行接收序列左循环移位i次；
(5)输出译码码字至此译码结束。
6.根据权利要求4所述的一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，所述利用动态规划比较p(k)和s(i)，计算序列之间的Levenstein距离并同步记录其对应的插入删节错误的位置具体为：
计算递归计算网格图中到达每一节点(i,j)的三种路径的最小距离值d(xi,yj)时，同步记录得到最小距离值所对应的路径，即上个节点的坐标值；
当i＝n,j＝n，终止递归，此时d(x,y)＝d(xn,yn)即为发送序列与接收序列的Levenshtein距离；
从节点(n,n)出发，利用pointer(i,j)中的信息，寻找一条到达节点(0,0)的最大似然路径并同步记录此路径上的插入、删节和替代错误；
若pointer(i,j)＝(i-1,j)表示发生删节；若pointer(i,j)＝(i-1,j-1)表示正确传输或替代，此时再判断xi和yj的值，若xi＝yj，则表示正确传输，若xi≠yj，则表示替代；若pointer(i,j)＝(i,j-1)，则表示发生插入。
7.根据权利要求6所述的一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，所述最小距离值d(xi,yj)具体为：
其中，表示删节的权重；若xi＝yj，w(xi,yj)＝0，若xi≠yj，w(xi,yj)表示替代错误的权重；表示插入的权重。
8.根据权利要求7所述的一种高鲁棒性DNA测序用条形码生成和读取方法，其特征在于，
若则pointer(i,j)＝(i-1,j)；
若d(xi,yj)＝d(xi-1,yj-1)+w(xi,yj)，则pointer(i,j)＝(i-1,j-1)；
若则pointer(i,j)＝(i,j-1)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910252268.7A CN110060734B (zh) | 2019-03-29 | 2019-03-29 | 一种高鲁棒性dna测序用条形码生成和读取方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910252268.7A CN110060734B (zh) | 2019-03-29 | 2019-03-29 | 一种高鲁棒性dna测序用条形码生成和读取方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110060734A true CN110060734A (zh) | 2019-07-26 |
| CN110060734B CN110060734B (zh) | 2021-08-13 |
Family
ID=67318026
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910252268.7A Active CN110060734B (zh) | 2019-03-29 | 2019-03-29 | 一种高鲁棒性dna测序用条形码生成和读取方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110060734B (zh) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110929542A (zh) * | 2019-11-19 | 2020-03-27 | 天津大学 | 基于分组纠错码的测序条形码构造与软判决识别方法 |
| CN113345521A (zh) * | 2021-05-31 | 2021-09-03 | 天津大学 | 一种采用大片段dna存储的编码与恢复方法 |
| CN113343736A (zh) * | 2021-06-21 | 2021-09-03 | 天津大学合肥创新发展研究院 | 一种dna测序用条形码识别算法的硬件加速装置 |
| CN113380322A (zh) * | 2021-06-25 | 2021-09-10 | 倍生生物科技(深圳)有限公司 | 人工核酸序列水印编码系统、水印字符串及编码和解码方法 |
Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2002014553A2 (en) * | 2000-08-11 | 2002-02-21 | Favrille, Inc. | A molecular vector identification system |
| US20160046930A1 (en) * | 2013-02-20 | 2016-02-18 | The Johns Hopkins University | Methods of sequencing nucleic acids in mixtures and compositions related thereto |
| CN106504180A (zh) * | 2016-11-03 | 2017-03-15 | 大连大学 | 基于dna编码的图像加密纠错方法 |
| CN106788458A (zh) * | 2016-12-02 | 2017-05-31 | 天津大学 | 面向插入删节与替代错误的硬判决导向前后向估计方法 |
| CN106795568A (zh) * | 2014-10-10 | 2017-05-31 | 因维蒂公司 | 测序读段的de novo组装的方法、系统和过程 |
| WO2018114980A1 (en) * | 2016-12-19 | 2018-06-28 | Universiteit Gent | Polynucleotide shuffling method |
| CN108893464A (zh) * | 2018-07-13 | 2018-11-27 | 广州华银医学检验中心有限公司 | 一种免疫组库高通量测序文库的构建方法 |
| CN109416928A (zh) * | 2016-06-07 | 2019-03-01 | 伊路米纳有限公司 | 用于进行二级和/或三级处理的生物信息学系统、设备和方法 |
| CN109460822A (zh) * | 2018-11-19 | 2019-03-12 | 天津大学 | 基于dna的信息存储方法 |
-
2019
- 2019-03-29 CN CN201910252268.7A patent/CN110060734B/zh active Active
Patent Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2002014553A2 (en) * | 2000-08-11 | 2002-02-21 | Favrille, Inc. | A molecular vector identification system |
| US20160046930A1 (en) * | 2013-02-20 | 2016-02-18 | The Johns Hopkins University | Methods of sequencing nucleic acids in mixtures and compositions related thereto |
| CN106795568A (zh) * | 2014-10-10 | 2017-05-31 | 因维蒂公司 | 测序读段的de novo组装的方法、系统和过程 |
| CN109416928A (zh) * | 2016-06-07 | 2019-03-01 | 伊路米纳有限公司 | 用于进行二级和/或三级处理的生物信息学系统、设备和方法 |
| CN106504180A (zh) * | 2016-11-03 | 2017-03-15 | 大连大学 | 基于dna编码的图像加密纠错方法 |
| CN106788458A (zh) * | 2016-12-02 | 2017-05-31 | 天津大学 | 面向插入删节与替代错误的硬判决导向前后向估计方法 |
| WO2018114980A1 (en) * | 2016-12-19 | 2018-06-28 | Universiteit Gent | Polynucleotide shuffling method |
| CN108893464A (zh) * | 2018-07-13 | 2018-11-27 | 广州华银医学检验中心有限公司 | 一种免疫组库高通量测序文库的构建方法 |
| CN109460822A (zh) * | 2018-11-19 | 2019-03-12 | 天津大学 | 基于dna的信息存储方法 |
Non-Patent Citations (6)
| Title |
|---|
| JOAQUÍN EZPELETA等: "Designing robust watermark barcodes for multiplex long-read sequencing", 《BIOINFORMATICS》 * |
| JOHN A. HAWKINS等: "Indel-correcting DNA barcodes for highthroughput sequencing", 《PNAS》 * |
| TILO BUSCHMANN等: "Levenshtein error-correcting barcodes for multiplexed DNA sequencing", 《BUSCHMANN AND BYSTRYKH BMC BIOINFORMATICS》 * |
| 乔鑫: "蜘蛛牵引丝蛋白MaSp1表达载体构建及其在大肠杆菌中的表达", 《中国优秀硕士学位论文全文数据库_医药卫生科技辑》 * |
| 谢旻晟: "用于高通量DNA测序的编码方法研究", 《10.7666/D.Y2706686》 * |
| 陈为刚等: "基于DTMB中LDPC码与BCH码的乘积码构造", 《电视技术》 * |
Cited By (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110929542A (zh) * | 2019-11-19 | 2020-03-27 | 天津大学 | 基于分组纠错码的测序条形码构造与软判决识别方法 |
| CN110929542B (zh) * | 2019-11-19 | 2021-12-07 | 天津大学 | 基于分组纠错码的测序条形码构造与软判决识别方法 |
| CN113345521A (zh) * | 2021-05-31 | 2021-09-03 | 天津大学 | 一种采用大片段dna存储的编码与恢复方法 |
| CN113343736A (zh) * | 2021-06-21 | 2021-09-03 | 天津大学合肥创新发展研究院 | 一种dna测序用条形码识别算法的硬件加速装置 |
| CN113343736B (zh) * | 2021-06-21 | 2023-07-25 | 天津大学合肥创新发展研究院 | 一种dna测序用条形码识别算法的硬件加速装置 |
| CN113380322A (zh) * | 2021-06-25 | 2021-09-10 | 倍生生物科技(深圳)有限公司 | 人工核酸序列水印编码系统、水印字符串及编码和解码方法 |
| CN113380322B (zh) * | 2021-06-25 | 2023-10-24 | 倍生生物科技(深圳)有限公司 | 人工核酸序列水印编码系统、水印字符串及编码和解码方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110060734B (zh) | 2021-08-13 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110060734B (zh) | 一种高鲁棒性dna测序用条形码生成和读取方法 | |
| Kiah et al. | Codes for DNA sequence profiles | |
| Chandak et al. | Overcoming high nanopore basecaller error rates for DNA storage via basecaller-decoder integration and convolutional codes | |
| US11784666B2 (en) | Decoding signals by guessing noise | |
| KR940019106A (ko) | 신호 처리 방법 | |
| CN101288232B (zh) | 对数据进行编码和解码的方法以及设备 | |
| US20170109229A1 (en) | Data processing method and device for recovering valid code words from a corrupted code word sequence | |
| US6711711B2 (en) | Error correctible channel coding method | |
| CN111858507A (zh) | 基于dna的数据存储方法、解码方法、系统和装置 | |
| Shomorony et al. | Torn-paper coding | |
| CN112100982A (zh) | Dna存储方法、系统和存储介质 | |
| Chen et al. | Sequencing barcode construction and identification methods based on block error-correction codes | |
| Krishnan et al. | Barcodes for DNA sequencing with guaranteed error correction capability | |
| Goyal et al. | Sequence reconstruction problem for deletion channels: A complete asymptotic solution | |
| JP2008118327A (ja) | ビタビ復号方法 | |
| CN113687976A (zh) | 面向dna信息存储的编码和解码方法与装置 | |
| RU2327297C2 (ru) | Способ декодирования блоковых кодов со стираниями элементов | |
| Yin et al. | PREMIER—PRobabilistic error-correction using Markov inference in errored reads | |
| Shafir et al. | Sequence reconstruction under stutter noise in enzymatic DNA synthesis | |
| CN110929542B (zh) | 基于分组纠错码的测序条形码构造与软判决识别方法 | |
| Kracht et al. | Using the Davey-MacKay code construction for barcodes in DNA sequencing | |
| WO2022115159A1 (en) | Trellis based reconstruction algorithms and inner codes for dna data storage | |
| Shafir et al. | Sequence design and reconstruction under the repeat channel in enzymatic DNA synthesis | |
| Tang et al. | Correcting a substring edit error of bounded length | |
| Sabary et al. | The error probability of maximum-likelihood decoding over two deletion channels |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |