化学诱导后白木香转录组文库的构建与测序

吴宏清 王磊 陶美华 高晓霞 白玲 章卫民

白木香［Aquilaria sinensis（Lour.）Gilg］又称土沉香，瑞香科沉香属植物，是国产沉香的唯一植物资源［1］。当白木香树干受到物理、化学伤害或真菌侵染的情况下，可分泌出一种气味芬芳的防御性黑色树脂，即为我国传统名贵药材沉香。在自然界，沉香的形成需要几年至十几年、甚至数十年的时间。为了快速获得珍贵的沉香药材，人们通过各种方法人工造香。布兰切特和范贝克［2］利用亚硫酸氢钠、氯化钠、甲酸等化学物质诱导沉香属（Aquilaria spp.）植物产生沉香；Chen等［3］用氯化钠诱导成年白木香产生的沉香与天然沉香的化学成分极为相似。倍半萜类化合物是沉香的主要药效成分［4］，Kumeta和 Ito［5］及 Okudera 和 Ito［6］的研究结果表明，水杨酸和茉莉酸甲酯能诱导沉香属植物的悬浮细胞产生沉香倍半萜前体物质α-愈创木烯（α-guaiene）、α-蛇麻烯（α-humulene）和 δ-愈创木烯（δ-guaiene）。Xu等［7］从伤害诱导的白木香细胞中成功克隆到主产物为δ-愈创木烯的合成酶基因，但是从前体物质到沉香特征产物间的代谢途径仍然未知。

转录组测序是最近发展起来的利用深度测序进行转录组分析的技术［8］，目前已广泛应用于生物学、医学和临床研究及药物研发等。转录组测序结合表达谱分析的方法，可以在没有参考基因组的条件下对转录组进行全测序，以获得的转录组信息为参考，对不同样品的表达谱进行基因注释，通过比较获得不同表达谱间的差异表达基因，进而研究其基因功能。因此，利用转录组测序有助于白木香结香相关功能基因的发现，阐明沉香特征产物的代谢途径，揭示人工诱导白木香结香的分子机制。

本研究对5年树龄、化学诱导后1年的成年白木香植株进行总RNA的提取，获得的总RNA用于Illumina转录组测序，经测序文库的构建，上机测序，数据过滤，序列组装，旨在获得完整的白木香转录组信息，为后续的表达谱分析积累基础数据。

1 材料与方法

1.1 材料

1.1.1 试验材料 试验样品采自广东省信宜市珍稀沉香发展有限公司的白木香基地，参考王磊等［9］进行结香试验，其中用于转录组测序的白木香样品取自同一株5年树龄、甲酸处理后1年的已结香的白木香植株及另外一株未进行结香试验的白木香植株。样品采集时利用前端成钩状的刨刀，除去树皮后，在树干上从外到内依次刨下样品，外围未变色部分为白木样品（W样品），与白木相邻的一圈深棕色木材为结香样品（A样品），白木样品与结香样品间还有部分颜色为浅棕色的木质部，作为结香与未结香间的过渡样品（T样品），最内侧已腐烂的部分作为腐木样品（D样品），以及从另外一株未结香植株上采集的白木样品（C样品）。采集完样品后立即用锡箔纸包裹置于液氮中保存。

1.1.2 主要仪器与试剂 BioSpec-nano生命科学紫外/可见分光光度计；EPS 601电泳仪；GE ImageQuant 350凝胶成像系统；Hettich VNIVERSAL-32R台式冷冻离心机；Agilent 2100生物分析仪；Illumina HiSeqTM2000测序仪。焦碳酸二乙酯（DEPC）购自广州杰顺生物科技有限公司；改良异硫氰酸胍-CTAB提取液（38％水饱和酚，1mol/L异硫氰酸胍，2％ CTAB，100mmol/L NaAc-HAc pH5.2，2mol/L NaCl，2％ PVP），用前混匀；抽提液Ⅰ（水饱和酚∶氯仿∶异戊醇 = 25∶24∶1）；抽提液Ⅱ（氯仿∶异戊醇 = 24∶1）。

1.2 方法

1.2.1 白木香总RNA的提取 采用改良异硫氰酸胍-CTAB法分别提取白木香W、A、T、D和C样品总RNA：样品用液氮研磨后迅速分装到含有改良异硫氰酸胍-CTAB提取液的离心管中，剧烈震荡，室温静置5min；等体积的抽提液I抽提2次，离心取上清；等体积的抽提液II抽提1次，离心取上清；加入1/2体积的无水乙醇以及与上清等体积的4mol/L LiCl，颠倒混匀，-30℃静置过夜，离心取沉淀；沉淀溶于适量DEPC处理水中，加入1/10体积3mol/L NaAc-HAc，混匀后，加入3倍体积无水乙醇，-30℃静置30min，4℃离心10min；75％乙醇洗涤沉淀2次，溶于30μL DEPC处理水中低温保存。Agilent 2100生物分析仪对总RNA的RIN值及28S∶18S比值进行检测，验证RNA的完整性。等量合并各样品总RNA进入下一步操作。

