1 引言

在当今这个信息大爆炸的时代,全世界在近两年内产生的信息量比过去五年的信息总量还多,数字信息正以惊人的速度增长和积累,根据国际权威机构Statista的统计,全球数据量在2019年多达41 ZB。现阶段人们使用的存储设备,如磁盘和半导体等逐渐暴露出先天不足的缺点^[1-3],因此寻找新一代存储技术刻不容缓。

脱氧核糖核酸(DNA)是一种天然的信息载体,具有容量大、存储密集、并行存取、无磨损和寿命长等先天优势^[4]。此外,随着DNA合成及测序技术的快速发展,研究人员把新一代数据存储介质的目光投向DNA,利用DNA中A(腺嘌呤)、G(鸟嘌呤)、C(胞嘧啶)和T(胸腺嘧啶) 4个碱基对二进制数据信息进行编码,并结合DNA人工合成技术来存储文本文档、图片、音频和视频等数据信息。

DNA信息存储技术作为信息领域和合成生物学领域的一项交叉融合技术,其开辟了一种新的高效存储模式,对节约存储能源消耗及推进大数据存储发展有着重要的影响及作用。

早在上世纪70年代,就有国外学者提出了利用DNA的不同状态来表示信息的想法。近年来,随着微软、欧洲生物信息学研究所、哥伦比亚大学、华盛顿大学、哈佛大学和英国剑桥顾问公司等众多科研机构和公司的加入^[5-8],DNA的存储性能在信息存储领域中正飞速发展。2017年,哈佛大学医学院的Shipman等^[9]利用了DNA剪切技术中的CRISPR-Cas基因编辑系统,将带有黑白图像和短视频影像编码的DNA序列导入一群活大肠杆菌的基因组,以此实现信息的快速复制;2018年,Organick等^[10]编码和存储了超过1300万个DNA寡核苷酸中的35个不同文件,且可以使用随机访问的方法独立、无误地恢复每个文件;2018年6月美国五角大楼对外宣称,借助DNA全新存储策略来保护公民的大量敏感信息,美国高级情报研究计划署打造了一种能够在DNA等人造分子结构上写入数据的桌面设备。

但国内对该领域的研究尚处于起步阶段,目前华中科技大学、天津大学及国防科技大学均有研究团队在进行相关研究,我国对于DNA数据存储的系统研究在逐步开展,并在2018年度“合成生物学”重点专项中将DNA数据存储技术作为其中一个子项。

本文研究一种基于Raptor码的DNA信息存储技术。该技术拥有抗数据损坏的稳健性,能够以恒定不变的成本对信息进行编解码,保证所得编码符号的平均度值恒定,降低编码的复杂度,缩短编码时间,提高DNA信息存储的性能。

图 1. DNA-Raptor码技术流程图

Fig. 1. Flow chart of DNA-Raptor code technology

下载图片 查看所有图片

2 方法设计与原理

基于Raptor码的DNA信息存储技术流程如图1所示。先将需要存储的信息文件利用二进制转换器转为二进制码流,对二进制码流采用Raptor码进行编码使其含有A、T、C、G的碱基序列,在序列两端加入合成DNA所需的引物段(一段用于促进DNA合成的碱基序列),再将这些碱基序列使用生物技术合成DNA链,即完成信息的存储。当读取DNA中的信息时,先使用聚合酶链反应(PCR)技术对DNA链进行扩增复制,获得多段相同的DNA副本链,再对DNA副本链使用DNA测试技术将其转换为碱基序列,对碱基序列进行解码以恢复成二进制文件,将二进制文档输入转换器中即可获得所存储的信息。

整个信息的存储技术流程中最关键的部分为信息编码部分,该部分主要由编码、纠错和筛选三个环节构成。

2.1 编码

Raptor码是一种以LDPC(Low-Density Parity-Check)为内码,LT(Luby Transform)为外码的编解码算法。先使用LDPC对信息源符号进行编码以获得中间符号,再对中间符号进行LT编码以生成编码符号。由于预编码环节的加入,确保Raptor码能够在较低的编解码复杂度下,仍具备良好的解码性能。具体的编码步骤如下^[11]。

step1:生成校验矩阵G。假设ε_pre与δ_pre分别为预编码的译码开销与译码失败概率的上限,则预编码LDPC的码率为

R=1+εpre/21+εpre。(1)

根据ε_pre=0.1,得到编码的码率R=0.95。G选择的度值为4,输入的信息符号数k为240,输出的中间符号数L为256,实验采用的是(16,4,4)LDPC码,其中16表示LDPC码长,4表示LDPC码的每一列非零元素的个数,4表示LDPC码的每一行非零元素的个数。

