北沙柳株型候选基因表达量与株型相关性的验证

贺嵘，赵恺，贺玉娇，阿拉腾苏和，王爱君，宁静，韩若霜，孙贵荣，王瑞平，张国盛*

（1.内蒙古农业大学 林学院，内蒙古 呼和浩特 010019；2.鄂尔多斯市造林总场，内蒙古 树林召 014300；3.内蒙古林业科学研究院，内蒙古 呼和浩特 010013；4.内蒙古自治区林业和草原种苗总站，内蒙古 呼和浩特 010010；5.鄂尔多斯市林业和草原局，内蒙古 鄂尔多斯 017010）

植物株型的形成是指植物在生长和发育全过程中，与其有关的组织器官发生的改变。随着植株的生长，侧生分生细胞会不断地生长和分化，并逐渐形成自己的器官。这一过程中，通过一个或多个基因的调节，使得枝条呈现出不同的形态，进而影响植株的株型。植物株型是影响生物量和利用效率的关键性状之一，开展株型定向培育的重要意义，已经在农作物及经济树种育种中得到了充分的证明。植物株型多数是数量遗传性状，由多个彼此独立或协同的基因共同决定的，基因表达发生改变，性状也会变化。植物基因的调控功能是多方位的，综合研究分析对挖掘完善基因功能是十分必要的。李忠南等报道玉米（Zea mays）的分蘖率由4 对主基因与多对微效基因功能决定的［1］。WRKY家族参与调控植物的生长发育等生理过程，拟南芥（Arabidopsis thaliana）中WRKY34与WRKY2共同调控植物花粉的发育［2］。WRKY转录因子也可以通过调控HY5的表达来调控植物下胚轴的生长发育［3］。在水稻（Oryza sativa）理想型形成过程中克隆得到的IPA1其功能获得性突变体表现为分蘖数少、抗倒伏、茎秆粗壮［4］，IPA1不仅可以结合DEP1来调控水稻的株型与产量，还可以增强水稻的抗病性［5］。

北沙柳（Salixpsammophila）是“三北”地区土地荒漠化防治的首选树种［6］，也是木型材、纤维板、刨花板、造纸、柳编、活性炭等产业的重要原料。北沙柳的株型呈现直立型、开散型、中间型等特征。不同的利用方式对北沙柳株型的要求不同，作为木材产业化的原料，理想的株型是通直少（或无）侧枝的直立型；而多分枝开散生长的株型（开散型）对于防沙固土是理想株型。北沙柳纤维短，壁腔比小，细胞壁薄，易打浆成纸，冯文莘等［7］对北沙柳硫酸盐法纸浆造纸进行了可行性研究。周宝［8］发现北沙柳漂白硫酸盐浆和针叶木漂白盐浆混在一起可以增强纸的撕裂指数。聂勋载等［9］通过实验证明了，去皮的北沙柳枝条能生产出品质好的瓦楞纸和纸箱板。北沙柳平茬后更新速度快，原材料价格低，枝条直径小，利于纤维分离，可降低人造板材的成本。内蒙古伊金霍洛旗建立了第一个北沙柳刨花板厂［10］。随着人们对北沙柳的研究的深入，北沙柳的应用范围也越来越广。如何能够通过现代生物技术，结合北沙柳生物学特性，开展快速精准定向化育种，以满足沙产业发展的需要是北沙柳资源精细化利用核心技术之一。

本研究利用定量PCR 技术和转录组测序技术，通过分析预测群体、训练群体、验证群体基因表达量与冠高比的相关性，创建基于相关系数筛选北沙柳株型后备基因方法，为北沙柳株型定向育种提供理论基础，具有一定的现实意义。

1 材料与方法

1.1 试验材料

本研究试验材料全部来自于鄂尔多斯市造林总场沙柳国家林木种质资源库的无性系测定林（2017 年营造）。水培群体均于内蒙古农业大学智能温室进行培育（平均空气湿度为60%RH，平均空气温度为24 ℃，平均土壤湿度为44%RH，平均土壤温度为17 ℃，平均CO2浓度为550 g·m-3）。水培装置主要使用塑料盆等体积较大的容器作为培养容器，但过程中要更换水和营养液、人工供气或在营养液中加入过氧化氢来解决氧气供应问题。

1.2 研究方法

1.2.1 株型特征的调查