1.2.2 白木香转录组测序文库的构建 用带有Oligo（dT）的磁珠富集总RNA样品中的mRNA，加入fragmentation buffer将mRNA打断成短片段，以打断后的mRNA为模板，用六碱基随机引物合成第一条cDNA链。加入缓冲液、dNTPs、RNase H和DNA polymerase I合成第二条cDNA链，然后用QiaQuick PCR纯化试剂盒纯化产物，用EB缓冲液洗脱后做末端修复、加A并连接测序接头，然后用琼脂糖凝胶电泳进行片段大小选择，最后进行PCR扩增，建好的测序文库用Illumina HiSeqTM2000进行测序。

1.2.3 Illumina HiSeqTM2000上机测序［11］使用 Illumina HiSeqTM2000测序平台进行转录组文库的测序。样品为白木香W、A、T、D和C各样品合并后的转录组测序文库。测序得到的原始图像经base calling转化为序列数据Raw reads。

1.2.4 数据过滤 对测序所得的Raw reads进行过滤，滤去的数据包括含接头的reads，N的比例大于5％的reads，重复的和质量数较低的reads（质量值Q≤10的碱基数占整个read的20％以上），过滤后所得为Clean reads，后续分析都基于此Clean reads。

1.2.5 序列的 De novo 组装［12］使用 Trinity［10］软件对Clean reads做De novo组装。将具有一定长度overlap的reads连成更长的片段Contig，然后与Clean reads重新比对，通过paired-end reads确定Contig所属的转录本以及在转录本中的分布，Trinity软件能将这些Contig连在一起，得到两端不能再延长的序列。然后使用Tgicl对其进行去冗余和进一步拼接，并对其进行同源转录本聚类，得到最终的Unigene。聚类后Unigene分为两部分，一部分是clusters（以CL开头），另一部分是singletons（以Unigene开头）。

2 结果

2.1 白木香总RNA的获得

使用改良异硫氰酸胍-CTAB法提取白木香各组织总RNA，经Agilent 2100生物分析仪检测，RIN值最小为6.7，28S∶18S均大于1.0，RNA总量远大于20μg，满足转录组测序的需求，结果如表1所示。

表1 用于转录组测序的白木香总RNA质量

2.2 测序产量统计

如表2所示，Illumina HiSeqTM2000上机测序后获得Raw reads共58804828条，过滤后获得Clean reads共54685634条，总测序长度为4921707060nt，Q20值达97.45％，测序质量较高。

2.3 组装结果统计及质量评估

如表3所示，经初步组装后，共获得190109条Contigs，平均长度有324nt，N50值为549，进一步组装后，共获得83467条Unigenes，平均长度高达702nt，N50值较高，达1120，序列组装理想，使得白木香的转录组数据得到较好的保存。

表2 白木香转录组测序统计

表3 组装结果

图1为所获得Contigs的组装统计结果。序列长度大于等于500nt的Contigs高达26786条，占总Contigs的14.09％；其中，大于等于1000nt的Contigs达2518条，占总Contigs的1.32％；大于等于2000nt的Contigs达1429条，占总Contigs的0.75％；大于等于3000nt的Contigs有1012条，占总Contigs的0.53％。

图2表示Contigs进一步组装后获得Unigenes的组装统计结果。序列长度大于等于1000nt的Unigenes高达17155条，占总Unigenes的20.56％；其中，大于等于2000nt的Unigenes有5189条，占总Unigenes的6.22％；大于等于3000nt的Unigenes有1691条，占总Unigenes的2.03％。

3 讨论

转录组测序技术是建立在新一代高通量测序平台（如Roche GS FLX或Illumina HiseqTM2000）上的cDNA测序技术，自2008年Nature和Science上分别发表利用转录组测序技术研究裂殖酵母（Schizosaccharomyces pombe）、酿酒酵母（Saccharomyces cerevisiae）转录组的论文［13，14］后，转录组测序技术已成为研究转录组的革命性工具。对比转录组学研究的其他方法，转录组测序可检测任意物种的转录组，无需前提信息，无需克隆，理论上可检测所有转录信息。

图1 白木香转录组Contigs组装统计

