盐敏感型与耐盐型高粱对盐胁迫反应的转录组差异分析

范 娜，彭之东，白文斌

（山西农业大学高粱研究所 / 高粱遗传与种质创新山西省重点实验室，山西榆次 030600）

盐胁迫是一类重要的非生物胁迫，严重制约了农业的可持续发展[1]。盐胁迫条件下，土壤溶液渗透离子浓度增高，渗透势发生改变，导致细胞质膜透性改变、影响植物根系对矿质元素的吸收和植株形态的发育，抑制光合作用、呼吸作用以及物质代谢等途径[2]。盐胁迫对作物的危害是一个复杂的生理过程，而作物抵抗盐胁迫的能力主要由植物本身的遗传基因控制，不同的作物对盐胁迫的响应方式也有所差异[3]。高粱作为一种重要的粮食兼经济作物，具有较强耐盐能力。根据土壤盐分变化规律可知，高粱苗期对盐分胁迫最为敏感，苗期生长的好坏影响高粱的成熟度及产量[4]。目前国内外对高粱耐盐性的研究大都集中在形态和生理水平上，对高粱耐盐胁迫的转录组分析、高粱耐盐代谢途径和关键耐盐调控基因挖掘等方面研究较少[5–6]。转录组测序技术(RNA-seq)在农作物抗逆性研究中应用的越来越广泛，对植物抗逆基因发掘和机制的研究具有重要意义。

本研究以高粱作为试材，结合苗期生长指标的测定，对盐胁迫处理表型差异较大的高粱材料进行转录组测序，对测序数据进行功能注释，从转录组水平宏观分析盐胁迫下基因的表达特性，通过差异基因聚类分析、pathway分析和实时荧光定量RTPCR技术等，从分子水平研究植物盐胁迫反应的生理机制。发掘其中的关键调控基因，探求高粱在盐胁迫下的耐受机制，旨在为高粱育种和耐盐栽培提供一定的生理生化和分子理论基础。

1 材料和方法

1.1 试验材料及室内耐盐鉴定

前期采用室内盐分模拟和大田鉴定相结合的方法对600份高粱亲本材料进行筛选，通过测定发芽势、发芽率、发芽指数、根长、叶长、叶片酶活性、叶片Na+和K+含量、幼苗生长指标、大田农艺指标及产量指标，筛选出耐盐表型材料67B和盐敏感表型材料3560R[7–8]，作为本次试验的供试材料。

1.1.1 叶片酶活性测定 取盐处理(150 mmol/L NaCl溶液处理) 12 h后的叶片材料，剪取擦拭干净叶3～5片，3次重复，经液氮处理放入超低湿冰箱中保存待测，利用南京建成试剂盒测定超氧化物歧化酶(SOD)、过氧化氢酶(CAT)和丙二醛(MDA)含量。

1.1.2 叶片 Na+、K+含量测定 取盐处理 12 h 后叶片材料，采用火焰光度计法测定Na+、K+含量[7]。

1.1.3 幼苗生长指标测定 用直尺测定植株根长、芽长，用电子天平称量植株鲜重[8]。

1.2 转录组测序

分别取盐处理及对照组(加入等体积蒸馏水处理) 6、12 h叶片材料，前期试验数据表明该时间点是高粱响应盐胁迫最敏感的时期。为了排除外在因素影响，确保数据的准确性，设3次重复，并保证材料选取的平行性。对每个样品分别进行 RNA提取，用于后续试验[9]。

1.2.1 总 RNA 的提取 采用 TIANGEN RNA plant plus Reagent试剂盒提取总RNA。利用1%琼脂糖电泳、Nanodrop2000、Qubit 2.0 和 Aglient2100 方法检测RNA样品的纯度、浓度和完整性。

1.2.2 数字基因表达谱分析 (DGE)文库构建和转录组测序，测序数据生物信息学分析 cDNA文库的构建及转录组测序样品RNA由百迈客生物科技有限公司制备。RNA样品通过质量检测进入Illumina HiSeq 2000平台进行转录组测序。测序得到的所有原始reads (raw data)经过过滤得到纯净数据 (clean data)，纯净数据再与指定参考基因组(Sbicolor_v2.1)比较得到比对数据(mapped data)。用HISAT2软件将纯净数据与指定参考基因组进行比对，再用StringTie软件将比对上的reads进行组装，从而得到完整的比对数据 (mapped data)。

