脱水对栓皮栎种子基因表达的影响1)

李东兴 钱家连 许慧慧 张国伟 任俊杰 王茜 许洋 王利兵 于海燕 李迎超

(中国林业科学研究院林业研究所，北京，100091) (河北省洪崖山国有林场)(河北邢台市信都区农业农村局) (中国林业科学研究院林业研究所)

栓皮栎(Quercusvariabilis)为壳斗科(Fagaceae)栎属(Quercus)落叶乔木，耐贫瘠和干旱，是营造防护林、水源涵养林和防风林的优良树种，具有较高的经济和生态价值[1-2]。栓皮栎树皮也称软木，作为一种可再生的生物质材料，被广泛用于葡萄酒瓶塞和保温隔热等诸多领域[3]；栓皮栎种子中营养丰富，淀粉含量高，可用于酿酒和提取浆纱，也是生产燃料乙醇的“非粮”原料[2]。同时，栓皮栎壳斗也是制备活性炭和糖醛等的优质原料。

然而，绝大多数栎属植物无性繁殖困难[4-5]，其苗木繁育主要依靠播种[6]，种子质量决定了成活率及苗木质量。栎属植物种子成熟后期都不经历成熟干燥后获得脱水耐性的过程，这导致它们对脱水敏感[7]。栓皮栎种子是典型的脱水敏感顽拗性种子，其成熟脱落时含水量高，干燥后在低温(-18℃左右)环境中贮藏的常规方法会使种子迅速失去活力，进而死亡[8-9]。顽拗性种子脱水敏感性是植物重要的功能特征和由多基因共同响应的次级进化[10-11]。种子脱水过程中，细胞感受到水分流失信号，信号经传导后，在转录和翻译等不同水平上做出响应，以调控相关基因表达和蛋白合成，应对脱水造成的伤害。顽拗性种子脱水过程中，与转录因子(MYB1R1、ERF106、ABI3、RL40和TGA10)、激素信号传导(PP2C、JAZ1、JAZ2和JAZ3)和抗氧化酶(SUPEROXIDE DISMUTASE和CATALASE)相关的基因差异表达，它们可能在种子脱水响应中发挥重要调控作用，参与种子脱水敏感性的分子调控[13-15]。不同顽拗性种子脱水敏感性不尽相同[16-17]，其调控机制也不相同。虽然已有研究对栓皮栎顽拗性种子脱水敏感性进行探究，指出种子最佳贮藏温度为(0～2)℃[8]，但是种子脱水敏感性分子调控机制在很大程度上仍然是未知的。

因此，本研究以从同一栓皮栎优良单株上收集的种子为试验材料，进行脱水和萌发试验，然后对脱水敏感关键时期的种子进行转录组测序，对获得差异表达基因的功能和其参与的代谢通路进行分析，筛选种子脱水敏感性的候选基因，分析其可能的调控作用。以期为栓皮栎种子脱水敏感性分子调控机制的研究提供参考，也为其苗木繁育、种子贮藏与种质资源长期保存提供一定理论依据。

1 材料与方法

供试种子来源于陕西省周至县楼观台林场中的同一栓皮栎优良单株，2019年9月中旬获得种子后，首先通过目视法手工去除杂质和肉眼可见的劣质种子，然后用55 ℃温水浸泡水选，5 min后所有仍然漂浮在水面上的种子被丢弃[18]，剩下的种子经沥水且用滤纸擦干表面水分后，平铺在阴凉处阴干(12～15 h)。在此期间，每隔2～3 h翻动1次种子。阴干完成后，挑选健康、均匀且没有萌发的种子，在国家林业和草原局北方林木种子检验中心内进行试验。

1.1 脱水试验设计

参照程继铭等[19]的方法，采用硅胶快速脱水法进行脱水试验。脱水处理梯度为1、2、3、4、5、6、7、8、9、10、11、12、13、14、15 d(T1～T15)，共15个处理，以不经脱水处理的种子为对照(CK)。每个处理均随机挑选45粒种子，重复3次，共135粒种子，分别置于装有硅胶的大号自封袋中(m(种子)∶m(硅胶)=1∶3)，混合均匀后，在室温约25 ℃，相对湿度为40%～45%的实验室中进行脱水试验。到达预脱水时间后，每个处理分别随机取100粒种子进行萌发试验，30(每组10粒分3组)粒种子经液氮速冻后，迅速放置于冰箱里-80 ℃储存，用于转录组测序分析。

