CN110442472B - 一种dna数据存储混合错误纠正与数据恢复方法 - Google Patents

Abstract

本发明公开了一种DNA数据存储混合错误纠正与数据恢复方法，包括：根据双端读段重叠部分的编辑距离、标号部分的校验信息是否正确两个标准对测序读段进行筛选；根据已恢复的标号和文件号将读段进行分簇，将簇内读段分割为中间重合部分和不重合部分两部分；若中间重合部分或不重合部分的副本数量大于设定的门限值，采用聚类的方法确定中心序列，否则通过多序列合并确定；将每个读段对应数据部分的碱基分成若干个预设长度片段，每一个片段根据列序号的奇偶与前后对应的片段进行联合纠错；纠错采用多序列合并最终得到重复码长度为预设长度片段的可靠恢复。本发明主要解决测序读段中存在的插入/删节错误，以及采用低测序覆盖读段的合并。

Description

一种DNA数据存储混合错误纠正与数据恢复方法

技术领域

本发明涉及利用脱氧核糖核酸(DNA)的数据存储领域，尤其涉及一种DNA数据存储混合错误纠正与数据恢复方法。

背景技术

脱氧核糖核酸(DNA)是一种由脱氧核糖和四种含氮碱基(包括腺嘌呤A、胸腺嘧啶T、胞嘧啶C、鸟嘌呤G)组成的双链结构，是所有生命的遗传信息载体，控制着生命的发育、延续以及生命机能运作，是自然界天然的、最重要的信息存储载体。随着生物技术的发展，尤其是DNA合成与测序技术的发展，利用DNA序列作为数字化数据信息存储载体，已经具有技术可行性。DNA数字信息存储指的是把数字化信息存储于DNA的碱基序列之中，用不同的碱基或碱基组合表示数据。该技术利用寡核苷酸合成仪或高通量的芯片合成技术合成完全人工设计的包含数字信息的DNA序列来存储数据，利用DNA测序仪来读取所存储的信息。

DNA作为存储介质，与现有的磁带、硬盘、光盘以及固态Flash等存储媒介相比具有体积小、密度大，信息保持时间持久等特点。在密度方面，根据有关研究机构的分析，认为DNA存储的密度与磁盘、磁带以及固态存储等相对都有7个数量级的提升；在保持时间方面，如果采用合适的强化存放方法，DNA分子可以存储上百年，甚至上万年。DNA存储的缺点是：合成和读取DNA所存储的数字化信息成本较高，但DNA分子的日常存储成本相对便宜。值得期待的是，目前DNA合成与测序技术发展迅速，其发展速度甚至超过了半导体领域的摩尔定律，因此未来合成和读取用于存储的DNA序列的成本有望大幅度降低，成为未来重要的颠覆性存储方法。

2013年1月，欧洲分子生物学实验室的生物信息学研究所(EMBL-BMI)的NickGoldman教授团队在Nature发表的论文中将英文语言文本、扩展ASCII码对照表的文本、一篇pdf格式论文、JPEG图片和MP3格式的音频文件成功的存储并读取。在保证可靠性方面，该方案中采用了较为简单的重复编码方案。Goldman教授建议的读取方法的基本流程为：首先，根据同一双端测序数据中两条读段的中间重合部分的91个碱基之间的汉明距离和奇偶校验信息对测序读段进行筛选；然后，根据标号和文件号对测序读段进行分簇，并通过大数合并的方法得到各簇的中心序列；进一步，通过大数合并的方法得到多重覆盖区域的中心序列；最后，将所有片段进行拼接，并将拼接后的序列进行解码得到计算机多媒体文件。

该方法已经在中国申请专利，发明人为尼克·高曼，约翰·伯尼，名称为“DNA中数字信息的高容量存储”，申请时间为2013年5月31日，申请号为201380028511，本发明主要针对该项专利发明的采用重复码的DNA数据存储方法，提供一种优化的读取数据处理与纠错恢复方法。