step2:将信息符号D(x₁,x₂,…,x_k)与LDPC码的校验矩阵G^-1相乘,获得中间符号矩阵C= [c0,c1,…,cL-1]T,即

C=G-1D。(2)

step3:采用概率转移法对度分布函数进行改进,构建一个适用于短码长弱LT码的度分布函数μ(d),表达式

μ(d)=0.04d+0.495d2+0.167d3+0.08d4+0.07d5+0.039d8+0.025d9+0.035d19+0.049d25,(3)

式中:d为度值。

step4:对预编码获得的中间符号C进行LT编码。先根据μ(d)随机选择一个度值d(1≤d≤k),再从L个中间符号中随机选择d个符号,对这d个符号进行模二异或运算以生成一个编码符号y_i,重复step4即可产生无数个编码符号(y₁,y₂,…,y_i)。

接收端收到编码符号后,并不需要获得全部中间符号,仅需解码一定个数的中间符号,再利用LDPC的性质即可恢复原始符号信息。

DNA存储技术中的信息需经过传播、DNA合成、PCR复制和DNA测序等众多过程。当信息在信道中传输时,由于白噪声等一些信道噪声较大极易引入误差,为了确保DNA信息的存储质量,加入的纠错码起着至关重要的作用。RS(Reed-Solomon)纠错码由于其性能良好,目前被越来越多的研究人员应用在DNA信息存储技术中,以保障信息的存储质量。实验在RS纠错码的基础上,结合DNA自身的结构特点,提出一种四进制RS纠错码机制^[12],编码流程如图2所示。

图 2. 四进制纠错编码流程

Fig. 2. Flowchart of quaternary error correction coding process

下载图片 查看所有图片

设计RS纠错码的关键在于,需确定伽罗华域元素表及生成多项式g(x),其中x为对应的信息符号。四进制RS纠错码的伽罗华域GF[(2²)²]由GF(2²)基于生成多项式g(x)=Z²+Z+2形成,其中Z为对二进制RS纠错码转换成四进制RS纠错码的信息符号。四元复合域GF[(2²)²]与二元扩展域GF(2⁴)本质上为GF(16)的伽罗华域。

(n,k)RS纠错码的编码算法通过添加t=n-k个冗余符号,将k个信息符号扩展为n个,其中t为冗余符号,n为冗余符号和输入信息符号组成的总的符号。从信息空间到编码空间的映射满足

v(x)=r(x)+xn-km(x),(4)r(x)=xn-km(x)mod[g(x)],(5)

式中:v(x)为编码序列;r(x)为余数;m(x)为信息序列;g(x)为n-k次生成多项式。

以(15,13)RS纠错码为例,用来说明四进制RS纠错码的纠错过程。若n=15,k=13,则t=n-k=2,即每13位信息能纠正1处错误。已知生成多项式g(x)的根为-1和-β,可得

g(x)=(x+1)(x+β)=x2+(1+β)x+β。(6)

而(1+β)映射至四元复合域中为β⁴,所以g(x)可以表示为

g(x)=x2+β4x+β。(7)

因RS码具有线性特性,则余数多项式可写为

r(x)=m0{xn-kmod[g(x)]}+m1{xn-k+1mod[g(x)]}+…+mk-1{xn-1mod[g(x)]},(8)

式中:x^n-k+amod [g(x)],∀a∈[0,k-1]为给定的n,k及g(x)计算获得的常量;m₀、m₁、m_k-1为对应项的系数。由(8)式可知,信息序列中任何项发生变化都会独立影响最终余数,因此,可单独计算每个信息符号的余数,再将这些余数相加可获得整段信息序列的余数。四进制RS纠错码的编码步骤如下。

step1:将四进制信息序列进行两两分组以转换成β^j形式,j∈{0,k-1},再对转换获得的序列重新分段,确保每段含k个符号。

step2:计算x^n-k+amod [g(x)]的值。由于n-k=2,故获得的值为C₂x+C₁,其中C₁、C₂为一元方程的各项系数。

step3:根据(8)式,将step1中计算的值与信息多项式的对应符号相乘。

step4:将每个信息符号的余数进行相加,获得整段信息的余数。再代入(4)式计算v(x)。

step5:对每段信息的多项式执行step3和step4,直至编码完成,再将结果序列转换为四进制序列。

改进后,四进制RS纠错码对于每26个碱基能纠正1个错误,不仅与DNA的结构特点相匹配,也提高信息存储的正确性。

2.2 筛选