1.2.3 差异表达基因筛选 将盐胁迫处理的样品测序得到的基因表达数据与对照处理的样品测序得到的基因表达数据进行比较，得到差异表达基因。该检测过程中，采用Benjamini-Hochberg 的校正方法，筛选标准为差异倍数(fold change，FC)≥2 且错误发现率 (false discovery rate，FDR)<0.001。差异倍数(FC)=处理组某基因表达量/对照组该基因表达量。

1.2.4 差异基因的注释和分类 将差异表达基因分别与 NR (NCBI non-redundant protein)、Swiss-Prot(Swiss-Prot proteinsequence)、KEGG (Kyoto encyclopedia of genes and genomes)、COG (clusters of or thologousgroups of proteins)和 GO (gene ontology)数据库进行比对，获得与差异表达基因对应的蛋白功能注释信息。

1.2.5 实时荧光定量PCR分析(qRT-PCR) 为了验证高粱耐盐胁迫转录组结果的准确性，从差异表达基因数据库中随机挑选5个基因，分析实时荧光定量PCR结果是否与转录组结果一致。

2 结果与分析

2.1 耐盐鉴定

表1为盐胁迫处理15天后高粱幼苗生长状况，耐盐材料生长速率快，表现出较强的耐盐性，而盐敏感材料盐胁迫下的根短、鲜重轻和生长速率低，对盐胁迫较敏感。高粱苗期对盐胁迫的响应表现与后期高粱生长表现一致。

表1 高粱幼苗生长状况Table 1 Growth indices of sorghum seedlings

盐胁迫下导致细胞质中的Na+/K+增高，使液泡膜ATP酶和焦磷酸酶活性下降，从而抑制植物的正常生长，因此是否具有高Na+/K+是判断植物耐盐性的重要依据之一[10]。表2显示，盐胁迫条件下，耐盐材料Na+含量小于盐敏感材料，耐盐材料可以提高Na+的选择吸收及其在植株体内的积累与分配；耐盐材料和盐敏感材料高粱叶片K+含量差异不显著。

表2 盐胁迫对高粱叶片K+、Na+浓度和K+/Na+值的影响Table 2 Effects of salt stress on K+ and Na+ concentration and K+/Na+ in sorghum leaves

植物处于盐分逆境时，耐盐材料比盐敏感材料的SOD、CAT活性高，MDA含量低，说明植物对逆境抵抗能力的强弱与维持体内较高抗氧化酶活性的能力强弱有关。表3显示，耐盐材料抗氧化酶活性都较强，清除活性氧、过氧化氢的能力较强。在盐胁迫条件下，耐盐材料可以维持较高的过氧化氢酶活性，受到盐胁迫后该酶活性升高幅度相对较大，进而保持了较强的过氧化氢清除能力，能够及时清除过量积累的活性氧。

表3 盐胁迫对高粱叶片SOD、CAT、MDA含量的影响Table 3 Effects of salt stress on the contents of SOD, CAT and MDA in sorghum leaves

2.2 高粱耐盐转录组

2.2.1 转录组测序质量可行性分析 盐胁迫组与对照组所有样本(共12个)转录组测序结果经过测序质量控制筛选，共得到 90.07 Gb Clean Data，每个样本测序所得的数据量均大于 6 G (6.1～8.2 G)，测序所得的clean reads数分别为21608762和30477544条，12个样本的转录组序列GC (鸟嘌呤和胞嘧啶)含量介于55.19%到56.85%之间，略高于AT (腺嘌呤和胸腺嘧啶)含量，各样品Q30碱基百分比均不小于94.40%。将转录组数据与参考基因组进行比对，各样本的reads比对效率在93.07%～95.11%。综合以上几个测序reads质量评价指标，说明这12份样本的测序质量较高，满足DEG筛选及数据分析的要求。

2.2.2 与参考基因组比对结果统计 盐胁迫下两个品系共有5040个差异表达基因。对差异表达基因在所选参考基因组染色体不同位置上的reads的分布数量进行生物信息学分析，分析比对到的序列在10条染色体上的分布情况。将比对到的不同染色体上的reads进行位置分布统计，绘制比对到的序列(mapped reads)在所选参考基因组上的覆盖深度分布图(图1)。

如图1所示，测序结果在高粱的10条染色体上都有reads覆盖，且覆盖率差别不大，即在整个基因组中都有reads分布且分布均匀，其中覆盖比例最高的是基因NC_012870.2，其次NC_012871.2和NC_012872.2(NC—耐盐对照)。

图1 Mapped reads在参考基因组上的位置及覆盖深度分布图Fig.1 Map of the location and depth of coverage of mapped reads on the reference genome

