RNA-Seq定量分析盐肤木对铅胁迫的响应

2019-07-19夏丽丹胡华英周垂帆

农业环境科学学报 2019年7期

夏丽丹，张虹，胡华英，曹升，周垂帆 *

（1.福建农林大学林学院，福州 350002；2.海峡两岸水土保持协同创新中心，福州 350002；3.福建长汀红壤丘陵生态系统国家定位观测研究站，福州 350002）

铅是一种有毒的重金属，采矿、印染、电镀等行业排放的废弃物中含有大量的铅，这使得土壤环境污染严重，且铅具有强累积性，可通过食物链进入人体，与人体中的多种酶发生反应，损伤人的免疫、消化、神经与生殖系统，从而威胁人体健康[1]。

铅不是植物生长的必需元素，过量的铅会对植物产生毒害效应：叶片小、发育不良、生物量减少、叶绿素含量降低、光合作用下降等[2-4]。但植物在长期的进化过程中会形成特定的各种适应性的生长习性和生理、形态特征[5]，一些植物则会表现出较强的重金属耐性。揭示植物适应铅胁迫的分子机制，有利于抗性植物的培育，从而发挥植物修复物种的最大优势。目前的研究已发现并推测出众多参与重金属响应的基因，其中以与重金属运输相关的蛋白基因的研究最为广泛，主要包括HMAs（重金属ATP酶）家族，ABC（ATP结合盒）转运家族、CDF（阳离子扩散促进者）家族和ZIP（锌/铁调节转运蛋白）家族等，这些蛋白基因在提高植物对重金属的耐性上起到了一定的作用。此外，还有一些基因家族对重金属离子在细胞中的运输及提高植物的抗性等有重要作用，主要有YSL（黄色条纹样）蛋白家族、Nramps（天然抗性相关巨噬细胞蛋白）家族与CTR（铜离子转运蛋白）家族等[6-7]。现已有40多种木本植物完成全基因组测序工作，但大多集中于葡萄、石榴等可食用果树[8]，而对重金属耐性植物的研究较少。

盐肤木（Rhus chinensis）是我国主要经济树种之一，具有良好的经济与药用价值。调查研究发现，盐肤木是重金属污染严重矿区为数不多的能够自然定居的木本植物，生物量较高、生长快、适应性强、根系发达、根萌蘖性强，对铅、铬等重金属表现出较强的耐性，因此可作为南方重金属污染区生态修复的先锋植物[4，9-10]。但目前关于盐肤木重金属胁迫的研究主要集中于重金属污染胁迫下的生理学响应及光谱特征[11]、迁移机制及富集特征[4]等，而有关盐肤木重金属胁迫下基因响应机制的研究鲜有报道。这极大制约了我们对盐肤木蛋白基因在重金属耐性方面的认知。目前，RNA-Seq技术是以新一代高通量测序为平台的RNA测序技术，该技术能够在单核苷酸水平上对任一物种的整体转录情况进行检测，能够更精确、更广泛地提供物种在特定状态下的基因表达情况，是深入研究复杂转录组的强大工具[12-14]。

鉴于此，本研究利用Illumina HiSeqTM2000高通量转录组测序技术对盐肤木铅胁迫下的根系组织进行转录组测序，并对差异基因进行Gene Ontology（GO）与Pathway分析，试图揭示盐肤木响应铅胁迫的相关分子机理，以期为盐肤木应用于南方土壤重金属污染修复提供理论依据。

1 材料与方法

1.1材料的种植与处理

本研究以盐肤木为供试植物。取大小均一、饱满的盐肤木种子，在98%浓硫酸中浸泡消毒90～105 min，每隔10 min用玻璃棒搅拌一次。消毒完成后，用流水清洗种子，搓去种子表面黑色蜡质，并用纯水浸泡24 h，而后洗净置于湿滤纸上培养，每2 d更换一次滤纸，并移除发霉和空壳的种子。该发芽过程在20±2℃人工气候培养箱中避光进行。待芽长至约2 cm，将其移至营养土中栽培，待植株长至20 cm左右，取长势较为一致的幼苗进行土培胁迫试验。

胁迫试验土壤取自福建农林大学后山红黄壤，参考鲁如坤[15]的土壤元素测定方法，测得土壤背景值，结果见表1。试验用土1.5 kg·盆-1，每盆种植1株，重复3次。设置土壤铅浓度为 0（CK）、250（Pb250）、1000（Pb1000）mg·kg-1（实测值分别为0、278.843 、1 103.498 mg·kg-1），铅以 Pb（NO3）2溶液的形式加入土壤中，用纯水浇灌保持土壤湿度为40%～50%。土壤平衡2～3d后，将盐肤木植株移栽到供试土壤中（福建农林大学下安温室大棚），温度为25℃左右，胁迫一个月后取根进行测序试验。