四进制RS纠错码主要保障信道传输的准确无误性。此外,DNA的信息存储过程与经典的信息理论信道过程有所不同,DNA序列中的错误模式依赖于输入的碱基序列。因此,对DNA序列中常见的错误模式进行研究,发现均聚物及鸟嘌呤和胞嘧啶(GC)含量是DNA合成和测序出错的主要因素^[13-14],所以除了引入RS纠错码进行信道纠错外,还需对碱基序列进行GC含量及均聚物的筛选,降低合成和测序的难度及错误率。

为了定量分析,将DNA序列分为两类:有效类和无效类。当GC含量在40%~60%之间,且其最大均聚物长度不超过4个核苷酸时,该序列将被认为有效,反之视为无效丢弃。

实验提出了两种筛选方案,如图3所示。方案A:先对信息码流进行GC含量及均聚物的筛选,筛选后再将该码流转换为碱基序列并输出;方案B:先将码流转为碱基序列后,再对碱基序列进行筛选。

图 3. 筛选方案的流程

Fig. 3. Screening plan process

下载图片 查看所有图片

2.3 性能指标

选用编码效率及编码时间作为衡量标准,用以衡量DNA信息存储框架的优劣。其中编码效率E定义为每个碱基上可存储的信息比特数,表达式为

E=FO,(9)

式中:F为输入文件的比特数;O为输出的碱基个数。由于DNA序列共含有4个碱基,若将二进制码流按照四进制转换为碱基序列,则理论上每个碱基的编码效率达到2 bit。编码时间定义为整个编码过程从开始至结束的CPU时间。

3 实验结果及分析

实验过程:先对所提的四进制RS纠错方案进行仿真。输入符号的大小选为1×10⁵的二进制码流矩阵,信道为高斯白噪声,信噪比设置为10 dB,对加入RS纠错码后的编码符号分别引入0.02%,0.04%,0.10%,0.20%的出错率,再对其进行解码纠错,统计解码纠错后仍存在的错误符号个数,计算纠错后的错码概率。为了避免结果的偶然性,每组实验需重复进行100次,最终记录每组的平均值,实验结果如表1所示。从表1可以看到,错码概率越大,引入的错误也越多,纠错后出现错码的概率也越大,但对整体而言,四进制RS码能够解决大多数情况的出错问题。

接着对所提的两种筛选方案进行性能对比。对大小为9 KB的英文文件使用Raptor码进行编码后以获得编码文件,对编码文件分别使用两种筛选方案进行筛选后以转换生成DNA碱基序列。记录两种方案的编码效率及编码时间,所得结果如表2所示。从表2可以看到,两种编码方案的编码效率一致,但方案A的编码时间略少于方案B,这是因为方案A经过筛选后再进行碱基转换,比方案B少一环节,故后续分析均采用方案A。

图4为DNA-Raptor的编码流程,其中GGCA为编码获得的碱基序列。信息符号经Raptor编码后转换为四进制序列,对该序列加入四进制RS纠错码后再进行GC含量及均聚物的筛选,若通过筛选则转换为碱基序列,反之丢弃,进入下一轮编码。

为了评估所提的DNA-Raptor存储方案的性能,将文本、音频和图片共三种格式的数据文件分别作为输入文件输入至DNA-Raptor存储框架中,编码转换为碱基序列以实现存储。表3为不同文件格式的存储性能参数。

图 4. DNA-Raptor码的编码过程

Fig. 4. Coding process of DNA-Raptor code

下载图片 查看所有图片

为了验证所提方案的可行性,对DNA链进行生物实验。选择了9 KB的英文文本(Jane Austen的《Pride and Prejudice》部分节选)作为实验的测试文件,使用Raptor码对其进行编码,转换成碱基序列后再进行DNA合成,共获得18条DNA链。

生物实验主要包含三个步骤:根据获得的碱基序列合成DNA链,对每段存储信息的DNA链进行PCR扩增复制后获得副本并对所得副本的DNA链进行测序。合成DNA前需先使用相关软件设计各段碱基序列的引物,再通过菌落PCR获得靶序列,完成后通过酶鉴定检测来扩增结果,酶样品的电泳图显示获得的序列与预期序列相同。采用DNA测序技术中的Sanger测序法对DNA链进行测序,从而获得测序结果。最后对测序结果进行解码转换,成功恢复存储内容。

