高通量测序技术及其在贝类转录组研究中的应用

刘敏 叶晟 杨凤 吴立新

摘要 高通量测序是传统DNA测序的一次重大技术创新，它为贝类养殖物种的研究带来了新的机遇。着重对454、Solexa 和 SOLiD高通量测序平台的技术原理、数据分析及其在贝类转录组的研究应用进行综述，同时对高通量测序未来发展方向做了总结和展望。

关键词 高通量测序；转录组；贝类；应用

中图分类号 S188 文献标识码 A 文章编号 0517-6611（2017）07-0130-04

High-throughput Sequencing Technology and Its Application in the Study of Shellfish Transcriptome

LIU Min，YE Sheng，YANG Feng，WU Li-xin*

（College of Fisheries and Life Science， Dalian Ocean University，Dalian，Liaoning 116023）

Abstract High-throughput sequencing is a major technology innovation of traditional DNA sequencing， which brings new opportunities to the research of shellfish species. This review focused on their technical principles， data analysis of the 454， Solexa and SOLiD high-throughput sequencing platforms and their application in the study of shellfish transcriptome. At the same time， the future development direction of high-throughput sequencing was summarized and forecasted.

Key words High-throughput sequencing；Transcriptome；Shellfish；Application

核酸测序技术是一种重要的分子生物学分析方法，被广泛地应用在生物学研究中。20世纪70年代，Sanger提出了具有里程碑意义的双脱氧核苷酸末端终止法[1]，该方法在“人类基因组计划”实施中发挥了巨大作用。随着测序技术的不断发展，第一代测序（Sanger测序）因其成本高、通量低和速度慢等缺点而不能满足研究者的需求[2]。自2000年以后，高通量测序技术平台相继出现，主要包括罗氏的454测序平台、Illumina的Solexa测序平台 和 ABI的SOLiD测序平台，它们以其通量高、速度快、准确度高、成本低等特点而被广泛应用在非模式生物测序研究中。高通量测序技术在贝类物种的转录组研究中发挥了重要作用，主要应用在贝类分子标记筛选、功能基因的挖掘和非编码miRNA方面的研究中。笔者综述了高通量测序技术的原理、数据分析及其在贝类转录组的研究应用，并展望了高通量测序今后的研究发展方向。

1 高通量测序的原理及测序步骤

高通量测序（High-throughput sequencing）又叫作深度测序或下一代测序，该技术是对Sanger测序的一次重大技术创新，在一次反应中能对数百万条的核酸序列进行测定。目前，高通量测序技术主要包括Roche/454[3]，Illumina/Solexa[4]和ABI/SOLiD 测序平台[5]，下面主要针对这3种测序技术的原理进行介绍。

1.1 Roche/454

454测序平台是采用边合成边测序的原理，其测序步骤如下：①DNA文库构建，将所测核酸序列打断成几百bp的小片段后，在两端分别加上特异性的接头。②乳液PCR，上述步骤带接头的DNA片段会与一个磁珠结合，磁珠被扩增试剂乳化之后，形成油包水的微反应器系统。测序PCR反应就发生在液滴包裹的微反应器里，每一个油滴体系中只包含1个DNA模板。扩增以后，每一个DNA分子都能进行大量富集，最后构成一个克隆集落。454测序仪的测序通道体积狭窄，仅仅能容纳一个磁珠。③转移至PTP 载板（Pico Titer Plate），将每个磁珠移至GS FLX测序仪载样板的每个小孔内，为测序做准备。④焦磷酸测序，在PTP孔内添加被激活液处理的磁珠，测序系统会把底物的聚合与荧光信号释放偶联起来。当发生碱基配对时，释放的焦磷酸基团在酶的作用下产生级联反应，获得测序荧光信号。通过捕捉的荧光信号，就能够达到核酸测序的目的。与其他高通量测序技术相比，454测序平台具有测序读长长和速度快的优点，适合基因組的重测序、从头测序和宏基因组学等。然而，454测序平台也存在一些缺点，如成本高、准确率低，对相同碱基聚合区域难以处理。454 测序平台在性能上很难有升级和突破，所以罗氏公司从2016 年起逐步停止焦磷酸测序仪的生产[6]。