表1 土壤背景值Table 1 Soil background value

1.2 RNA提取与测序

总RNA的提取采用Trizol法[16]，由广州基迪奥生物公司完成RNA的提取、质控、文库构建与测序。

1.3 基因功能表达、注释及富集分析

1.3.1 Unigene注释

Unigene基本功能注释信息给出Unigene的蛋白功能注释、COG/KOG功能注释、GO功能注释、Path⁃way注释等。通过blastx将Unigene序列比对到蛋白数据库Nr、Swissprot、KEGG（系统分析基因产物在细胞中的代谢途径以及这些基因产物功能的数据库）和COG/KOG（e-value＜0.000 01）（基因产物进行直系同源分类的数据库），得到与给定Unigene具有最高序列相似性的蛋白，从而得到该Unigene的蛋白功能注释信息。

1.3.2 GO功能显著性富集分析

GO是一个国际标准化的基因功能分类体系，提供了一套动态更新的标准词汇表来全面描述生物体中基因和基因产物的属性。GO功能分析一方面给出差异表达基因的GO功能分类注释；另一方面给出差异表达基因的GO功能显著性富集分析，通过GO功能显著性富集分析能确定差异表达基因行使的主要生物学功能。计算得到的p-value通过FDR[17]校正之后，以q-value≤0.05为阈值，满足此条件的GO term定义为在铅胁迫下样品间差异表达基因中显著富集的GO term。

1.3.3 Pathway功能显著性富集分析

在生物体内，不同基因相互协调行使其生物学特性，基于Pathway的分析有助于更进一步了解基因的生物学功能。KEGG是有关Pathway的主要公共数据库[18]。Pathway显著性富集分析以KEGG Pathway为单位，应用超几何检验，找出与整个基因组背景相比，在差异表达基因中显著性富集的Pathway。通过Pathway显著性富集分析能确定差异表达基因参与的最主要生化代谢途径和信号转导途径。通过计算以q-value≤0.05为阈值筛选出样品间差异基因的主要代谢路径。

2 结果与分析

2.1 样品测序质量与组装结果

使用短reads比对软件Bowtie 2[19]将高质量clean reads比对参考基因序列（表2），CK、Pb250、Pb1000与参考基因的匹配率分别为85.69%、84.58%、84.56%，表明测序数据可用性较高。使用短reads组装软件Trinity[20]做转录组从头组装（表3）发现，转录本拼接长度大多集中在200～299 nt，转录本长度介于201～15 581 nt，N50值为1566，拼装组装效果符合要求。

表2 各样品测序质量统计Table 2 Sampling quality statistics of each sample

表3 组装结果统计Table 3 Assembly result statistics

2.2 Unigene基本注释统计

2.2.1 数据库注释统计

如图1所示，盐肤木70 459个Unigene与KEGG、KOG、Nr与Swissprot四大数据库分别比对得到20 363、23 655、41 047、34 476个 Unigene，共有 43 006个Unigene得到注释，27 453个Unigene未得到注释，说明盐肤木转录组中还存在大量未知基因，而有14 087个Unigene在四大数据库中均得到注释，占比32.76%。

2.2.2 物种分布统计

利用blastx将组装出来的Unigene序列与Nr数据库进行比对，有9927个Unigene与柑橘、3509个Unigene与可可、1563个Unigene与麻风树等植物的序列同源（图2）。这些基因为本研究转录组的注释提供了丰富的参考序列，并且从这些数据中可以推断出盐肤木与柑橘的进化关系较近。

图1 四大数据库注释维恩图Figure 1 Four big database annotations Venn diagram

图2 物种分布统计图Figure 2 Species distribution chart

2.3差异基因数量统计

对3个样品间的差异基因进行两两比较，用FDR与log2FC来筛选差异基因，筛选条件为FDR＜0.05且|log2FC|＞1，结果见图3。

图3 样品间差异基因数量统计Figure 3 Quantitative statistics of differential genes between samples

在本研究中，CK-VS-Pb250、CK-VS-Pb1000的上调基因表达量分别为7446、7641个，下调基因分别为8602、8134个。相对于CK而言，Pb1000的上调基因比Pb250增加了195个，表明盐肤木在较高铅浓度下一些平常不表达的基因开始表达，这与其抵御铅胁迫的耐性相关[21]；Pb1000的下调基因比Pb250降低了468个，这可能是因为铅胁迫对盐肤木的生理系统造成了伤害，使其活性降低，基因表达量下降[22]。Pb250-VS-Pb1000的上调、下调基因分别为11 286、9683个，皆大于前两者，表明盐肤木在遭受高浓度铅胁迫时，虽然大量基因表达量降低，但其自身会调用更多不常用的基因来应对逆境。

