整合转录组与共表达网络分析揭示玉米苞叶数关键模块与枢纽基因

《Genes》：Transcriptome and Co-Expression Network Analysis Reveals Key Modules and Hub Genes Associated with Husk Number in Maize

　　

摘要 背景：苞叶数是影响玉米籽粒脱水速率和产量的关键农艺性状，但其调控机制仍 largely unknown。本研究旨在通过整合转录组测序、层次聚类和加权基因共表达网络分析（WGCNA），鉴定与苞叶数相关的基因模块和枢纽基因。方法：研究人员对两个苞叶数极端差异的玉米自交系（LJ311，约20片苞叶；PHG35，约5片苞叶）在四个发育阶段（V7–V10）的穗组织进行了RNA-seq分析。共鉴定出12,605个假定苞叶相关基因。结合差异表达分析、层次聚类和WGCNA对候选基因进行优先排序。结果：层次聚类发现了一个晚期特异表达的基因簇（Cluster 3，44个基因），在高苞叶系中于V8–V10阶段显著上调。WGCNA鉴定出三个与苞叶数显著相关的模块：MEblue和MEbrown（正相关）以及MEturquoise（负相关）。Cluster 3富集在三羧酸（TCA）循环、DNA复制和甘油磷脂代谢通路中。将枢纽基因（|kME| > 0.9）与高置信度差异表达基因（DEGs）（|log2FoldChange| > 2）取交集，得到两组优先候选基因：上调基因富集于细胞周期、油菜素内酯（BR）信号和磷脂代谢；下调基因富集于丝裂原活化蛋白激酶（MAPK）信号、糖酵解和细胞壁降解。结论：这些结果表明苞叶数可能受一个潜在的BR–MAPK–细胞周期调控轴调节，为未来功能验证和分子育种提供了调控网络和候选基因。

研究背景与问题

玉米是全球最重要的粮食作物之一，其机械化收获是实现全程机械化的最后瓶颈。当前中国大部分玉米品种在收获时的籽粒含水量高达30%–40%，远高于机械化收获适宜标准的15%–25%，导致收获效率低、干燥储存成本增加、籽粒破损率高，并在温暖潮湿条件下易引发霉变、真菌感染和穗发芽等问题。因此，提高后期籽粒脱水速率、降低收获期含水量成为现代玉米育种的核心目标之一。苞叶作为包裹玉米果穗的变态叶，其形态特征显著影响籽粒脱水速率：苞叶面积与收获期籽粒含水量呈显著正相关，苞叶层数少、叶片薄、包裹松散及早衰等性状均有助于加速脱水。然而，尽管已有数量性状位点定位和全基因组关联分析研究鉴定了若干与苞叶数相关的位点和候选基因，但由于该性状具有复杂的数量遗传特征，其全局基因调控网络和核心功能模块仍未阐明。为此，研究人员利用转录组测序结合计算生物学方法，以期系统解析苞叶数决定的分子基础。

研究内容与结论

研究人员选取了两个苞叶数极端差异的玉米自交系——LJ311（约20片苞叶）和PHG35（约5片苞叶），在其四个关键发育阶段（V7、V8、V9、V10）采集穗组织进行RNA-seq分析。通过整合差异表达分析、无监督层次聚类和加权基因共表达网络分析（WGCNA），系统鉴定了与苞叶数变异相关的共表达模块和枢纽基因，构建了涵盖潜在促进和抑制机制的双向调控模型。研究揭示了苞叶数可能受一个BR–MAPK–细胞周期调控轴调节的分子网络，为理解苞叶数的分子基础提供了新见解，并为培育适于机械化收获的玉米品种提供了有价值的基因资源。

主要技术方法

1. 1.
   转录组测序与数据处理：对LJ311和PHG35两个自交系在V7至V10四个发育阶段的穗组织（包含苞叶）进行RNA-seq，每个基因型每阶段设置三个生物学重复。原始数据经fastp质控后，使用HISAT2比对至玉米B73参考基因组（RefGen_v5），StringTie定量基因表达水平。通过与已发表的玉米穗时空表达数据集比较，去除主要在穗内部结构表达的基因，获得12,605个假定苞叶相关基因用于后续分析。
2. 2.
   差异表达分析：采用DESeq2建立~ Stage + Genotype模型，以|log₂FoldChange| > 1和padj < 0.05为标准鉴定差异表达基因（DEGs）。同时拟合含交互项的模型进行似然比检验，鉴定基因型×阶段互作效应显著的基因。
3. 3.
   无监督层次聚类：从12,605个基因中选取方差最大的前200个基因，基于1-Pearson相关系数为距离、complete linkage法进行层次聚类，将基因划分为三个簇（Cluster 1–3）。
4. 4.
   加权基因共表达网络分析（WGCNA）：选取方差最大的前8000个基因构建带符号共表达网络，软阈值设为10。采用blockwiseModules函数识别模块，计算模块特征基因（ME）并与苞叶数和基因型进行Pearson相关性分析。定义|kME| > 0.9且|GS| > 0.7的基因为严格枢纽基因。
5. 5.
   候选基因筛选与验证：将正相关模块（蓝色和棕色）的枢纽基因与上调DEGs（|log₂FoldChange| > 2）取交集得到上调候选基因；将负相关模块（青绿色）的枢纽基因与下调DEGs取交集得到下调候选基因。选取5个代表性候选基因进行qRT-PCR验证。

研究结果