2020 年9 月和2021 年9 月分别对无性系测定林的株高、冠幅（东西向和南北向）、地径、萌生枝条数量、叶片特征等指标进行了调查，计算株型特征值冠高比（平均冠幅与高度之比）。2 年计算获得的冠高比呈现相同的趋势，既能体现无性系间的差异，又具有稳定性，且具有正态分布趋势，是体现株型性状较理想的特征值（图1）。

图1 无性系群体冠高比变化趋势Fig.1 Trends of crown-to-height ratio in training groups

1.2.2 预测群体的构建及转录组测序分析

2021 年5 月筛选9 个北沙柳无性系，即直立型、中间型、开散型各设置3 个生物学重复（直立型12-27，11-30，5-25；中间型7-35，17-12，17-13；开散型13-15，16-36，19-49），分别采集枝条，在内蒙古农业大学智能温室中进行水培构成预测群体，水培2个月采集嫩茎。

采用TIANGEN RNAprep Pure Plant Plus Kit试剂盒（天根生化科技有限公司，北京）提取RNA。将预测群体嫩茎的RNA 送到基迪奥生物科技有限公司进行转录组测序。为了保证数据的质量，在使用信息前要对原始数据进行过滤，最终得到clean reads，并将clean reads 进行碱基分析［11］。共检测出40 049 个基因，已知基因有37 865 个，未知基因有2184 个。并利用RSEM 计算每个样本中所有基因的表达量，表达量用FPKM 展示［12］。利用Excel 中CORREL 的平方计算所有基因序列FPKM 值与冠高比的相关系数作为预测值，将预测值按由大到小的顺序进行排列。不同基因的预测值不同，预测值的变化范围是0.0～0.9。

将基因reads count 的数据使用DESeq2［13］进行差异分析：首先对reads count进行标准化，然后计算假设检验的概率，最后通过多重假设检验进行校正，得到FDR 值。基于差异分析的结果，筛选FDR＜ 0.05 且|log2FC|＞ 1 的基因为差异显著基因。为了更清晰的进行数据分析，可使用火山图分析，越靠近两端的基因，差异程度越大。GO 功能分析可以给出基因的GO 功能分类注释。KEGG 是有关Pathway 的主要公共数据库。通过Pathway 显著性富集能确定差异基因参与的主要生化代谢和信号转导途径。蛋白互作关系是通过差异基因蛋白互作网络进行分析的。

本试验是通过直立型与中间型相比，将差异基因进行GO 功能分析与KEGG 通路分析进而筛选北沙柳株型相关候选基因。

1.2.3 目标基因的选择

为了验证目标基因在训练群体中的表达特征，探讨预测值与北沙柳的株型是否相关，主要从以下几个方面选择目标基因。首先选择已在北沙柳中克隆并验证与分枝有关的SpsLAZY1b［14］和SpsTAC2［15］基因。然后将预测值与转录组测序结果相结合，选择与调控株型有关的差异基因ZFP4［16］、TB1［17］、SPA2［18］、ABF2［19］、PYL1［20］。最后选择预测值较大的ATX1［21］、FHY1［22］和RFK1［23］基因。

将目标基因登录号在数据库Phytozome 13（https：//phytozome-next. jgi. doe. gov/）中进行比对，获得对应的桤柳CDs 序列，将CDs 序列在Primer3Plus（https：//www. bioinformatics. nl/cgi-bin/primer3plus/primer3plus.cgi）设计引物（表1），然后由擎科公司（北京）合成。

表1 引物序列Table 1 Primer sequence

1.2.4 训练群体、验证群体的构建

2021 年7 月按照生长旺盛、病虫害少且株型特征有明显变化的标准筛选100 个无性系，标记后组成北沙柳训练群体，并采集嫩茎（液氮低温保存待用）。

将训练群体按照“冠高比的平均值±0.5 倍标准差”分为6 个株型组，每组抽取4 个无性系，共24个无性系作为验证群体。2022 年7月采取枝条于内蒙古农业大学智能温室进行水培，作为水培验证群体；同时2022 年7 月进行野外采样，作为野外验证群体，均每隔10 d 进行一次采样，共采5 次。将采集嫩茎用液氮低温保存。

1.2.5 基因表达量的测定