2.4 差异基因的GO分析

GO共有3个本体，分别描述基因的分子功能（Molecular function）、细胞组分（Cellular component）及参与的生物过程（Biological process）。盐肤木铅胁迫转录组基因注释结果如图4和图5所示。

CK-VS-Pb250 GO分类统计见图4：在分子功能分类中，共有12个差异表达基因显著富集，分别是结构分子活性、2类RNA结合相关基因、4类水解酶相关基因、2类磷酸酶相关基因、2类聚合酶相关基因和西格玛因子活动；在生物过程分类中，共有34个差异表达基因显著富集，分别是基因表达、电子传递链、氧化还原过程、光合电子传递链、RNA修饰、细胞离子稳态、核糖体、17类代谢过程相关基因、5类合成过程相关基因、3类复合物相关基因和2类光合作用相关基因；在细胞组分分类中，共有24个差异表达基因显著富集，分别是膜蛋白复合物、3类核糖体亚基相关基因、4类细胞器相关基因、3类细胞质相关基因、4类细胞及细胞器相关基因和9类光合作用相关基因。

CK-VS-Pb1000 GO分类统计见图5：在分子功能分类中，共有11个差异表达基因显著富集，分别是结构分子活性、四吡咯结合、6类氧化还原酶相关基因、单加氧酶活性、脂肪酶激活剂活性和rRNA结合；在生物过程分类中，共有2个差异表达基因显著富集，分别是基因表达和碳水化合物衍生物分解代谢过程；在细胞组分分类中，共有16个差异表达基因显著富集，分别是3类复合物相关基因、4类核糖体及核糖体亚基相关基因、2类细胞器相关基因、4类细胞相关基因和3类光合作用相关基因。

图4 CK-VS-Pb250差异基因GO功能注释Figure 4 CK-VS-Pb250 differential gene GO function annotation

由图4和图5可得，在分子功能上，Pb250胁迫下磷酸酶与水解酶相关基因表达量较高，而Pb1000胁迫下则是氧化还原酶，此外，四吡咯结合相关基因也有较高的基因表达量。在生物过程上，Pb250胁迫下富集了34个差异基因，而Pb1000仅有2个，表明植物受损严重，基因表达量骤降；Pb250富集量最多的差异基因是代谢过程，主要是有机物质的代谢，其中核糖体蛋白相关基因居多，如40S核糖体S3-3样蛋白、60S核糖体蛋白L7、核糖体蛋白L11等。在细胞组分上，Pb250胁迫下主要是细胞结构相关基因富集较多，也存在较少数量的核糖体亚基，而Pb1000下核糖体相关基因出现大量富集，主要为核糖体相关亚基与核糖核蛋白复合物。

2.5 差异基因的KEGG分析

由表4可得，核糖体是CK-VS-Pb250、CK-VSPb1000与Pb250-VS-Pb1000差异样品中共有的代谢通路，表明核糖体是盐肤木适应铅胁迫的主要代谢通路。在CK-VS-Pb250与CK-VS-Pb1000中分别有1211、1026个核糖体DEGs，表明高浓度铅胁迫使得核糖体相关基因表达量减少。光合作用代谢通路仅存在于CK-VS-Pb250中，表明轻度铅胁迫下，盐肤木仍能有效进行光合作用，而重度铅胁迫使得植物光合作用受到抑制，植物生长受到影响。此外，光合天线蛋白代谢通路的差异表达基因主要有叶绿素a/b结合蛋白Ⅱ型前体、叶绿素a/b结合13、叶绿素a/b结合蛋白3、预测叶绿素a/b结合蛋白CP24 10B等与叶绿素相关的基因。与CK相比，Pb1000较Pb250的光合天线蛋白通路减少了2个DEGs，分别是叶绿素a/b结合蛋白和预测叶绿素a/b结合蛋白1D样。

图5 CK-VS-Pb1000差异基因GO功能注释Figure 5 CK-VS-Pb1000 differential gene GO function annotation

表4 差异表达基因的KEGG代谢途径分析Table 4 Analysis of KEGG metabolic pathways of differentially expressed genes

3 讨论