参考基因组的reads分为外显子、内含子和基因间区。理论上，比对到外显子区的mRNA的reads其基因组注释较为清晰，比对到内含子是由于mRNA前体和发生可变剪切的内含子保留，而比对到基因间区则表明基因组注释不完善。如图2所示，12个高粱样本比对到reads在外显子、内含子与基因区间的比例分别为91.19%，3.92%，4.89%。

图2 基因组reads在3个区域的分布比例Fig.2 Percentage of the genome reads in the three regions

12个高粱样本进行转录组测序及比对拼接共得到31272个基因序列，其中在数据库中有功能注释的有29901个。这些基因在KEGG和KOG注释比例仅为31.1%和49.7%，说明虽然高粱基因组已经全测序完毕，但大部分基因的功能及参与调控途径并未研究清楚，科研工作者们还有大量的研究工作要做。

2.2.3 差异表达基因分析 在相同盐胁迫条件下，盐敏感材料和耐盐材料的差异表达基因数是不同的(图3)。盐敏感材料中共有2772个差异表达基因，其中上调基因有1443个，下调基因有1329个。耐盐材料中共有2268个差异表达基因，其中上调基因有1217个，下调基因有1051个。盐敏感材料的差异表达基因数目略多于耐盐材料，两品系中上调基因均略多于下调基因，差异表达基因的差异倍数主要集中在2～5倍。

图3 盐胁迫处理后高粱盐敏感和耐盐材料差异表达基因表达量变化Fig.3 Volcanic diagram of gene expression levels of salt-sensitive and salt-tolerant sorghum materials after salt stress

将耐盐材料和盐敏感材料的差异表达基因做韦恩图分析，如图4所示，盐敏感材料和耐盐材料差异表达基因有1009个是相同的，这些基因可能是高粱响应盐胁迫的共同表达基因，盐敏感材料中有1763个基因特异表达，耐盐材料中有1259个基因特异表达，这些特异表达的基因可能决定着两个品系在同一盐胁迫处理条件下的不同表型，说明盐胁迫造成了两个品系的基因表达变化。

图4 两个品系差异表达基因韦恩图Fig.4 Venn diagram of differentially expressed genes in two sorghum materials

根据表4和图5可知，N-1-1_N-1-2_N-1-3 (耐盐材料处理 6 h) vs N-2-1_N-2-2_N-2-3 (耐盐材料处理12 h)差异表达基因2268个，其中上调表达基因1217个(N-1-1表示耐盐材料盐胁迫处理|重复|，其他类推)；Y-1-1_Y-1-2_Y-1-3 (盐敏感材料处理 6 h)vs Y-2-2_Y-2-1_Y-2-3 (盐敏感材料处理 12 h)差异表达基因352个，其中上调表达基因178个(Y-1-1表示盐敏感材料盐胁迫处理|重复|，其他类推)；N-1-1_N-1-2_N-1-3 vs Y-1-1_Y-1-2_Y-1-3 差异表达基因1387个，其中上调表达基因612个；N-2-1_N-2-2_N-2-3 vs Y-2-1_Y-2-2_Y-2-3 差异表达基因 1241个，其中上调表达基因536个；Y-1-1_Y-1-2_Y-1-3 vs Y-2-1_Y-2-2_Y-2-3差异表达基因2772个，其中上调表达基因1443个。

图5 两个品系盐胁迫12 h与盐胁迫6 h相比差异表达基因的韦恩图Fig.5 Venn diagram of differently expressed genes between salt stress of 12 h and 6 h of the two sorghum cultivars

表4 不同高粱材料盐胁迫12 h与胁迫6 h相比差异表达基因的表达数Table 4 Number of differently expressed genes (DEGs) between sorghum cultivars at salt stress from 6 h to 12 h

耐盐材料和盐敏感材料中的差异表达基因，按照盐胁迫后的变化幅度进行聚类(图6)。盐敏感材料中的差异表达基因可分为4类：第一类为表达上调幅度较大的基因，共295个(占10.6%)，第二类为表达上调幅度较小的基因，共249个(占9.0%)，第三类为表达下调幅度较大的基因，共377个(占13.6%)，第四类为表达量虽有差异但变化幅度不大的基因，共1851个(占66.8%)。耐盐材料中的差异表达基因可分为5类：第一类为表达下调幅度较小的基因，共53个(占2.3%)，第二类为表达下调幅度较大的基因，共283个(占12.5%)，第三类为表达上调幅度较大的基因共222个(占9.8%)，第四类为表达上调幅度较小的基因共51个(占2.2%)，第五类为表达变化幅度不大的基因共1659个(占73.1%)。表达上调幅度较小的基因比例盐敏感材料中(9.0%)高于耐盐材料中(2.2%)。