1.2 Illumina /Solexa

Solexa测序平台的原理是边合成边测序，它是基于桥式PCR和可逆终止子的反应，测序步骤如下：①构建文库，首先把DNA序列随机打断成小片段，并在两端加上特定的接头。②锚定桥接，将第1步获得的DNA 片段变性为单链，并和接头引物结合，形成后续扩增的桥状结构。放大反应在带有8个通道 的flow cell玻璃管内进行，每一个通道的内表面有許多个固定的单链接头。③预扩增，添加Taq 酶和不带标记的底物dNTP进行桥式PCR 扩增。在变性时，单链桥型片段锚定到附近的固相表面。经过连续的扩增反应，上百万条成簇的双链DNA待测片段将会在flow cell 的固相表面上。④单碱基延伸测序，在flow cell玻璃管中添加DNA 聚合酶、4种不同的荧光标记dNTP和引物进行扩增，当互补链合成时，就会释放出带有荧光标记dNTP的荧光，系统根据捕捉的荧光信号就可以推断出待测的序列信息。

1.3 ABI/SOLiD

SOLiD测序平台是利用连接酶在连接过程中进行测序的原理，即边连接边测序。 SOLiD测序文库的构建和微乳滴 PCR 扩增与454技术类似，只是SOLiD的微珠更小，只有1 μm。它的特点是采用双碱基编码技术（two-base encoding）进行误差校正，其测序步骤如下：①制备单链DNA文库，将DNA待测序列打断成很小的片段，并在两头添加不同的接头，连接载体，获得单链DNA文库。②微乳化PCR（emulsion PCR），将带有接头的单链DNA固定在磁珠表面，进行 PCR 扩增，并对扩增产物进行 3′端修饰。③转移至上样玻片，将磁珠放置在SOLiD上样玻片（slide），准备测序。④进行双碱基编码测序，在体系中加入通用引物、DNA连接酶和带有 8 个碱基单链荧光探针混合物（3′-XXnnnzzz-5′），此探针3′端的第1位和第2位的碱基是确定的，5′端是带有荧光标记的（根据种类不同在6～8位zzz上加不同的荧光）。SOLiD测序包括5轮测序反应，每一轮反应又由多个连接反应构成。第1轮的第1次测序反应中，当探针3′端的第1和2位在连接酶的作用下进行配对连接时，就会发出荧光信号。记录荧光颜色信息后，除去6～8 位碱基和淬灭5′末端荧光信号，也就是说，实际上每一次测序的位置都相差5个碱基，第1轮得到第1和第2位的颜色信息。以此类推，第2次测得第6和7位的颜色信息，第3次测得第11和12位的颜色信息，第4次测得第16和17位的颜色信息……第1轮测序之后，将新合成的链变性、洗脱，引物重置，进行下一轮测序。由于之后的每一轮引物都要比前一轮引物向前移动1位，所以第2轮将会测得第0～1位、5～6位、10～11位和15～16位等的荧光颜色信息。5轮反应结束后，可以得到全部位置的荧光颜色信息，同时每个位点都会被检测2次，准确率达99.94%。SOLiD测序仪每次循环可以测75～110 bp，与Solexa技术类似，后续的序列拼接工作需要足够运算能力的计算机（主要指CPU和内存）和具有一定生物信息学操作经验的处理人员[7]。