1.2 种子含水量测定

参照GB 2772—1999《林木种子检验规程》，采用稍加改进后的低恒温烘干法测定种子初始含水量(CIM)。首先将样品盒和盖烘干至恒质量后，把种子迅速切成片状放入样品盒中，然后在103 ℃的烘箱中烘17 h(达到预设温度后开始计时)。3次重复，每次重复10粒种子，CIM按照以下公式计算：

CIM=((L-N)/(L-P))×100%。

式中：L为样品盒和盖以及样品的烘前质量(g)；N为样品盒和盖以及样品的烘后质量(g)；P为样品盒和盖的质量(g)。

不同脱水处理后种子含水量用种子初始含水量减去脱水量计算得到。

1.3 萌发试验

选择直径为12 cm的培养皿，双层滤纸加一层纱布做为发芽床，然后置于光暗周期为30 ℃光照8 h，20 ℃黑暗16 h的光照培养箱中进行萌发试验[20]，5次重复，每个重复20粒种子。以胚根露出至少2 mm作为萌发标志[21]，每天观察并记录1次萌发情况，萌发持续时间为28 d[22]。种子萌发率(PG)按照以下公式计算：

式中：Dt表示在t日时的萌发数量；N表示种子总数。

1.4 种子总RNA提取、转录组测序和分析

用Trizol试剂(Invitrogen,CA,USA)提取栓皮栎种子总RNA后，分别用Agilent 2100 Bioanalyzer和RNA 6000 Nano LabChip Kit(Agilent,CA,USA)检测所提样品RNA的质量和纯度。样品合格后在杭州联川生物技术股份有限公司Iilumina测序平台(Illumina HiSeqTM4000)进行高通量测序，得到转录组数据。然后使用Cutadapt软件对测序产生的原始数据进行预处理，过滤掉不合格的序列后得到有效数据，接着用HISAT2将获得的有效数据与参考基因组[23]比对。使用edgeR对String Tie组装和定量完的基因进行差异分析，获得差异表达基因(DEGs)。对差异表达基因进行GO和KEGG富集分析，并以P<0.05作为标准，筛选差异表达基因显著富集到的GO类别和KEGG代谢途径。

1.5 实时荧光定量PCR

为验证RNA-seq结果的准确性，选择4个差异表达基因进行qRT-PCR验证。使用PrimeScriptTMReagent Kit with gDNA Eraser试剂盒将RNA反转录成第一条链的cDNA。然后按照制造商的说明，用KAPA SYBR FAST qPCR Master Mix(KapaBiosystems,USA)进行qRT-PCR试验。使用NCBI上的Primer-BLAST设计引物，引物序列如表1所示，以PP2A[24]为内参基因，采用2-ΔΔCt计算相对表达量[25]。

表1 qRT-PCR引物序列

2 结果与分析

2.1 种子萌发率变化

随着脱水时间的延长，栓皮栎种子含水量和萌发率呈现出一直降低的趋势(图1)。未经脱水处理种子的萌发率及含水量分别为84.00%和34.27%，脱水12 d(T12)后，种子萌发率降为0，此时种子含水量为17.17%。脱水初期(CK～T2)种子萌发率下降最快；经过一个较慢的下降阶段后(T2～T4)，萌发率又迅速下降(T4～T6)，最后缓慢降至0。此外，Pearson相关性分析结果表明栓皮栎种子含水量与萌发率间存在极显著的正相关(P<0.01,R=0.976)，这说明其对水分流失高度敏感。此外，临界含水量和致死含水量分别指50%左右和全部种子无法存活时的含水量，能够作为顽拗性种子脱水敏感性的标准[26-27]。本研究中，当栓皮栎种子含水量低于28.20%(T2)时，近一半种子失去生活力，即栓皮栎种子临界含水量可能在28.20%左右；当种子含水量低于17.79%(T11)时，种子基本全部失去生活力，即种子致死含水量可能在17.79%左右。因此，为探究栓皮栎种子脱水过程中的基因表达变化，我们对未脱水、脱水2 d和11 d的栓皮栎种子样品进行转录组测序分析。

竖线表示均值加减标准误

2.2 转录组测序质量分析