DNA扩增过程中随机错误除替代错误外，还包括插入(Insertion)和删节(Deletion)错误。尼克·高曼等提出的读取方法在数据筛选时使用汉明距离作为筛选条件，可能会将仅发生少量插入或者删节错误的读段舍弃，减少了用于数据恢复的样本数量，浪费了合成与测序资源。另一方面，在数据量较少时，直接使用大数合并的方法恢复数据，会造成数据恢复的可靠性降低；如果考虑读段中的插入或者删节错误，该合并方法无法有效工作。该方案采用了四重重复码方法，重复码的合并也存在类似问题。本发明在数据筛选时以编辑距离作为筛选条件，提高了数据的利用率；同时，在数据恢复时，使用聚类或多序列合并的方法，对编辑错误进行了纠正，保证了数据恢复的可靠性。

发明内容

本发明提供了一种DNA数据存储混合错误纠正与数据恢复方法，本发明主要解决测序读段中存在的插入/删节错误，以及采用低测序覆盖读段的合并，详见下文描述：

一种DNA数据存储混合错误纠正与数据恢复方法，所述方法包括：

(1)对双端读取序列中的一条读段进行反转互补，获取两条有重叠的读段；

(2)根据双端读段重叠部分的编辑距离、标号部分的校验信息是否正确两个标准对测序读段进行筛选；

(3)根据已恢复的标号和文件号将读段进行分簇，将簇内读段分割为中间重合部分和不重合部分两部分；

(4)若中间重合部分或不重合部分的副本数量大于设定的门限值，使用聚类的方法确定中心序列，否则通过多序列合并确定；

(5)将每个读段对应数据部分的碱基分成若干个预设长度片段，每一个片段根据序列序号的奇偶与前后对应的片段进行联合纠错；纠错采用多序列合并最终得到重复码长度为预设长度片段的可靠恢复。

其中，所述步骤(2)具体为：

计算双端读段中间重合部分之间的编辑距离，保留编辑距离小于设定门限值的双端读段；

根据标号部分序列恢复标号并对其进行校验，保留校验正确的双端读段。

进一步地，所述根据标号部分序列恢复标号并对其进行校验，保留校验正确的双端读段具体为：

根据标号部分序列恢复标号；检验双端测序中标号部分的奇偶校验信息，若校验正确，保留当前双端读段，并记录该双端读段的标号及文件号；否则，丢弃该双端读段。

其中，所述步骤(4)具体为：

计算簇内不重合部分副本数量，若副本数量大于设定的门限值，则使用聚类的方法确定中心序列；否则，针对编辑错误，通过多序列合并确定中心序列；

计算簇内中间重合部分副本数量，若副本数量大于设定的门限值，则使用聚类的方法确定中心序列；否则，针对编辑错误，通过多序列合并得到中心序列；

将不重合部分和中间重合部分的中心序列合并，得到当前簇的中心序列。

进一步地，所述步骤(5)具体为：

(5.1)依次选取片段所对应的副本；

(5.2)判断当前副本所在读段对应的文件号的奇偶，若文件号为奇数，则将对应的副本进行反转互补；

(5.3)判断片段对应的副本是否全部选取完毕，若存在未选取的副本，则返回执行步骤(5.1)；否则，针对插入、删节与替代错误，通过多序列合并实现重复码片段的纠错恢复。

其中，所述多序列合并具体为：

(6.1)分别计算当前段中各副本与其它副本之间的编辑距离之和；将编辑距离之和最小的副本作为参考序列，其它副本作为待处理片段；

(6.2)将与参考序列之间编辑距离大于等于设定删除门限的待处理片段删除，分别计算各待处理片段与参考序列的编辑距离，若编辑距离小于设定的门限值，则保留该待处理片段；否则在当前段中删除；

(6.3)计算当前段中副本数量，若副本数量小于3，则将参考序列作为中心序列，流程结束；否则，执行步骤(6.4)；