为了全方位评价所提的DNA-Raptor信息存储方案的性能,将所提方案与近几年主流的DNA信息存储方案^[4-7,15-16]进行对比,性能参数如表4所示。从表4可以看到,基于Raptor码的DNA存储框架的每个碱基的平均编码效率为1.49 bit,所提方案引入存储信息中的冗余占比较低,编码效率普遍优于大部分其他方案,仅略低于文献[ 7]提出的DNA喷泉码方案,故更详细的对这两种方法进行实验对比。

实验主要利用Raptor码对信息进行存储,而文献[ 7]利用的是LT码,所以先对比这两种编码的方式。分别使用LT码与Raptor码编码同一个图片文件,统计所得编码信息的度值,获得的度值-频数曲线如图5(a)所示。从图5(a)可以看到,LT码的度值较分散,并且由于鲁棒孤波分布的影响,在d=17处出现了一个波峰,最大度值为46,整个编码复杂度较高;Raptor码的度值大部分集中在5以内,最大度值为20,峰值消失,小度值的出现使频数明显增加,降低编码的复杂度,提高解码的成功率。

图 5. LT码与Raptor码的编解码性能。(a)度值与频数的变化曲线;(b)冗余与解码成功率的变化曲线

Fig. 5. Encoding and decoding performance of LT code and Raptor code. (a) Change curves of degree value and frequency; (b) change curves of redundancy and decoding success rate

下载图片 查看所有图片

图5(b)为译码成功率的变化曲线,其中横坐标为LT编码中的信息冗余量ε,不包括预编码引入的其他冗余。从图5(b)可以看到,解码成功率随着冗余量的增加而上升。此外Raptor码的解码性能明显优于LT码,当冗余量为0.06时,Raptor码的成功解码率已达到0.90。

与LT码最大不同,Raptor码能够保证编码符号的平均度值恒定,从而降低编解码的复杂度,提高存储性能,缩短时间。分别将大小为9 KB的文本文件作为DNA喷泉码及DNA Raptor码两个存储方案的输入文件,设置ε均为0.05,对LT编码过程中的度值及所选信息符号段数进行统计。DNA喷泉码的平均度值为4.38,DNA-Raptor码由于具有预编码环节,故当利用LT编码时可选用度数更低的弱LT码进行编码,平均度值为3.00。图6为两种方案中每个信息符号参与编码的频率统计。从图6可以看到,每个符号均参与编码,并且度值的出现频率集中在0.01~0.055之间,各信息符号编码的覆盖率良好。

图 6. 两种方案的信息符号频率。(a) DNA-Raptor码;(b) DNA喷泉码

Fig. 6. Information symbol frequency of two plans. (a) DNA-Raptor code; (b) DNA fountain code

下载图片 查看所有图片

接着对两种方案的编码性能进行统计,结果如表5所示。从表5可以看到,DNA-Raptor码的编码效率略低于DNA喷泉码,这是因为Raptor码在预编码环节引入部分信息的冗余。对于编码时间来看,DNA-Raptor码明显缩短编码时间,编码的复杂度更低。

此外,在DNA喷泉码中将平均度值降至2.73,获得的信息符号频率如图7所示。从图7可以看到,第11个信息符号的频率为0,说明信息符号未编码至编码符号中,后期解码将会失败。DNA-Raptor码解码时不需要恢复全部中间符号,仅需解码一定比例的中间符号再结合预编码的性质即可解码以获取源文件,降低解码的失败概率。

图 7. DNA喷泉码的信息符号频率

Fig. 7. Frequency of information symbols in DNA fountain code

下载图片 查看所有图片

与DNA喷泉码相比,所提的DNA-Raptor码信息存储框架最大的优势在于编解码复杂度低,能够缩短信息编解码的时间,提高解码的成功率。此外就纠错机制而言,DNA喷泉码中每528个碱基只能纠正8个错误,而所提的四进制RS纠错码每528个碱基能纠正21个错误,极大地保证了信息传输的准确性。

4 结论

提出一种基于Raptor码的DNA信息存储方案,将Raptor码中可生成无数编码符号及编解码复杂度低的特性应用于DNA信息存储领域,克服其他编码方式易丢失信息、可延展性差和编解码复杂度高等缺点。此外针对DNA含有4个碱基的特点,提出了四进制RS纠错机制,提高编码效率的同时,也提高信息的有效性及抗干扰性,并且针对GC含量及均聚物对DNA信息存储的影响提出筛选方案。最后使用三种不同格式的文件验证所提的DNA-Raptor存储方案的编码效率及编码时间的性能,并进行生物实验验证整个信息存储过程的可行性,同时也将所提存储方案与DNA喷泉码存储方案进行对比。实验结果表明:所提方案的编码时间短,复杂度低,适用于存储各种格式的文件。