454、SOliD和Solexa技术的测序原理、运行成本和数据输出量均不相同，但其共同之处是不需要克隆，可以直接进行高通量测序。此外，这3种测序技术也都要进行文库构建、片段扩增、测序、信号捕获等环节。研究者可以根据研究目的、试验经费、测序深度来选择测序平台。近年来，利用高通量测序对贝类进行测序研究刚刚兴起，加之传统的Sanger测序成本昂貴，使得贝类的基因组资源缺乏。对于分析贝类这样的非模式生物，需要找出邻近物种的基因组数据或较完整的转录组信息，若是没有搜寻到相关信息，在经费允许的情况下454测序可作为首要选择，从而为研究获得较长的序列片段。如果试验经费有限，Solexa测序也可以在被考虑范围之内，主要体现在2个方面：一方面是可以产生很大的数据量，从而增加测跨区域片段的概率；另一方面是功能强大的软件不断地被开发出来，可能使拼接序列的读长变得稍长些。从这两方面来看，Solexa测序在某种程度上也削减了短序列带来的缺陷。

2 高通量测序数据分析步骤

2.1 数据处理

测序获得的原始读段（raw reads），结果以fastq格式储存，并被提交到NCBI SAR数据库中。数据获得后要对测序读段进行统计，通过统计测序质量值分布和碱基分布情况，评估数据量是否满足后续分析要求。之后对原始数据进行质量评价，过滤掉低质量序列、接头序列和核糖体，获得高质量的数据进行接下来的组装分析。

2.2 数据组装

转录组测序序列组装包括有参基因组装和无参考基因组装（从头组装）。

有参基因组装：利用Scripture或Cufflinks软件，将测序读段比对到基因组上，再通过重叠信息得到graph，将graph中的信息进行转换，最后获得转录本。有参基因组装具有内存需求小、灵敏度高、污染小和拼接的转录本序列较完整等优点，但也存在一些缺点，如依赖参考序列和软件错误比对等。

无参考基因组装：又称从头组装，常用软件为Trinity。以Trinity软件对样本组装情况为例，根据序列之间的重叠信息组装获得contig，再通过序列的contig和双末端信息之间的相似性比对，然后进行contig聚类，而后在局部组装得到transcript[8]。通常把最主要的transcript当成Unigene。获取的Unigene进行后续的注释、定量、基因结构分析、基因的表达量差异和差异表达分析。

与有参基因组装相比，从头组装不依赖参考序列、比对软件，能较好地重建。但是，同样存在缺点，如要有较大内存和更高的测序深度、敏感性强、高相似度的转录本可能被合并等。

2.3 比对及注释

利用比對软件将组装得到的Unigene序列和公共数据库COG[9]、Swiss-Prot[10] 、 NR[11]、KEGG[12]进行相似性搜查，若是Unigene序列被比对到高优先级的数据库，就不用进行下一轮比对，否则仍会与之后的database比对完。在BLAST结果中，取比对到相似性最高的蛋白，最后得到这个Unigene的功能注释信息。将NR比对的结果输入到Blast2GO软件[13]中，就可以获得Unigene序列的GO注释信息。获得GO注释后，通过WEGO软件[14]可以统计GO功能分类，从而了解物种的基因功能。目前，采用的蛋白质数据库包括PDB、TrEMBL 和Swiss-Prot等。

3 高通量测序在贝类转录组研究中的应用

转录组是指细胞或特定组织在某一状态下转录出来的所有RNA集合。以真核生物基因组研究为例，转录组测序比全基因组测序针对性更强，仅仅研究被转录的部分，既减小了研究范围又降低了成本。此外，通过基因表达谱量的变化确定是否有某个特定基因，转录组是作为衡量功能基因组研究的一个重要指标。近些年，高通量测序技术被广泛应用在贝类转录组研究中，同时大量的贝类转录组研究文章被发表，取得了很大成果。

3.1 高通量测序在贝类分子标记开发中的应用