图6 盐胁迫处理下两个品系差异表达基因的聚类分析Fig.6 Cluster analysis of differently expressed genes of the salt-tolerant and salt-sensitive cultivars after salt stress treatment

2.2.4 差异表达基因的GO分类与富集分析 GO数据库是一个结构化的标准生物学注释系统，旨在建立基因及其产物知识的标准注释体系。GO分析有助于我们理解基因的功能，GO条目即基因的注释条目。

对盐敏感材料和耐盐材料的所有差异表达基因进行GO功能富集分析(图7)，覆盖了绝大多数GO条目。盐胁迫后，盐敏感材料中的差异表达基因涉及了52个GO分类条目，其中与生物过程相关的有22个条目，上调表达基因522个，细胞过程上调表达基因430个，单组织过程上调表达基因343个；与分子功能相关的有15个条目，基因数目较多的是催化活性上调表达基因509个，及结合功能上调表达基因509个。

图7 盐敏感材料和耐盐材料中差异表达基因GO分类Fig.7 GO classification of differently expressed genes in salt-sensitive and salt-tolerant materials

耐盐材料在盐胁迫处理后的差异表达基因涉及了三大一级条目中的52个二级GO分类条目，其中与生物过程相关的有22个条目，基因数目较多的是代谢过程上调表达基因494个，细胞过程上调表达基因373个，单组织过程上调表达基因330个；与分子功能相关的有13个条目，基因数目较多的是结合功能上调表达基因433个与催化活性上调表达基因475个。

2.2.5 差异表达基因的KEGG分类及Pathway显著性富集分析 从盐胁迫下盐敏感材料的差异表达基因KEGG分类(图8)可知，差异表达的基因分布在110条pathway中，其中差异基因数目最多的前5条pathway分别是植物激素信号转导(25条，6.53%)，碳代谢(24条，6.27%)，苯丙烷类合成(24条，6.27%)，淀粉与蔗糖代谢(23条，6.01%)及氨基酸生物合成(21条，5.48%)，这些pathway都与代谢通路和信号转导有关。

图8 盐敏感材料盐胁迫下差异表达基因的KEGG分类Fig.8 KEGG classification of differently expressed genes in salt-sensitive materials under salt stress

由耐盐材料的差异表达基因的KEGG分类(图9)可知，差异表达的基因分布在104条pathway中，差异表达基因数目最多的前5条pathway分别是苯丙烷类合成(20条，6.73%)，植物激素信号转导(16条，5.39%)，DNA复制(15条，5.05%)，植物病原菌互作(15条，5.05%)及碳代谢(13条，4.38%)，这些途径与信号转导通路和代谢通路有关。前5个pathways中，盐敏感材料的差异表达基因在代谢通路中的占比高，而耐盐材料的差异表达基因在传导信号中的占比较高。

图9 耐盐材料盐胁迫后的差异基因KEGG分类Fig.9 KEGG classification of differently expressed genes in salt-tolerant material after salt stress

对差异表达基因参与的pathway进行富集分析。盐敏感材料富集在前5位的pathway分别是硫辛酸代谢，二苯乙烯、双苯庚烷与姜辣素合成，赖氨酸生物合成，有机含硒化合物代谢，DNA复制(图10)。耐盐材料富集在前5位的pathway分别是硫辛酸代谢、亚麻酸代谢、光合作用、硫胺代谢、DNA复制(图11)。盐敏感材料中差异表达基因主要集中在基础代谢和次生物质 (basic metabolism and secondary substances)合成途径，应该是盐胁迫引起了植株细胞内许多物质的降解引起的，即主要是盐胁迫伤害所致。耐盐材料中差异表达基因参与的途径除了次生产物合成途径外，还包括信号转导途径和防御(signal transduction pathways and defense)途径。

图10 盐敏感材料的差异表达基因KEGG富集图Fig.10 KEGG enrichment of differently expressed genes in salt-sensitive material

图11 耐盐材料的差异表达基因KEGG富集图Fig.11 KEGG enrichment diagram of differently expressed genes in salt-tolerant material