(6.4)依次选取待处理片段，递归计算当前待处理片段与参考序列在各节点的编辑距离d(xⁱ,y^j)，并记录上一节点的位置pointer(i,j)；

(6.5)回溯，处理编辑错误，从节点(n,n)出发，通过利用pointer(i,j)中的信息，寻找一条到达节点(0,0)的最大似然路径并同步处理此路径上的插入、删节错误，若pointer(i,j)＝(i-1,j)，则表示待处理片段发生删节，从当前碱基连同后续碱基右移一个碱基位置，空出的位置填零补充；若pointer(i,j)＝(i,j-1)，则表示发生插入，将当前碱基在片段中删除，后续碱基左移一个碱基位置；由此，得到一条当前待处理片段的新的序列副本；

(6.6)判断待处理片段是否全部处理完毕，若存在未处理的片段，则返回执行步骤(6.4)；否则，将各待处理片段处理编辑错误后新生成的副本和参考序列进行大数合并得到中心序列。

本发明提供的技术方案的有益效果是：

(1)高测序数据利用率：本发明在进行数据筛选时根据双端测序文件中两条读段中间重合部分之间的编辑距离选择读段，与基于汉明距离选择读段的方法相比，用于文件恢复的读段数量获得提升，提高了测序数据利用率；

(2)高可靠性：本发明在纠错方面，根据段内副本数量的不同，采用了不同纠错机制。当副本数量较多时，使用聚类的方法可快速得到较可靠的中心序列；当副本数量较少时，采用多序列合并的方法得到中心序列，可保证序列恢复的可靠性，进而保证了整个文件恢复的可靠性，尤其是在低测序覆盖下有较明显的优势。

附图说明

图1为一种DNA数据存储混合错误纠正与数据恢复方法的流程图；

图2为本发明根据编辑距离对测序读段进行筛选的流程图；

图3为本发明对簇内片段进行纠错并得到中心序列的流程图；

图4为本发明采用多序列合并方法实现纠错的流程图；

图5为本发明中多序列合并的流程图；

图6为本发明中用于存储数据的DNA序列组成的流程图；

图7为本发明中不同数据覆盖深度下测序数据利用率的对比图；

图8为本发明中重复码译码后未恢复片段数量的对比图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面对本发明实施方式作进一步地详细描述。

参见图1，一种DNA数据存储混合错误纠正与数据恢复方法，该方法包括以下步骤：

(1)对测序数据进行预处理，具体是：针对双端读取的序列，将其中一条读段进行反转互补，从而得到两条有重叠的读段；

(2)根据双端读段重叠部分的编辑距离、标号部分的校验信息是否正确两个标准对测序读段进行筛选；

(3)根据已恢复的标号和文件号将读段进行分簇，考虑到标号部分已处理，并将簇内读段分割为中间重合部分和不重合部分两部分；

(4)若中间重合部分或不重合部分的副本数量Num大于设定的门限值copy_num，使用聚类的方法确定中心序列，否则通过多序列合并确定；

(5)考虑插入、删节与替代错误的重复码译码，具体是将每个读段对应数据部分的100个碱基分成四个25nt的片段，每一个片段根据序列序号的奇偶，按照正常的顺序或者反转互补的顺序，与前后对应的片段进行联合纠错；纠错采用任意多序列合并的方法或者本发明提出的多序列合并方法实现，从而得到重复码长度为25nt片段的可靠恢复。

所述根据双端读段重叠部分的编辑距离、标号部分的校验信息是否正确两个标准对测序读段进行筛选，如图2所示，具体步骤为：

计算双端读段中间重合部分之间的编辑距离，保留编辑距离小于设定门限值ed_filter的双端读段；

根据标号部分序列恢复标号并对其进行校验，保留校验正确的双端读段。

进一步地，所述根据标号部分序列恢复标号并对其进行校验，保留校验正确的双端读段具体为：

根据标号部分序列恢复标号；检验双端测序中标号部分的奇偶校验信息，若校验正确，保留当前双端读段，并记录该双端读段的标号及文件号；否则，丢弃该双端读段。