T11、T2和CK3个时期栓皮栎种子转录组测序结果表明共得到了475 939 682条原始序列，过滤掉带接头(Adaptor)、含有N(N表示无法确定碱基信息)的比例大于5%和低质量(质量值Q≤10的碱基数占整个序列的20%以上)的读取片段后，共得到450 383 524条有效序列。它们的GC含量和Q30的值分别是43.50%～44.00%和98.05%～98.28%(表2)。与参考基因组(GCF_002906115.1)比对后，86.73%～88.33%的测序片段能够比对到参考基因组上，而且51.03%～52.61%的片段能唯一比对到参考基因组一个位置上。这些结果表明了本次测序数据质量较好，能够进行后续差异表达分析。

表2 转录组测序数据统计

2.3 差异表达基因的鉴定

对栓皮栎种子转录组测序结果进行差异表达基因分析。最终，以基因差异表达倍数≥1并且P<0.05为筛选标准，在T2与CK比较组中获得4 405个差异表达基因，包括2 724个上调和1 681个下调基因；在T11与T2比较组里鉴定出1 765个上调和2 301个下调，共4 066个差异表达基因；有3 208个上调和2 699个下调，共5 907个差异表达基因，在比较组T11与CK中被获得(表3)。这些差异表达基因可能在栓皮栎顽拗性种子脱水过程中发挥重要作用，与其脱水敏感性密切相关。

表3 不同比较组差异表达基因数

2.4 差异表达基因的GO和KEGG分析

基因本体论(GO)分析能够全面描述生物体中基因属性，其总共有3个主要类别，分别描述基因的分子功能、细胞组分和生物学过程。如图2所示，在3个比较组中，上调和下调的差异表达基因均主要富集在“防御反应”、“转录调控，DNA模板”和“氧化还原过程”等生物学过程类别；“细胞核”、“细胞质”和“质膜”等细胞组分类别以及“蛋白结合”、“ATP结合”和“DNA结合转录因子活性”等分子功能类别中。

图2 差异表达基因的GO富集分析

生物体内，不同基因相互协调行使其生物学功能，基于途经的分析有助于更进一步了解基因的生物学功能。本研究对获得的差异表达基因进行了KEGG富集分析，结果发现在比较组T2与CK中，上调的差异表达基因显著富集在“蛋白质在内质网中的加工”和“植物激素信号传导”等途径；“异黄酮生物合成”、“谷胱甘肽生物合成”和“植物激素信号传导”等是T11与T2中显著富集的通路；“谷胱甘肽生物合成”和“植物激素信号传导”等代谢途径在T11与CK中显著富集(图3)。比较组T2与CK中，下调的差异表达基因在“植物激素信号传导”和“MAPK信号级联”等途径中显著富集；“其它多糖降解”和“柠檬烯和蒎烯降解”等途径是T11与T2中显著富集的通路；T11与CK中，“植物激素信号传导”和“淀粉和蔗糖代谢”等是差异表达基因显著富集的通路(图3)。这些代谢通路有助于了解栓皮栎种子应对脱水时的代谢信息，以便更好地探究其脱水敏感性的潜在调控机制。

2.5 转录因子对脱水的响应

本研究中，我们在比较组T2与CK、T11与T2和T11与CK分别发现117、110和183个差异表达的转录因子基因(表4)。其中，MYB108(LOC112035993)、MYB1R1(LOC111996955)、WRKY28(LOC111998523)、WRKY72(LOC111987575)和ERF1B(LOC111996147)等14个转录因子表达持续上调；UNE10(LOC112020213)和CRF4(LOC112020213)等5个转录因子表达持续下调；ABI3(LOC112020213)、NGA1(LOC112036666)和RAP2-11(LOC112003881)等10个转录因子表达先上调后下调；转录因子MYB102(LOC111989495)和SRM1(LOC112023214)表达先下调后上调(表5、图4)。

表4 与转录因子和Ca2+途径相关的差异表达基因

2.6 Ca2+途径对脱水的响应

本研究中，我们在比较组T2与CK、T11与T2和T11与CK分别发现12、4和13个与Ca2+途径相关的差异表达基因(表4)。其中，有3个基因在3个比较组中均差异表达。钙调蛋白11(LOC111999077)表达持续上调；钙调蛋白3(LOC112017550)表达先上调后下调；钙调蛋白(LOC111983139)表达先下调后上调(图5)。