高通量测序技术作为一种快速开发分子标记的有效方法，越来越受到研究人员的青睐。通过高通量测序开发的分子标记主要包括简单重复序列（SSR）和单核苷酸多态性（SNP）。简单重复序列是指由2～6个核苷酸为重复基序的串联的DNA序列组成，具有高多态性、稳定性强和共显性遗传等特点，主要用于遗传分析、遗传图谱构建、亲缘关系鉴定和分子标记辅助育种。单核苷酸多态性是指由单个核苷酸变异所引起的DNA序列多态性，通常用于鉴别生物体的遗传变异。对SSR标记而言，一般采用MISA软件进行SSR标记筛查。MISA软件是一个可以大量识别和定位SSR序列的筛查工具，同时，它还提供了一个和引物设计Primer3的接口工具，利用这个工具可以将识别的SSR转成Primer3的格式，为设计引物提供了便利条件。Hou等[15]利用454 GS FLX测序系统对虾夷扇贝（Patinopecten yessoensis）进行转录组测序分析，检测到了超过49 000 SNPs和2 700 SSRs。Wang等[16]对栉孔扇贝（Chlamys farreri）进行高通量测序，并与之前研究的虾夷扇贝测序进行比较分析，38个假定的快速演变基因和大量的SNP被辨别。Niu等[17]在对缢蛏（Sinonovacula constricta）的转录组研究中也检测到大量的SSR和SNP。Franchini等[18]对中间鲍（Haliotis midae）进行了转录组测序，研究获得了数以千计的分子标记SSR（4 707个contigs鉴定出7 381）和SNP（4 380个contigs鉴定出11 934），同时为大量的SSR设计出高质量的PCR引物。Chen等[19]利用Illumina测序平台对河蚬（Corbicula fluminea）进行了转录组测序，鉴定出2 151条SSR序列。此外，还有其他贝类物种利用高通量测序技术去开发分子标记，如虾夷扇贝（Patinopecten yessoensis）[20]、大珠母贝（Pinctada maxima）[21]和褶纹冠蚌（Cristaria plicata）[22]等。从这些例子看出，利用454和Solexa测序技术可以大量挖掘分子标记SSR和SNP，它成为快速开发分子标记的一条有效途径。同时，这些分子标记非常有助于进一步分析群体遗传结构和构建遗传图谱。

3.2 高通量测序在挖掘贝类功能基因或新基因中的应用

根据试验方案：①取不同发育时期或经过胁迫（重金属、菌、温度或盐度胁迫等）处理的生物个体组织，然后提取RNA。②转录组测序后，利用组装的序列和数据库的序列进行相似性比对，研究者可以寻找到功能基因、新基因或差异表达基因等。

Moreira等[23]利用454测序技术对实施体外刺激的菲律宾蛤仔（Ruditapes philippinarum）进行转录组测序，产生了平均读长为250 bp的974 976个高质量读段，利用这些读段组装得到51 265个contigs，其中有22 915个contigs被注释。在注释的结果中，发现35个contigs包括大量的免疫相关基因，并且详尽地分析展示出几个免疫通路和进程的假定成员。同时，在补体途径中发现了双壳类从未被描述的分子序列。和454测序技术相比，尽管Solexa序列要比454所测读长短，但是它的应用性和测序深度要远大于454技术。Meng等[24]利用Illumina测序平台对镉暴露的虾夷扇贝（Mizuhopecten yessoensis）进行测序，14 240條序列被注释。通过GO注释和KEGG通路分析，辨别出720个基因涉及刺激反应和302个基因涉及免疫反应通路。此外，利用数字基因表达图谱DGE系统对鳃和消化腺检测转录组变化，在这2个组织分别检测到了7 556和3 002个差异表达基因。此外，Meng等[25]为了阐明铜暴露在免疫系统上的毒理學机制，以之前研究的虾夷扇贝（Mizuhopecten yessoensis）为材料，利用Illumina HiSeq 2000研究其鳃组织中的差异表达基因和转录本丰度。总共鉴定了1 312个上调基因和2 237个下调基因。同时，显著富集分析确定了9个GO术语和38个涉及铜暴露反应的途径。Zhang等[26]利用Illumina测序技术对美洲牡蛎（Crassostrea virginica）进行转录组分析揭示其先天免疫相关基因的丰富性。测序组装得到66 229个contigs，它涵盖了89.4%的已发表的EST和97.9%的线粒体基因，并且39 978个contigs和unigenes被辨别和注释。仍有26 251个contigs没有出现在比对的数据库中，说明Illumina测序可以胜任转录组深度测序的工作。