所述若中间重合部分或不重合部分的副本数量Num大于设定的门限值copy_num，使用聚类的方法确定中心序列，否则通过多序列合并确定，如图3所示，具体步骤为：

计算簇内不重合部分副本数量Num，若副本数量Num大于设定的门限值copy_num，则使用聚类的方法确定中心序列；否则，针对编辑错误，通过多序列合并确定中心序列；

计算簇内中间重合部分副本数量Num，若副本数量Num大于设定的门限值copy_num，则使用聚类的方法确定中心序列；否则，针对编辑错误，通过多序列合并得到中心序列；

将不重合部分和中间重合部分的中心序列合并，得到当前簇的中心序列。

所述考虑插入、删节与替代错误的重复码译码，具体是将每个读段对应数据部分的100个碱基分成四个25nt的片段，每一个片段根据序列序号的奇偶，按照正常的顺序或者反转互补的顺序，与前后对应的片段进行联合纠错；纠错采用任意多序列合并的方法或者本发明提出的多序列合并方法实现，从而得到重复码长度为25nt片段的可靠恢复，如图4所示，具体步骤为：

(5.1)依次选取片段所对应的副本；

(5.2)判断当前副本所在读段对应的文件号的奇偶，若文件号为奇数，则将对应的副本进行反转互补；

(5.3)判断片段对应的副本是否全部选取完毕，若存在未选取的副本，则返回执行步骤(5.1)；否则，针对插入、删节与替代错误，通过多序列合并的方法实现重复码片段的纠错恢复。

所述的简化的多序列合并的方法，如图5所示，具体步骤为：

令序列x＝(x₁,x₂,...,x_n)和y＝(y₁,y₂,...,y_n)分别表示长度为n的参考序列和待处理片段；令d(xⁱ,y^j)表示节点(i,j)的编辑距离(xⁱ和y^j分别表示序列x与y的前i个和前j个元素)；

(6.1)分别计算当前段中各副本与其它副本之间的编辑距离之和；将编辑距离之和最小的副本作为参考序列，其它副本作为待处理片段；

(6.2)将与参考序列之间编辑距离大于等于设定删除门限ed_erasure的待处理片段删除，分别计算各待处理片段与参考序列的编辑距离，若编辑距离小于设定的门限值ed_erasure，则保留该待处理片段；否则在当前段中删除；

(6.3)计算当前段中副本数量，若副本数量小于3，则将参考序列作为中心序列，流程结束；否则，执行步骤(6.4)；

(6.4)依次选取待处理片段，初始化当前待处理片段和参考序列的的编辑距离，即令d(x⁰,y^j)＝j，d(xⁱ,y⁰)＝i(0≤i,j≤n)，递归计算当前待处理片段与参考序列在各节点的编辑距离d(xⁱ,y^j)，并记录上一节点的位置pointer(i,j)，节点(i,j)处的编辑距离为：

其中，指示函数

当xⁱ＝y^j时，值为0；否则，值为1；上一节点位置为：

(6.5)回溯，处理编辑错误，从节点(n,n)出发，通过利用pointer(i,j)中的信息，寻找一条到达节点(0,0)的最大似然路径并同步处理此路径上的插入、删节错误，若pointer(i,j)＝(i-1,j)，则表示待处理片段发生删节，从当前碱基连同后续碱基右移一个碱基位置，空出的位置填零补充；若pointer(i,j)＝(i,j-1)，则表示发生插入，将当前碱基在片段中删除，后续碱基左移一个碱基位置；由此，得到一条当前待处理片段的新的序列副本；

(6.6)判断待处理片段是否全部处理完毕，若存在未处理的片段，则返回执行步骤(6.4)；否则，将各待处理片段处理编辑错误后的新副本和参考序列进行大数合并得到中心序列。