图3 差异表达基因的KEGG富集分析

表5 与转录因子相关的差异表达基因

图4 与转录因子相关的差异表达基因热图

图5 与Ca2+途径相关的差异表达基因热图

2.7 与植物激素信号传导相关的DEGs

脱水2 d后，我们发现与茉莉酸(JA)信号传导相关，编码茉莉酸—氨基酸合成酶的1个JAR1、茉莉酮酸酯ZIM结构域蛋白的2个JAZ、髓细胞组织增生蛋白(MYCs)的4个基因表达均上调；在比较组T11与T2中，2个JAZ基因表达上调，1个JAZ和5个MYC2基因表达下调(图6A；图7A和B)。此外，我们同样发现脱水2 d后，与水杨酸(SA)信号传导相关，编码TGACG基序结合因子的4个基因表达上调，1个基因表达下调；在比较组T11与T2中，1个TGA基因表达下调，1个编码病程相关蛋白1的基因表达上调(图6B；图7A和B)。比较组T11与CK中，与JA和SA信号传导相关的DEGs最多(图7C)。

A显示了茉莉酸信号传导的途径；B显示了水杨酸信号传导的途径。圆圈代表代谢产物，方框代表调节基因。红色表示基因表达上调，蓝色表示下调，橘色表示该基因家族中既有上调基因又有下调基因。

A代表T2与CK;B代表T11与T2;C代表T11与CK。

2.8 qRT-PCR验证

为验证RNA-Seq数据的准确性，我们对4个DEGs进行了qRT-PCR验证。尽管表达量并非完全一样，但是它们的表达趋势基本一致(表6)。这说明本研究RNA-Seq结果是可靠的，能够进行后续的差异表达分析。

表6 对4个DEGs的qRT-PCR验证

3 结论与讨论

顽拗性种子脱水后，测定含水量和萌发率是研究种子顽拗性的常用方法[19]。与我们的预期一样，栓皮栎种子萌发率随着含水量的下降而下降，并且含水量与萌发率之间呈现极显著的正相关(P<0.01,R=0.976)，相同的结果在其它顽拗性种子的研究中也被得到证实[12-13，20]。本研究中，栓皮栎种子含水量由34.27%(CK)下降到28.20%(T2)时，近一半的种子已丧失生命力，此时种子已受到严重的伤害。此外，当栓皮栎种子含水量从34.27%(CK)下降到17.17%(T11)时，种子几乎完全失去生命力，这与Yu et al.[28]对顽拗性七叶树(Aesculuschinensis)种子的研究类似。

植物中参与调控胁迫诱导型基因表达的各类转录因子，尤其是MYB、WRKY、AP2-EREBP和bHLH等转录因子家族，是信号传递网络中的关键调控因子，在植物对非生物胁迫的响应中发挥重要作用[29-30]。山茶(Camelliasinensis)种子脱水过程中，差异表达的MYB、WRKY和bHLH等转录因子家族基因可能参与种子脱水敏感性调控[12]。本研究中，4个MYB、5个WRKY、3个ERF和2个HSF类转录因子在栓皮栎种子脱水过程中持续上调。转录因子MYB108是ABA诱导细胞死亡的负调控因子，能够调节活性氧(ROS)的产生，还可能在转录水平上参与JA的信号传导，调控多种胁迫反应的信号传导[31-32]。顽拗性种子脱水过程中，ROS的积累是其不耐脱水的原因之一[26]。栓皮栎种子脱水过程中，持续上调的MYB108可能增强了ROS的产生，使种子受到脱水伤害。与顽拗性桑寄生(Taxillusichinensis)种子[13]不同，栓皮栎种子脱水过程中MYB1R1表达持续上调，这可能是因为不同顽拗性种子对脱水响应不同。转录因子ERF1B能够调控拟南芥(Arabidopsisthaliana)中JA响应基因，并作用于JA信号传导途径下游[33]，参与植物对胁迫的响应。栓皮栎种子脱水过程中，持续上调的ERF1B转录因子可能通过调控激素信号传导中相关基因的表达，在种子对脱水的响应过程中发挥作用。此外，我们同样发现2个ERF类转录因子在栓皮栎种子脱水过程中持续下调。CRF转录因子能够调控胁迫信号传导途径中的基因表达[34]，持续下调的CRF4转录因子可能干扰了栓皮栎种子中胁迫信号的传递，细胞中防御机制无法启动，受损程度不断加深，进而导致种子对脱水敏感。bHLH类转录因子能增强细胞清除ROS的能力，降低氧化应激对植物造成的伤害，提高对非生物胁迫耐受性[35]。拟南芥中过表达的AtUNE12基因使其能快速清除ROS，减少胁迫导致的细胞损伤[36]。我们发现UNE10、HEC1-like和bHLH94-like，3个bHLH转录因子基因在栓皮栎种子脱水过程中持续下调，推测这可能使种子清除ROS的能力下降，过量积累的ROS使细胞损伤程度不断加深，进而导致种子对脱水敏感。