采用TIANGEN FastKing gDNA Dispelling RT SuperMix 试剂盒（天根生化科技有限公司，北京）将RNA 反转录成cDNA，在-80 ℃冰箱中进行保存。以训练群体与验证群体的cDNA 作为模板，UBQ基因稳定性好，在植物中核苷酸序列的同源性在88% 以上，且袁梦如等人已在北沙柳中选择UBQ作为内参基因进行试验结果良好，所以选择UBQ作为内参基因［24］，用Roche LightCycler 480 Ⅱ仪器进行定量PCR 测定。反应体系为20 μL，程序设置为95 ℃ 5 min；95 ℃ 10 s，60 ℃ 15 s，72 ℃ 15 s，45 个循环；95 ℃ 10 s，60 ℃ 1 min；40 ℃ 10 s；3 个重复，4 ℃保存。

1.2.6 调控株型最优基因组合的筛选

利用SPSS 中秩和分析计算目标基因在株型组中的差异显著性、平均值、标准偏差和极值。选择达到统计显著水平（P＜0.05）的基因，然后将基因进行随机组合统计，组合方式为Cmn（m≤n，其中n为显著水平P＜0.05 的基因数量，并从n 中随机抽取m 个基因）。利用多元回归分析计算每个组合中基因表达量与冠高比的相关性系数。并以“基因个数最少，相关系数最高”为原则，初步筛选最佳与调控株型有关的基因组合。

1.2.7 高相关目标基因蛋白互作网络分析

利用STRING（https：//cn. string-db. org/）对相关系数高的目标基因分别进行蛋白互作分析。通过分析PPI 图中蛋白结构信息、蛋白与蛋白之间的相互作用以及蛋白与蛋白间的连接程度，来预测蛋白质的功能。

2 结果与分析

2.1 转录组测序结果分析

通过Pathway 显著性富集能确定差异基因参与的最主要生化代谢途径和信号转导途径。差异基因GO 功能富集分析得出差异基因主要富集在生物过程和分子功能中。在生物过程中，差异基因在新陈代谢过程、单一生物过程、细胞过程等过程中富集；在分子功能中，差异基因在催化作用、结合性、核酸结合转录因子的活性等功能中富集。且下调的差异基因数量比上调的多。基于差异分析结果，以P＜0.05 且|log2FC|＞1 为基础选择显著差异基因。嫩茎中与调控株型相关的差异基因有CCD、MYB4、WRKY75、RGA1、TCP12、RFK1、SPA、ATX1、ZFP4、PYL1、ABF、TB1等。

2.2 训练群体中目标基因相关系数变化特征

采集训练群体的嫩茎进行定量PCR 试验，验证目标基因表达量与北沙柳冠高比的关系。在嫩茎中，预测群体的相关系数高于训练群体，但基因表达量具有相同的变化趋势。预测值的变化范围为0.04～0.79，训练群体中相关系数的变化范围为0.05～0.60。ATX1的相关系数在两个群体中均相对较高为0.595 1，TAC2相对较低为0.051 7。FHY1、TAC2表达量与冠高比呈正相关，其余都呈现负相关（图2）。

图2 嫩茎中目标基因表达量趋势Fig.2 Expression trends of target genes in the axillaries in the hydroponic and training populations

2.3 验证群体中目标基因相关系数变化特征

野外验证群体与水培验证群体中目标基因表达量与冠高比的相关系数如表2 和表3 所示。在不同采样时间下，LAZY1b、ZFP4、TB1、SPA2、ABF2、PYL1、ATX1、RFK1基因与冠高比成负相关，TAC2、FHY1基因与冠高比成正相关，这与预测群体、训练群体结果一致。预测值较高的ATX1、FHY1、RFK1基因在验证群体中，相关系数仍然高于其它目标基因。

表2 野外验证群体相关系数Table 2 Field validation of population correlation coefficients

表3 水培验证群体相关系数Table 3 Hydroponic validation of population correlation coefficients

2.4 基因组合法

大群体下（训练群体），对目标基因进行秩和分析。结果表明，达到显著差异（P≤0.05）的有ATX1、TB1、ABF2、FHY1、RFK1、PYL1。将差异基因按照的组合方式进行组合。

ATX1、FHY1在6 个株型组中的差异显著性如图3 和图4 所示。6 个株型组的冠高比依次增大，ATX1基因表达量随着冠高比的增大而降低，且1组，2、3、4 组与5、6 组基因表达量之间有显著差异（P＜0.05）；FHY1基因表达量随着冠高比的增大而增大，且在1 组，3、4 组与5、6 组基因表达量有明显差异（P＜0.05）。