综上所述，本发明实施例根据双端测序数据中两条读段中间重合部分之间的编辑距离是否小于设定的门限值选择读段，提高测序数据利用率；然后将所有读段按照标号和文件号分簇后，在各簇内将读段分割为中间重合部分和不重合部分，并根据中间重合部分和不重合部分中副本数量是否大于设定的门限值，分别采用聚类或通过多序列合并的方法计算得到中心序列，进而组合得到簇中心序列。重复码的片段译码也采用多序列合并的方法。本发明纠错恢复方法，能对抗插入与删节错误，可提高测序数据利用率，实现采用少量测序数据的可靠恢复。

下面给出一个具体的实施例，说明本发明实施例给出的DNA数据存储的数据译码恢复方法的可行性。本实施例，使用引证文献使用数据，将本发明实施例中提出的混合错误纠错与数据恢复方法，与Goldman在论文与专利申请材料(引证文献)给出的恢复方法在不同数据覆盖深度下的性能进行了对比。其中，原始数据采用Goldman分析的数据，总共79,564,267组测序数据。该数据对应的原始DNA序列共153,305条，包括5个计算机多媒体文件。编码数据的具体格式如图6所示，每条DNA序列长度为117符号，除两端添加的1nt用于区分读取方向的标识信息外，该序列依序包含100nt的重叠编码数据和15nt的标引信息，其中标引信息分别包括2nt的文件编号、12nt的文件内编号以及1nt奇偶校验信息。原始数据中每条读段长度为104nt，在双端测序文件同一双端测序数据中，正向读取文件中读段的后91nt和反向读取文件中读段的前91nt为中间重合部分。为验证对比纠错与数据恢复方法，在双端测序文件中分别取75万、150万、300万和7900万条读段进行性能验证。

该实施例中，用于数据筛选的编辑距离预设门限值ed_filter设置为7，判断采用聚类或多序列合并的方法的副本数量门限值copy_num设置为9，处理编辑错误进行大数合并法中片段删除门限ed_erasure设置为7。具体中应该根据样本情况、测序方法与测序质量等参数动态选择上述参数，这里给出的参数仅针对Goldman提供的测序数据，其方法详见引证文献。

图7给出了不同数据覆盖深度下两种数据处理方法的测序数据利用率对比图，从图7中可以看出，本发明在进行数据筛选时根据双端测序文件中同一双端测序数据内的两条读段中间重合部分之间的编辑距离进行筛选，与基于汉明距离进行数据筛选的办法相比，可以利用的读段数量获得了提升，提高了测序数据利用率，对提高恢复质量，尤其是低测序覆盖度下的恢复质量有一定意义。

图8给出了重复码译码后未恢复片段对比图。从图8中可以看出，本发明实施例可保证中心序列恢复的可靠性，进而保证了整个文件恢复的可靠性，尤其是在低覆盖深度下，本方法有较明显的优势。为衡量整理的纠正混合错误与重复码纠正混合错误的整体性能，在仿真中，主要衡量碱基的错误率与缺失率，可以看出本方法在这两个指标上均取得了一定优势，尤其是采用的测序读段的数量为150万、300万时，有明显的改善，如果不考虑重复码的译码，仅仅考虑利用多副本的纠错，图8也给出了以25nt的段为单位衡量的错误率，可以看出，本方法均优于引证文献中的恢复方法，从而说明本发明，考虑存在插入与删节错误，并纠正插入、删节、替代混合错误的方法的有效性。

本领域技术人员可以理解附图只是一个优选实施例的示意图，上述本发明实施例序号仅仅为了描述，不代表实施例的优劣。

以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (6)

1.一种DNA数据存储混合错误纠正与数据恢复方法，其特征在于，所述方法包括：

步骤(1)对双端读取序列中的一条读段进行反转互补，获取两条有重叠的读段；

步骤(2)根据双端读段重叠部分的编辑距离、标号部分的校验信息是否正确两个标准对测序读段进行筛选；

步骤(3)根据已恢复的标号和文件号将读段进行分簇，将簇内读段分割为中间重合部分和不重合部分；