Qin等[27]利用454测序技术在8个不同的早期发育阶段对葡萄牙牡蛎（Crassostrea angulata）进行de novo转录组测序，提供了深度覆盖的EST数据库，并注释了大量的序列数据，同时鉴定和量化了6个独特序列涉及生长和早期发育的一些功能基因，例如肾上腺素能受体、多巴胺受体、表皮生长因子受体（EGFR）、胰岛素样生长因子1受体（IGFR），胰岛素诱导的蛋白质，卵泡抑素前体。Huan等[28]对文蛤（Meretrix meretrix）的4个不同幼虫期进行转录组测序，总共获得704 671个读数，并组装成124 737个独特序列（35 205个重叠群和89 532个单独序列）。进一步分析显示，这些序列中的94.66%（118 075个）是低表达水平的转录本。通过对UniProt蛋白质数据库进行搜索，15 000 215个（12.20%）序列被注释。通过分析每个独特序列的深度，进行初步定量分析，大量序列显示了4个幼虫阶段之间的转录差异，其与发育、生长、壳形成和免疫应答等相关。Song等[29]利用Hiseq 2500测序平台对脉红螺（Rapana venosa）的6个不同发育时期进行转录组测序和分析，对6个涉及神经分泌功能序列进行了量化和鉴定。对unigenes涉及的各种生物过程中的功能注释可以刺激对该物种早期发育机制的研究，同时了解早期发育和变态的机制将有利于该物种的水产养殖。

研究者利用高通量测序技术在对其他双壳贝类转录组研究中也发现了功能基因，如栉孔扇贝[16]、合浦珠母贝（Pinctada martensii）[30-31]、河蚬[19]、南极帽贝（Laternula elliptica）[32]、中国蛤蜊（Mactra chinensis）[33]和紫贻贝（Mytilus galloprovincialis）[34]等。

3.3 高通量测序在贝类MicroRNA发现和功能研究中的应用

微小RNA（miRNA）是一类长度为20～25个核苷酸（nt）的小非编码RNA，其在动物和植物中进行转录后的基因表达调节[35]。随着更多miRNA的调节靶基因和其功能被发现，miRNA指导的基因调节的广度和重要性开始备受关注。近年来，利用高通量测序有很多关于贝类的miRNA研究，其研究表明miRNA控制很多细胞功能，包括代谢、增殖、凋亡和免疫等，控制了不同生物代谢通路多个基因的表达，在贝类生长和发育中扮演重要角色。