图2 白木香转录组Unigenes组装统计

为了研究白木香未结香组织和结香组织间的表达差异，获得与结香相关的功能基因，可以采用数字基因表达谱分析的方法，筛选未结香与结香组织间的差异表达基因。然而，如今NCBI上未有白木香的参考基因组或参考转录组信息，因此必须先进行转录组测序，获得总的转录组信息，在获得转录组序列的基础上进行表达谱分析及差异基因的研究。

张争等［15］利用454测序平台，对机械伤害后的白木香茎的转录组进行测序，共获得22095条平均长度为314nt的Unigenes。本研究对化学诱导后白木香样品的转录组测序、组装，共获得83467条平均长度为702nt的Unigenes，转录组信息保存较完整，为化学诱导白木香结香机理的研究提供大量的基础数据。为了便于后续的表达谱分析，防止因个体差异引起的无关信息过多，本研究选择来自同一株白木香不同部位的W、T、A和D 4个样品用于后续的分析，可消除因来自不同植株的样品对差异基因筛选的干扰。将W样品作为对照组，T、A和D 3个样品分别作为试验组，获取各对照的差异表达基因，研究其基因功能，进而揭示化学诱导白木香结香的分子机理，获得与化学诱导相关的基因和转录因子。

4 结论

采用改良异硫氰酸胍-CTAB法提取白木香各组织总RNA，构建转录组测序文库后进行Illumina HiSeqTM2000上机测序，共获得54685634条Clean reads，总长度为4921707060nt，经多次组装获得83467条Unigenes，平均长度为702nt，N50值为1120，大于等于3000nt的Unigenes有1691条，占总Unigenes的2.03％，测序和组装质量较高。

［1］中国科学院中国植物志编辑委员会. 中国植物志［M］. 北京：科学出版社, 1999：290.

［2］布兰切特RA, 范贝克HH. 栽培的沉香木：US, ZL 02810500.1［P］. 2004-10-27.

［3］Chen HQ, Yang Y, Xue J, et al. Comparison of compositions and antimicrobial activities of essential oils from chemically stimulated agarwood, wild agarwood and healthyAquilaria sinensis（Lour.）Gilg trees［J］. Molecules, 2011, 16（6）：4884-4896.

［4］杨俊山. 沉香化学成分的研究概况［J］. 天然产物研究与开发,1998, 10（1）：99-103.

［5］Kumeta Y, Ito M. Characterization of δ-guaiene synthases from cultured cells ofAquilaria,responsible for the formation of the sesquiterpenes in agarwood［J］. Plant Physiol, 2010, 154（4）：1998-2007.

［6］Okudera Y, Ito M. Production of agarwood fragrant constituents inAquilaria calliand cell suspension［J］. Plant Biotechnology, 2009,26（3）：307-315.

［7］Xu YH, Zhang Z, Wang MX, et al. Identification of genes related to agarwood formation：transcriptome analysis of healthy and wounded tissues ofAquilaria sinensis［J］. BMC Genomics, 2013, 14（1）：227.

［8］祁云霞, 刘永斌, 荣威恒. 转录组研究新技术：RNA-Seq及其应用［J］. 遗传 , 2011, 33（11）：1191-1202.

［9］王磊, 章卫民, 高晓霞, 等. 一种人工诱导白木香产生沉香的方法 ：中国，CN102302041A［P］. 2012-01-04.

［10］Grabherr MG, Haas BJ, Yassour M, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome［J］.Nat Biotechnol, 2011, 29（7）：644-652.

［11］Feng C, Chen M, Xu CJ, et al. Transcriptomic analysis of Chinese bayberry（Myrica rubra）fruit development and ripening using RNA-Seq［J］. BMC Genomics, 2012, 13：19-33.

［12］Xiang LX, He D, Dong WR, et al. Deep sequencing-based transcriptome profiling analysis of bacteria-challengedLateolabrax japonicusreveals insight into the immune relevant genes in marine fish［J］. BMC Genomics, 2010, 11：472-492.

［13］Wilhelm BT, Marguerat S, Watt S, et al. Dynamic repertoire of a eukaryotic transcriptome surveyed at single-nucleotide resolution［J］. Nature, 2008, 453（7199）：1239-1243.

［14］NagalakshmiU, Wang Z, Waern K, et al. The transcriptional landscape of the yeast genome defined by RNA sequencing［J］.Science, 2008, 320（5881）：1344-1349.

［15］张争, 高志晖, 魏建和, 等. 三年生白木香机械伤害转录组学研究［J］. 药学学报 , 2012, 47（8）：1106-1110.