图3 ATX1 显著性分析Fig.3 ATX1 significance analysis

图4 FHY1 显著性分析Fig.4 FHY1 significance analysis

2.5 蛋白互作网络分析

将相关系数较高的ATX1、FHY1、RFK1基因进行蛋白互作分析，预测其功能，置信度为0.7（高等），结果如图5所示。

图5 高相关基因蛋白互作网络图（a）ATX1 蛋白（b）FHY1 蛋白Fig.5 Highly relevant gene protein interactions network： （a） ATX1 protein and （b） FHY1 protein

在ATX1蛋白互作网络中，与ATX1直接相关的基因有HMA5、RAN1、COPT5和HMA6。COPT1-6为铜转运蛋白，ATX1、CCH、CCS为铜伴侣蛋白。铜蛋白CCS可以提供铜离子，ATX1和CCH可将铜离子传递给HMA5、RAN1，RAN1再将铜离子整合到乙烯受体里，产生乙烯信号，通过调控乙烯的合成来影响植物的生长。

在FHY1蛋白互作网络中，与FHY1直接相关的是PHYA基因。COP1、PHYA、CRY1基因可以通过调控光的形态建成来影响植物的生长。AUX1为生长素受体，FHL为FHY1互为同源蛋白。PIFs是一种光敏色素蛋白，可以与赤霉素信号通路中的DELLA蛋白相互作用来调控植物的生长发育。

RFK1被报道是与植物育性恢复有关的基因，但与调控植物生长发育方面的基因互作较少，数据支撑较少，但仍具有一定的参考价值，可能在调控北沙柳生长发育方面有着重要作用。

3 讨论

3.1 北沙柳株型与目标基因相关系数的关系

IGT基因家族是与调节分枝角度有关的一类基因家族，包括LAZY、TAC等基因。在水稻中LAZY1基因调节茎向地性反应调控水稻分蘖角度［25］；在拟南芥中，lazy1突变体可增加分枝角度［26］。在桃树（Prunus persica）中，TAC1 突变后可导致分枝角度变小［27］。张磊［14］、袁梦如等［15］在北沙柳中分别克隆了SpsLAZY1b、SpsTAC2基因，并验证了SpsLAZY1b、SpsTAC2与北沙柳的分枝有关，但实验结果表明基因表达量与北沙柳的冠高比相关性较小。这可能是因为SpsLAZY1b、SpsTAC2在北沙柳株型中为非主效基因，是通过调节其他基因间接来影响株型；也可能受内源激素的调控影响基因的表达。

ZFP4、TB1、SPA2、PYL1、ABF2是转录组测序结果中通过GO 分析及KEGG 富集通路分析筛选出的差异基因，同样也是在其他物种中已被证明是与植物株型发育有关的基因。在北沙柳中可找到对应的同源基因。在麻风树（Jatropha curcas）中赤霉素与细胞分裂素共同作用对TB1起着负调控的作用，促进侧芽的萌发［28］；在拟南芥中SPA被认定为红光和远红光受体光敏色素信号传导过程中的下游负调控因子，光照强度或红光、远红光降低时，植物的分枝会受到抑制［29］；秦琳琳发现白桦（Betula platyphylla）中ZFP4基因含有多个脱落酸响应顺式作用元件，通过调控脱落酸信号途径来调控植物的生长［30］。PYL1与ABF2则可以通过调控脱落酸的合成来影响植物株型的发育［31］，这些基因都是通过调控下游基因或其他途径间接影响株型。在北沙柳中，水培群体预测值、训练群体以及验证群体实验结果的趋势一致，在大群体下（训练群体）相关系数在0.0～0.2 范围之间为低相关，但有一定的线性关系，初步预测这些基因在北沙柳中同样可以可调控植物的生长发育。