步骤(4)若中间重合部分或不重合部分的副本数量大于设定的门限值，使用聚类的方法确定中心序列，否则通过多序列合并确定；

步骤(5)将每个读段对应数据部分的碱基分成若干个预设长度片段，每一个片段根据序列序号的奇偶与前后对应的片段进行联合纠错；纠错采用多序列合并最终得到重复码长度为预设长度片段的可靠恢复。

2.根据权利要求1所述的一种DNA数据存储混合错误纠正与数据恢复方法，其特征在于，所述步骤(2)具体为：

计算双端读段中间重合部分之间的编辑距离，保留编辑距离小于设定门限值的双端读段；

根据标号部分序列恢复标号并对其进行校验，保留校验正确的双端读段。

3.根据权利要求2所述的一种DNA数据存储混合错误纠正与数据恢复方法，其特征在于，所述根据标号部分序列恢复标号并对其进行校验，保留校验正确的双端读段具体为：

根据标号部分序列恢复标号；检验双端测序中标号部分的奇偶校验信息，若校验正确，保留当前双端读段，并记录该双端读段的标号及文件号；否则，丢弃该双端读段。

4.根据权利要求1所述的一种DNA数据存储混合错误纠正与数据恢复方法，其特征在于，所述步骤(4)具体为：

计算簇内不重合部分副本数量，若副本数量大于设定的门限值，则使用聚类的方法确定中心序列；否则，针对编辑错误，通过多序列合并确定中心序列；

计算簇内中间重合部分副本数量，若副本数量大于设定的门限值，则使用聚类的方法确定中心序列；否则，针对编辑错误，通过多序列合并得到中心序列；

将不重合部分和中间重合部分的中心序列合并，得到当前簇的中心序列。

5.根据权利要求1所述的一种DNA数据存储混合错误纠正与数据恢复方法，其特征在于，所述步骤(5)具体为：

(5.1)依次选取片段所对应的副本；

(5.2)判断当前副本所在读段对应的文件号的奇偶，若文件号为奇数，则将对应的副本进行反转互补；

(5.3)判断片段对应的副本是否全部选取完毕，若存在未选取的副本，则返回执行步骤(5.1)；否则，针对插入、删节与替代错误，通过多序列合并实现重复码片段的纠错恢复。

6.根据权利要求1或4所述的一种DNA数据存储混合错误纠正与数据恢复方法，其特征在于，所述多序列合并具体为：

(6.1)分别计算当前段中各副本与其它副本之间的编辑距离之和；将编辑距离之和最小的副本作为参考序列，其它副本作为待处理片段；

(6.2)将与参考序列之间编辑距离大于等于设定删除门限的待处理片段删除，分别计算各待处理片段与参考序列的编辑距离，若编辑距离小于设定的门限值，则保留该待处理片段；否则在当前段中删除；

(6.3)计算当前段中副本数量，若副本数量小于3，则将参考序列作为中心序列，流程结束；否则，将参考序列作为中心序列，初始化其他序列与中心序列的编辑距离，执行步骤(6.4)；

(6.4)依次选取待处理片段，递归计算当前待处理片段与参考序列在各节点的编辑距离d(xⁱ,y^j)，并记录上一节点的位置pointer(i,j)；

(6.5)回溯，处理编辑错误，从节点(n,n)出发，通过利用pointer(i,j)中的信息，寻找一条到达节点(0,0)的最大似然路径并同步处理此路径上的插入、删节错误，若pointer(i,j)＝(i-1,j)，则表示待处理片段发生删节，从当前碱基连同后续碱基右移一个碱基位置，空出的位置填零补充；若pointer(i,j)＝(i,j-1)，则表示发生插入，将当前碱基在片段中删除，后续碱基左移一个碱基位置；由此，得到一条当前待处理片段的新的序列副本；

(6.6)判断待处理片段是否全部处理完毕，若存在未处理的片段，则返回执行步骤(6.4)；否则，将各待处理片段处理编辑错误后新生成的副本和参考序列进行大数合并得到中心序列。

Images