2.1 苞叶数表型差异与全局转录组特征：LJ311成熟苞叶数平均为22.56±1.61，PHG35为5.34±0.43。主成分分析显示第一主成分（PC1）解释48.3%方差，将两个基因型明显分开；第二主成分（PC2）解释12.3%方差，反映发育阶段效应。同一基因型-阶段组合的三个生物学重复紧密聚集，组内相关系数均高于0.95，表明转录组数据质量高，适合后续定量分析。

2.2 差异表达基因鉴定：采用阶段校正模型共鉴定到3764个DEGs，其中1900个在LJ311中上调，1864个下调。阶段特异性比较显示各阶段既有共享也有特有DEGs。交互项分析发现1947个基因表现出显著的基因型×阶段互作效应，其中39个Cluster 3基因（88.6%）呈现显著互作，表明其晚期上调是基因型-阶段依赖性的而非简单的基因型主效应。

2.3 无监督层次聚类鉴定晚期特异表达基因簇：将方差最大的前200个基因聚为三簇。Cluster 3（44个基因）表现出独特的表达模式：在LJ311中于V8–V10阶段强烈特异性上调，而在PHG35及LJ311的V7阶段表达极低。这种晚期特异性激活提示Cluster 3可能是一个功能协调的基因集，与苞叶数决定相关。

2.4 WGCNA鉴定与苞叶数显著相关的共表达模块：共获得十个共表达模块。蓝色模块（1371个基因，r=0.88，p<0.001）和棕色模块（1183个基因，r=0.85，p<0.001）与苞叶数和LJ311基因型正相关；青绿色模块（1963个基因，r=-0.97，p<0.001）与苞叶数强负相关。蓝色和棕色模块在LJ311中持续高表达，蓝色模块在V8–V10阶段进一步升高；青绿色模块主要在PHG35中表达。

2.5 层次聚类与WGCNA结果比较：Cluster 3与蓝色模块仅有7个基因重叠，与棕色模块无重叠，在单个基因层面的一致性有限。但Cluster 3所有44个基因的功能富集分析显示其显著富集于TCA循环、氧化磷酸化、DNA复制与修复、甘油磷脂代谢及类胡萝卜素/脱落酸代谢等通路，指向能量生产、细胞增殖和膜生物合成功能。

2.6 枢纽基因鉴定与优先候选基因精炼：在蓝色、棕色和青绿色模块中分别鉴定到123、231和282个严格枢纽基因。将正相关模块枢纽基因与上调DEGs（|log₂FoldChange| > 2）取交集得到324个上调优先候选基因；将青绿色模块枢纽基因与下调DEGs取交集得到266个下调优先候选基因。上调候选基因富集于DNA复制、错配修复、同源重组、油菜素内酯（BR）生物合成、植物MAPK信号通路、核糖体和甘油磷脂代谢等与细胞周期进程、激素信号和膜生物合成相关的通路。下调候选基因富集于糖酵解、果糖和甘露糖代谢、木聚糖降解、MAPK信号和磷脂酰肌醇信号系统等通路。

2.7 最终候选基因的表达模式与初步验证：热图清晰地将候选基因分为两大簇，完美匹配其调控方向。五株代表性候选基因的qRT-PCR验证结果与RNA-seq的倍数变化高度一致，证实了表达模式的可靠性。

讨论与结论

讨论部分指出，单一方法难以捕获复杂性状的完整调控网络，而整合差异表达分析、层次聚类和WGCNA的多策略方法可提高候选基因发现效率。Cluster 3的晚期特异性激活恰逢苞叶原基快速分化和伸长的关键时期，其富集的TCA循环、氧化磷酸化、DNA复制和甘油磷脂代谢共同构成支持苞叶后期发育的协调分子网络。三个关键模块的功能分化表明，正调控可能需要细胞增殖机器与细胞骨架-激素信号系统的协同作用，而负调控模块中内质网应激的富集提示快速苞叶生长造成的蛋白质合成压力可能通过未折叠蛋白反应限制器官发育。上调候选基因富集于BR信号和细胞周期，下调候选基因富集于MAPK信号和细胞壁降解，共同定义了苞叶数促进与抑制之间的分子边界，提示可能存在一个“早期决定-晚期执行”的调控模型。研究也指出了局限性：仅使用两个自交系导致基因型和苞叶数完全混杂；采样组织为整个幼穗；样本量较小；qRT-PCR仅验证表达模式未提供功能证据；部分表达差异可能反映发育不同步性。

研究结论：本研究通过整合差异表达分析、无监督层次聚类和WGCNA，系统探索了玉米苞叶数变异的分子基础。鉴定出三个与苞叶数显著相关的共表达模块（蓝色、棕色和青绿色），分别代表生长促进和生长抑制网络。层次聚类发现一个小的独特基因簇（Cluster 3）在高苞叶系LJ311的晚期发育阶段（V8–V10）强烈激活，富集于TCA循环、氧化磷酸化、DNA复制和膜脂代谢。通过将模块枢纽基因与高幅度差异表达基因取交集，获得一组优先候选基因。上调候选基因富集于细胞周期、油菜素内酯（BR）信号和磷脂代谢，而下调候选基因富集于MAPK信号、糖酵解和细胞壁多糖降解。综上所述，这些结果表明苞叶数可能与一个汇聚于BR–MAPK–细胞周期轴的潜在双向调控程序相关。本工作鉴定的候选基因为进一步的功能验证和分子育种提供了有价值的目标。