对于小RNA（miRNA）测序而言，Illumina和SOLiD测序平台要优于454测序平台，因前两者有更大的单次数据产出量，可以提供一定的测序深度，从而提高其测序覆盖度。Bao等[36]对泥酣（Tegillarca granosa）血细胞miRNA进行测序，2个小RNA文库被构建（对照组和Cd处理组），经过与miRBase 数据库进行比对，总共鉴定出215个保守的miRNA，同时在没有泥酣参考基因组序列的情况下，利用生物分析发现39个泥酣特有的miRNA。在Cd刺激的文库中有5个miRNA表达显著上调，11个miRNA表达显著下调，上调miRNA主要涉及发育和细胞进程，下调miRNA在胚胎发育、炎症反应和免疫应答中发挥重要作用。Zhao 等[37]利用2种牡蛎物种香港牡蛎（Crassostrea hongkongensis）和长牡蛎（C.gigas）进行低盐度胁迫，并对鳃组织构建的4个小RNA文库进行高通量测序。在长牡蛎中鉴定出137个保守的和65个新的miRNA，同时在香港牡蛎中也鉴定85个保守的和2个新的miRNA。此外，在香港牡蛎物种中，显示相反表达类型的5个miRNA-mRNA互作对被辨别。对这些辨别的miRNA进一步研究，可能会揭示其在应答胁迫中的重要作用，从而为阐明牡蛎耐盐的分子机理提供线索。Jiao等[38]利用Illumina/Solexa深度测序对合浦珠母贝MicroRNA进行描述和辨别，研究鉴定出258个MicroRNA，其中，205个是跨物种保守的，53个是其特有的。一些pm-miRNA被预测参与生物矿化，例如miR-2305和miR-0046。Zhou等[39]利用Illumina平台对来自牡蛎血细胞的3个小RNA文库进行测序，以研究细菌攻击和热胁迫后miRNA的潜在免疫调节作用。对数据进行分析，鉴定出199个牡蛎miRNA（71种已知的和128种新型的）。研究显示免疫攻击可以诱导表达miRNA，其可以调节免疫应答，例如吞噬作用和凋亡。此外，一些免疫相关的表达miRNA也可以通过热应激来调节改善牡蛎对环境的适应。Chen等[40]也对刺激后的牡蛎血细胞miRNA进行了免疫调节研究。关于miRNA的研究还有中间鲍[41] 和长牡蛎[42]。这些在某阶段特异表达的miRNA可能会在一些生物进程中发挥重要作用，因此，miRNA的测序分析将会促进靶基因功能研究。

4 展望

近些年，高通量測序技术被广泛应用在贝类转录组的功能基因研究中。对于贝类这样的非模式生物，大多是无参考基因组序列的物种，通过高通量测序可以在短时间内产生海量数据，这些数据信息将会丰富贝类物种的遗传数据库，并且有利于这些物種进行深入的分子生物学研究。

高通量测序虽然在转录组学和分子生物学的研究中发挥了重要作用，但是在大量序列的拼接上仍然有不足之处。后续的拼接（序列的准确性和完整性）工作对于分子标记的开发和获取有意义的信息非常重要。在分子标记开发方面，错误拼接会导致标记扩增失败，同时太短序列的拼接也会增添设计引物的困难。研究者可以着手对拼接方法和测序工艺进行改进，从而增加序列的拼接长度。面对这些庞大的数据分析，如何得到有用的生物信息学信息将成为之后研究的热点。今后高通量测序技术会在各领域有更广泛的应用，为分子标记和功能基因的挖掘提供一条便利的途径，从而极大地促进遗传分析和分子育种方面的研究。

参考文献

[1] SANGER F，NICKLEN S，COULSON A R.DNA sequencing with chain-terminating inhibitors[J].Proceedings of the national academy of sciences，1977，74（12）：5463-5467.

[2] SHENDURE J，JI H.Next-generation DNA sequencing[J].Nature biotechnology，2008，26（10）：1135-1145.

[3] EMRICH S J，BARBAZUK W B，LI L，et al.Gene discovery and annotation using LCM-454 transcriptome sequencing[J].Genome research，2007，17（1）：69-73.

[4] MORIN R D，BAINBRIDGE M，FEJES A，et al.Profiling the HeLa S3 transcriptome using randomly primed cDNA and massively parallel short-read sequencing[J].Biotechniques，2008，45（1）：81-94.

[5] CLOONAN N，FORREST A R R，KOLLE G，et al.Stem cell transcriptome profiling via massive-scale mRNA sequencing[J].Nature methods，2008，5（7）：613-619.

[6] VAN DIJK E L，AUGER H，JASZCZYSZYN Y，et al.Ten years of next-generation sequencing technology[J].Trends in genetics，2014，30（9）：418-426.

[7] STRICKLER S R，BOMBARELY A，MUELLER L A.Designing a transcriptome next-generation sequencing project for a nonmodel plant species[J].American journal of botany，2012，99（2）：257-266.

[8] GRABHERR M G，HAAS B J，YASSOUR M，et al.Full-length transcriptome assembly from RNA-Seq data without a reference genome[J].Nature biotechnology，2011，29（7）：644-652.