水解酶（XTH）是植物细胞壁重构过程中的关键酶之一，不仅具有松弛细胞壁、促进细胞生长的作用，而且也具有强化细胞壁和维持细胞壁完整性的作用[23-24]。在轻度铅胁迫下有4类水解酶相关基因出现显著富集，表明盐肤木在轻度铅胁迫下通过提高细胞中水解酶活性修复受损细胞壁。同时，在拟南芥和番茄中发现，超量表达CaXTH3能够提高转基因植株抗旱性和耐盐性[25-26]。核苷三磷酸酶通过提供能量来促进核内成熟mRNA穿过核孔复合体进入胞浆，是细胞核mRNA出核转运的主要限速酶[27]；焦磷酸盐是控制细胞内焦磷酸浓度的关键酶，而焦磷酸是RNA、DNA、蛋白质及糖原等生物大分子生物学合成过程中的副产物，其浓度会影响细胞内部分生理反应的平衡[28]。在轻度铅胁迫下这两类磷酸酶相关基因出现显著富集，表明盐肤木在轻度铅胁迫下通过提高相关磷酸酶活性来调节细胞内被扰乱的平衡，以维持细胞内的正常运作。同时，Park等[29]将拟南芥H+-PPase基因转入番茄中发现，转基因植株与野生型植株相比生长速度快、根系发达，且有较强的抗旱能力。在高浓度铅胁迫下，盐肤木的各类氧化还原酶活性显著增强，氧化还原酶能够进行大量的催化反应，生物合成结构复杂的化合物[30]，以抵御铅侵害。此外，四吡咯化合物是植物光合作用、呼吸作用等生物学过程不可或缺的重要组分，维系着植物的生长和发育，当四吡咯合成途径受阻时，一些光敏型中间代谢物积累，则会产生氧化胁迫，诱导细胞死亡，严重阻碍植物生长[31]。在重度铅胁迫下四吡咯结合相关基因出现显著富集，表明高浓度铅胁迫已对四吡咯化合物的结合产生阻碍，影响了植物正常生长。施翔等[4]和Souza等[32]发现，高浓度重金属胁迫下植物生物量降低。综上，通过不同浓度铅胁迫发现，细胞内水解酶、磷酸酶等相关基因在轻度铅胁迫下主要通过提高自身活性来增强细胞抗性，而在重度铅胁迫下，则是氧化还原酶相关基因通过合成各类化合物以抵御逆境。

中心体是非膜细胞器之一，是细胞中的微管组织中心，在提供细胞器定向运输支架和细胞运动调节中发挥重要作用[33]。此外，核糖体也属于非膜细胞器，本研究发现，无论是在轻度还是重度铅胁迫下，非膜细胞器相关基因都出现大量富集，表明盐肤木在铅胁迫下细胞受损严重，而其自身首先通过细胞器调节相关的细胞运动以应对逆境。细胞质是生命活动的主要场所，绝大多数的化学反应都在细胞质中进行，主要成分为核糖体、多种酶类和中间代谢物、各种营养物等。在轻度铅胁迫下有3类细胞质相关基因出现显著富集，表明盐肤木在遭受轻度铅胁迫时，通过调节自身基因抵御侵害，如提高酶类基因表达、增强代谢等。Yao等[34]发现在番茄与苹果愈伤组织中过量表达细胞质苹果脱氢酸酶基因，不仅能够促进细胞发育和植物生长，促进质子泵相关基因的表达和三磷酸腺苷酶的生成，而且能够调节有机酸转运和溶质势，减少ROS的生成，从而提高转基因作物对寒冷和盐胁迫的抗性。综上，铅的侵入使得盐肤木根系细胞受损严重，而盐肤木会通过调节细胞内细胞器、细胞质等相关的细胞运动以应对逆境。

从差异表达基因的GO与KEGG的分析中均发现，核糖体相关基因是盐肤木应对铅胁迫的主要调节基因。核糖体蛋白不仅参与了rRNA的加工、折叠、核糖体亚基组装和转运过程，还在亚基结构的稳定性、核糖体与各种翻译因子的相互作用和新生肽的折叠与定位等过程中发挥作用，甚至还可能承担着核糖体外的生物学功能[35-37]。核糖体由rRNA和蛋白质构成，是细胞内蛋白质合成的分子机器，如果核糖体不起作用，细胞就无法分裂，植物生长将停止。与CK相比，Pb1000较Pb250核糖体通路减少了185个DEGs，因此推测高浓度铅胁迫影响盐肤木根系的正常代谢，抑制植株正常生长，通过对根系生理特性如SOD、CAT等抗氧化酶和MDA的测定分析也印证了这一点[38]。此外，光合作用代谢通路中的捕光叶绿素a/b结合蛋白等能够促进叶绿体合成辅酶并提高光能利用率[39]，表明在高浓度铅胁迫下，光合作用受到一定的抑制，植物生长受到影响，笔者在对叶绿素含量的研究分析中也发现铅显著抑制了盐肤木叶绿素含量[38]。姚广等[40]也发现，铅胁迫显著抑制了玉米地上部分和地下部分的生长、降低了叶片光合色素含量，杜连彩[41]也有相似的研究结果。综上，核糖体相关基因是盐肤木应对铅胁迫的主要应答与调节基因。