2.2.6 实时荧光定量PCR验证 为了验证转录组测序结果的可靠性，从差异表达基因的数据库中随机选取5个基因进行实时荧光定量PCR验证，以ACTIN-1为内参，利用2–△△Ct计算各个样品及各组中的基因相对表达量，5个基因的RPKM值(RPKM是 reads per kilobase per million mapped reads 的缩写，代表每百万reads中来自于某基因每千碱基长度的reads数。本试验以5个基因对照的RPKM值作为基准进行比较，RPKM值与转录组分析数据差异均不显著，且变化趋势基本一致，表明9RT-PCR分析结果与测序结果吻合，转录组测序结果准确(图12)。

图12 目的基因各样品中的相对表达量Fig.12 Relative expression of target gene in each sample

3 讨论

转录组测序通过分析基因在盐胁迫处理前后表达的模式，发掘与逆境相关的基因以及代谢通路，是研究植物盐胁迫分子机制的重要手段[11]。

3.1 盐胁迫对高粱幼苗生长的影响

室内试验表明，植物处于盐分逆境时，耐盐材料生长速率快，表现出较强的耐盐性，Na+含量和Na+/K+明显增加，而K+含量明显下降。耐盐材料Na+上升的幅度明显小于盐敏感材料，耐盐材料SOD、CAT活性高，说明植物对逆境抵抗能力的强弱与维持体内较高抗氧化酶活性的能力强弱有关。耐盐材料抗氧化酶活性都较强，清除活性氧、过氧化氢的能力较强；在盐胁迫条件下，耐盐材料可以维持较高的过氧化氢酶活性，受到盐胁迫后该酶活性升高幅度相对较大，进而保持了较强的过氧化氢清除能力，能够及时清除过量积累的活性氧，细胞膜质过氧化程度较低。总之盐胁迫下，耐盐材料生长较好、生长速率快，表现出较强的耐盐性；盐敏感材料相对生长缓慢，耐盐性弱，研究结果与董明等[5]、宝力格[12]研究一致。

3.2 差异表达基因表达水平聚类分析

转录调节是植物应答逆境的关键步骤，盐胁迫下同一物种不同品系耐盐性不同其基因转录有不同的应答模式[12]。本试验表明，盐胁迫处理后，两品系材料的差异表达上调基因数量略多于下调基因数量，但是盐敏感材料的差异表达基因数目略多于耐盐材料。

将盐敏感材料中的2772个差异表达基因和耐盐材料中的2268个差异表达基因做层次聚类分析，盐敏感材料中的差异表达基因按照表达模式的不同可分为4类，耐盐材料中的差异表达基因按照表达模式的不同可分为5类。

3.3 差异表达基因分析

GO富集分析显示激素响应、光合作用以及碳代谢相关途径，旨在建立基因及其产物知识的标准注释体系。GO分析有助于我们理解基因的功能，GO条目即基因的注释条目[13]。KEGG分析还发现差异表达基因涉及激素信号传导、光合作用以及类黄酮代谢途径，这些代谢途径在应对盐胁迫过程中起着重要作用。

本研究表明，盐敏感材料和耐盐材料的差异表达基因GO分类在各注释条目上的分布趋势大致相同，各条目上的差异表达基因数目相差不大，富集情况基本一致，说明盐敏感材料和耐盐材料对盐胁迫的响应途径是相同的。盐敏感材料中各条目上的上调表达基因数略多于下调表达基因数，而耐盐材料中各条目上的上调表达基因数目显著多于下调表达基因数目，说明耐盐材料是通过增加上调表达基因的数量来抵抗逆境的。在盐胁迫中起主要作用的是调控结合功能和催化活性的基因，参与的主要生物过程有细胞过程、单组织过程和代谢过程，这与董明等[5]的研究结果一致。

盐敏感材料和耐盐材料的差异表达基因在KEGG各pathway中的分布趋势差别很大，差异表达基因数排名前5的条目有3条相同，两条不同。盐敏感材料中两条不同的条目分别是淀粉与蔗糖代谢及氨基酸生物合成，与基础代谢有关。耐盐材料中两条不同的条目分别是植物病菌互作和DNA修复，推测与其耐盐性有关。盐胁迫下，耐盐材料和盐敏感材料代谢通路有所不同。本研究结果与张飞等人的研究结果[14–16]相符，耐盐和盐敏感材料在对照和盐渍胁迫条件下在苯丙烷类生物合成、苯丙氨酸代谢、类黄酮生物合成3个途径中差异表达基因较多，可能是造成两个材料耐盐性差异的重要原因。

4 结论

高粱盐敏感材料和耐盐材料耐盐基因表达差异显著。耐盐材料调控基因涉及植物激素信号转导、碳代谢、苯丙烷类合成、淀粉与蔗糖代谢及氨基酸生物合成。盐敏感材料上调表达基因包括：细胞过程上调表达基因430个，单组织过程上调表达基因343个，催化活性上调表达基因509个，结合功能上调表达基因509个；耐盐材料上调表达基因包括：代谢过程上调表达基因494个，细胞过程上调表达基因373个，单组织过程上调表达基因330个，结合功能上调表达基因433个，催化活性上调表达基因475个。