[9] TATUSOV R L，GALPERIN M Y，NATALE D A，et al.The COG database：A tool for genome-scale analysis of protein functions and evolution[J].Nucleic acids research，2000，28（1）：33-36.

[10] APWEILER R，BAIROCH A，WU C H，et al.UniProt：The universal protein knowledgebase[J].Nucleic acids research，2004，32（S1）：115-119.

[11] 鄧泱泱，荔建琦，吴松锋.nr数据库分析及其本地化[J].計算机工程，2006，32（5）：71-74.

[12] KANEHISA M，GOTO S，KAWASHIMA S，et al.The KEGG resource for deciphering the genome[J].Nucleic acids research，2004，32（S1）：277-280.

[13] CONESA A，GTZ S，GARCA-GMEZ J M，et al.Blast2GO：A universal tool for annotation，visualization and analysis in functional genomics research[J].Bioinformatics，2005，21（18）：3674-3676.

[14] YE J，FANG L，ZHENG H K，et al.WEGO：A web tool for plotting GO annotations[J].Nucleic acids research，2006，34（S2）：293-297.

[15] HOU R，BAO Z M，WANG S，et al.Transcriptome sequencing and de novo analysis for Yesso scallop（Patinopecten yessoensis）using 454 GS FLX[J].PloS One，2011，6（6）：21560.

[16] WANG S，HOU R，BAO Z M，et al.Transcriptome sequencing of Zhikong scallop（Chlamys farreri）and comparative transcriptomic analysis with Yesso scallop（Patinopecten yessoensis）[J].PloS One，2013，8（5）：63927.

[17] NIU D H，WANG L，SUN F Y，et al.Development of molecular resources for an intertidal clam，Sinonovacula constricta，using 454 transcriptome sequencing[J].PloS One，2013，8（7）：67456.

[18] FRANCHINI P，VAN DER MERWE M，ROODT-WILDING R.Transcriptome characterization of the South African abalone Haliotis midae using sequencing-by-synthesis[J].BMC Research Notes，2011，4：59.

[19] CHEN H H，ZHA J L，LIANG X F，et al.Sequencing and de novo assembly of the Asian clam（Corbicula fluminea）transcriptome using the Illumina GAIIx method[J].PloS One，2013，8（11）：79516.

[20] DING J，ZHAO L，CHANG Y Q，et al.Transcriptome sequencing and characterization of Japanese scallop Patinopecten yessoensis from different shell color lines[J].PloS One，2015，10（2）：116406.

[21] DENG Y W，LEI Q N，TIAN Q L，et al.De novo assembly，gene annotation，and simple sequence repeat marker development using Illumina paired-end transcriptome sequences in the pearl oyster Pinctada maxima[J].Bioscience，biotechnology，and biochemistry，2014，78（10）：1685-1692.

[22] PATNAIK B B，WANG T H，KANG S W，et al.Sequencing，De Novo assembly，and annotation of the transcriptome of the endangered freshwater pearl bivalve，Cristaria plicata， provides novel insights into functional genes and marker discovery[J].PloS One，2016，11（2）：148622.

[23] MOREIRA R，BALSEIRO P，PLANAS J V，et al.Transcriptomics of in vitro immune-stimulated hemocytes from the Manila clam Ruditapes philippinarum using high-throughput sequencing[J].PloS One，2012，7（4）：35009.

[24] MENG X L，LIU M，JIANG K Y，et al.De novo characterization of Japanese scallop Mizuhopecten yessoensis transcriptome and analysis of its gene expression following cadmium exposure[J].PloS One，2013，8（5）：64485.

[25] MENG X L，TIAN X，LIU M，et al.The transcriptomic response to copper exposure by the gill tissue of Japanese scallops（Mizuhopecten yessoensis）using deep-sequencing technology[J].Fish & shellfish immunology，2014，38（2）：287-293.