FHY1、ATX1、RFK1是与北沙柳株型相关性高的基因。在大群体（训练群体）中ATX1相关性最高为0.595 1，RFK1相关性为0.229 3，相关系数在0.2～0.6 范围之间为中等程度的相关。FHY转录因子现已证明是远红光信号通路的重要组成部分，在拟南芥中，FHY转录因子在调节开花时间、枝条分生组织的发育等重要生命过程中都具有非常重要的功能［32］。通过蛋白互作网络分析可知，FHY1可能通过与光敏色素的相互作用来影响北沙柳的株型，也可能在生长素或赤霉素的影响下调节北沙柳的株型；铜是植物生长发育所需的微量元素，ATX1为铜伴侣蛋白，其作用就是供应植物生长所需的金属离子，将过量的金属离子排出从而达到解除重金属毒害的作用。通过蛋白互作网络分析可知，ATX1可能通过转运铜离子调节铜的浓度来影响北沙柳的生长发育，也可能通过调节乙烯信号来影响北沙柳的生长发育；RFK1目前在北沙柳上未见报道，在本研究中嫩茎部分的RFK1基因与北沙柳的株型相关性较高，这说明RFK1基因作为北沙柳株型的后备基因作用是不可忽视的。

3.2 水培试验对北沙柳株型的影响

北沙柳是一种沙生植物，耐干旱。本试验中，将水培后的北沙柳作为预测群体进行转录组测序，筛选目标基因。通过训练群体与验证群体的定量PCR 结果可知，水培群体与野外群体中目标基因与冠高比的相关系数差异不大，具有相同的规律。这说明水培实验不会影响北沙柳的株型。

在水培试验中，水培北沙柳的根系洁白光滑，新生不定根较直，二级根较短，木质化程度低，因为水环境中存在低氧胁迫，激发北沙柳特有的生理生化过程，保证根部缺氧部位获得氧气。虽然生理生化过程有所改变，但通过试验发现，水培群体与野外群体变化规律一致。水培技术，是无土栽培技术中的一种，其原理是运用现代生物技术，在不改变植物遗传物质的基础上，通过外源激素或者环境胁迫等方式激发原有生根基因的表达，促使其形成新根。作为一种新型的农业栽培方式，具有以下几方面的优点：（1）植物长势好，产量高，品质好。水培能提供并稳定维持作物生长所需的光、温、水、气、肥等环境条件，可充分发挥作物的生产潜力。相比土培，水培植株生长速度快、长势强。（2）节省资源。水培系统避免了水资源与营养液的流失以及被土壤中的微生物吸收固定，节约资源的同时又提高了肥水的利用率。水培能摆脱土壤的制约，使植物生长不再依附于土壤的限制，扩大种植范围和空间，节省土地。（3）病虫害少。水培技术是室内栽培方式，属于相对封闭的环境，在一定程度上可避免外界环境的伤害，如土壤微生物、害虫、化学有害物、重金属等对植物的伤害。

3.3 预测群体的构建

本试验是基于转录组测序结果，利用性状与相关基因建立北沙柳株型后备基因的筛选。构建预测群体的目的是为了初步得到植物株型与基因表达量之间的相关性，即以使用较少的样本，获得更准确的数据为基本原则，预测功能基因与性状之间的相关性。

北沙柳开散程度不同使其株型多种多样。本实验选择9 个样本（直立型、中间型、开散型并各设置3 个生物学重复）进行了转录组测序。结果表明，FPKM 值与冠高比的预测值差异非常大，高的可达到0.9 以上，低的在0.000 1 以下。目标基因（SpsLAZY1b、SpsTAC2基因预测值较小，FHY1、ATX1、RFK1基因预测值较大，ZFP4、TB1、SPA2、ABF2、PYL1基因的预测值趋中）在9 个样本中的预测值与在大群体中的相关系数变化趋势一致，预测值大的相关系数也大，预测值小的相关系数也小，且在不同采样时间点下相关性均与预测值一致，验证了转录组测序结果的准确性，证明了9 个北沙柳无性系可以用来构建预测群体。

4 结论

本研究以不同株型北沙柳作为研究材料，构建预测群体、训练群体及验证群体。选择ZFP4、TB1、SPA2、PYL1、ABF2、SpsLAZY1b、SpsTAC2、FHY1、ATX1、RFK1为目标基因。试验证明目标基因表达量与北沙柳冠高比的相关性在预测群体、训练群体以及验证群体中具有相同的趋势，可以利用转录组测序结果筛选北沙柳株型候选基因，基因表达量与目标性状相关系数越高，功能基因与性状的关系越密切。且水培试验不会改变北沙柳的株型。按照“基因个数最少，相关系数最高”的原则，选择（ATX1+FHY1）为鉴别北沙柳株型的最佳基因组合。