植物对各种胁迫因子的响应中，Ca2+作为关键的第二信使，传递胁迫信号，引发细胞的防御反应，Ca2+结合蛋白(CMLs)对Ca2+信号传递至关重要[37]。CML37能通过影响茉莉酸—异亮氨酸缀合物的合成，进而将Ca2+和JA信号连接起来，CML37基因突变体中，JAR1表达降低，JA信号传导途径受到干扰[38]。本研究中，CML11在种子脱水过程中持续上调，这可能影响了Ca2+和JA信号传导，导致种子不耐脱水。

植物感受到环境胁迫信号后，体内会大量合成JA。JA在腺苷酸形成酶的催化下形成高生物活性的JA-Ile，它能与JA受体COIl特异性结合为COIl-JA-Ile复合物，然后通过招募JAZ形成COI1-JA-Ile-JAZ三元复合体，解除对下游MYC2等转录因子的抑制作用，激活JA响应基因的表达[39-40]。本研究中，脱水2 d后，JAR1,JAZ和MYC2基因表达上调，这与Luo et al.[41]对棉花(Gossypiumspp.)的研究有些类似。JAR1基因的上调可能促进了COI1-JA-Ile-JAZ复合物的形成，解除JAZ对MYC2的抑制作用，上调的MYC2通过调控下游基因的表达响应种子脱水。拟南芥脱水过程中。MYC2基因的显著下调可能是JA信号传导途径中，下游基因不表达的原因[42]。当种子继续脱水时，编码JAZ的基因表达失衡，MYC2基因表达下调。这可能导致某些下游胁迫响应基因在栓皮栎种子继续脱水时不表达，使种子受到严重伤害。JA和SA能相互串扰，共同调控多种信号网络，进而调控植物防御反应[43]。植物感受到胁迫信号时，SA能和NON-EXPRESSER OF PATHOGENESIS-RELATED GENES 1(NPR1)结合[44]，影响NPR1结构域变化，诱导下游基因转录。NPR1与bZIP类转录因子TGA家族和热休克转录因子HSFA1间均存在互作关系，进而在植物对SA信号的响应中发挥重要调控作用[45-46]。栓皮栎种子脱水过程中，水分流失诱导SA信号通路中TGA相关基因表达失衡，既有上调也有下调，这表明TGA对栓皮栎种子响应脱水中SA信号的调控作用不同。拟南芥II类TGA，即tga2、tga5和tga6突变体中，SA下游的PR1基因表达量比野生型高，其负调控SA下游基因，但该突变体中SA诱导PR1表达能力丧失，植物不能产生SAR[46]。当水分继续流失时，TGA基因表达下调，PR-1基因表达上调，其可能通过调控下游胁迫相关基因，在栓皮栎种子应对严重水分流失过程中发挥作用。因此，我们推测脱水诱导栓皮栎顽拗性种子中JA和SA信号传导途径的变化是其脱水敏感性的潜在调控机制。

综上所述，栓皮栎顽拗性种子脱水过程中，MYB108、MYB1R1、ERF1B、UNE10、CRF4、CML11、JAZ、MYC2和TGA等转录因子及Ca2+、JA和SA信号相关的基因差异表达。它们可能参与栓皮栎种子脱水过程中ROS的产生和清除，以及胁迫信号的传导过程，进而在种子脱水敏感性的调控中发挥作用。本研究为栓皮栎种子脱水敏感性分子调控机制的研究提供理论依据，有助于种子贮藏和种质资源的长期保存，同时也能为其它顽拗性种子的研究提供一定参考。