[26] ZHANG L L，LI L，ZHU Y B，et al.Transcriptome analysis reveals a rich gene set related to innate immunity in the Eastern oyster（Crassostrea virginica）[J].Marine biotechnology，2014，16（1）：17-33.

[27] QIN J，HUANG Z X，CHEN J，et al.Sequencing and de novo analysis of Crassostrea angulata （Fujian oyster）from 8 different developing phases using 454 GSFlx[J].PloS One，2012，7（8）：43653.

[28] HUAN P，WANG H X，LIU B Z.Transcriptomic analysis of the clam Meretrix meretrix on different larval stages[J].Marine biotechnology，2011，14（1）：69-78.

[29] SONG H，YU Z L，SUN L N，et al.De novo transcriptome sequencing and analysis of Rapana venosa from six different developmental stages using Hi-seq 2500[J].Comparative biochemistry and physiology，Part D：Genomics and proteomics，2016，17：48-57.

[30] ZHAO X X，WANG Q S，JIAO Y，et al.Identification of genes potentially related to biomineralization and immunity by transcriptome analysis of pearl sac in pearl oyster Pinctada martensii[J].Marine biotechnology，2012，14（6）：730-739.

[31] SHI Y H，YU C C，GU Z F，et al.Characterization of the pearl oyster（Pinctada martensii） mantle transcriptome unravels biomineralization genes[J].Marine biotechnology，2013，15（2）：175-187.

[32] TRUEBANO M，BURNS G，THORNE M A S，et al.Transcriptional response to heat stress in the Antarctic bivalve Laternula elliptica[J].Journal of experimental marine biology and ecology，2010，391（1/2）：65-72.

[33] ZHANG J J，LI H J，QIN Y J，et al.Identification of functional genes involved in Cd2+ response of Chinese surf clam（Mactra chinensis）through transcriptome sequencing[J].Environmental toxicology and pharmacology，2016，41：113-120.

[34] GERDOL M，DE MORO G，MANFRIN C，et al.RNA sequencing and de novo assembly of the digestive gland transcriptome in Mytilus galloprovincialis fed with toxinogenic and non-toxic strains of Alexandrium minutum[J].BMC Research Notes，2014，7（1）：722.

[35] OBERNOSTERER G，LEUSCHNER P J F，ALENIUS M，et al.Post-transcriptional regulation of microRNA expression[J].RNA，2006，12（7）：1161-1167.

[36] BAO Y B，ZHANG L L，DONG Y H，et al.Identification and comparative analysis of the Tegillarca granosa haemocytes microRNA transcriptome in response to Cd using a deep sequencing approach[J] .PloS One，2014，9（4）：93619.

[37] ZHAO X L，YU H，KONG L F，et al.High throughput sequencing of small RNAs transcriptomes in two Crassostrea oysters identifies microRNAs involved in osmotic stress response[J].Scientific reports，2016，6：22687.

[38] JIAO Y，ZHENG Z，DU X D，et al.Identification and characterization of microRNAs in pearl oyster Pinctada martensii by Solexa deep sequencing[J].Marine biotechnology，2014，16（1）：54-62.

[39] ZHOU Z，WANG L L，SONG L S，et al.The identification and characteristics of immune-related microRNAs in haemocytes of oyster Crassostrea gigas[J].PloS One，2014，9（2）：88397.

[40] CHEN H，WANG L L，ZHOU Z，et al.The comprehensive immunomodulation of NeurimmiRs in haemocytes of oyster Crassostrea gigas after acetylcholine and norepinephrine stimulation[J].BMC Genomics，2015，16（1）：942.

[41] PICONE B，RHODE C，ROODT-WILDING R.Identification and characterization of miRNAs transcriptome in the South African abalone，Haliotis midae[J].Marine genomics，2017，31：9-12.

[42] XU F，WANG X T，FENG Y，et al.Identification of conserved and novel microRNAs in the pacific oyster Crassostrea gigas by deep sequencing[J].PloS One，2014，9（8）：104371.