《Genome》：Molecular network analysis of gene co-expression driving seed oil content in soybean
编辑推荐：
大豆种子成分受蛋白质与油脂含量间负相关关系的调控，这种关系又受到环境因素的复杂影响。在加拿大，西部种植的大豆种子蛋白质含量持续低于东部，这一格局由基因型与环境的互作（genotype by environment interaction, G×E）驱动，进而调
大豆种子成分受蛋白质与油脂含量间负相关关系的调控，这种关系又受到环境因素的复杂影响。在加拿大，西部种植的大豆种子蛋白质含量持续低于东部，这一格局由基因型与环境的互作（genotype by environment interaction, G×E）驱动，进而调节蛋白质-油脂平衡。为揭示这些差异背后的转录网络，研究人员对种植于东部与西部对比环境中的10个大豆基因型的叶片转录组进行了多年加权基因共表达网络分析（weighted gene co-expression network analysis, WGCNA）。分析鉴定出在不同年份及环境中均保守存在的油脂相关共表达模块，其特征基因（eigengene）在每年均与种子油含量呈现显著且一致的相关性。一组核心的重复出现的枢纽基因（hub gene）构成了该网络的骨架，其中许多基因与主要种子成分数量性状位点（quantitative trait locus, QTL）共定位。值得注意的是，一个关键枢纽基因位于第20号染色体上蛋白质-油脂权衡调控因子的邻近区域，而另一些枢纽基因则定位在第10号染色体上与更有利的蛋白质-油脂平衡相关的一个QTL区间内。模块表达与油含量之间的关系表现出基因型与环境的互作效应，为东部与西部区域间表型差异提供了分子基础。这些枢纽转录因子是培育种子成分改良大豆品种的首选候选基因。
大豆种子油含量调控网络：跨环境保守枢纽基因与蛋白-油脂权衡的分子解析
一、研究背景与科学问题
大豆（Glycine max (L.) Merr.）是全球最重要的经济作物之一，其种子兼具高蛋白（约40%干重）与高油（约20%干重）的独特成分特征，使其成为植物油的最大单一来源，同时脱脂豆粕广泛用于动物饲料。然而，大豆育种面临一个根本性挑战：种子蛋白质含量与油含量之间存在广泛记载的负相关关系。这种权衡在代谢层面根源于种子发育期间光合产物前体（主要是碳骨架）的竞争性分配——用于氨基酸合成进而积累蛋白质的碳骨架，无法同时用于脂肪酸从头合成以储存油脂，由此形成生化瓶颈。尽管已鉴定出数百个与大豆油和蛋白质含量相关的QTL，且转录组研究也定位了大量候选基因，但油含量是由众多基因协同作用控制的复杂数量性状，传统聚焦单个基因的研究方法难以捕捉调控代谢通量的复杂调控架构，许多QTL区间内仍含有未表征基因。
此外，种子成分显著受环境因素及基因型与环境互作（G×E）影响。在加拿大存在明显的地理格局：东部（安大略省）生长季较长且温暖，大豆种子蛋白质较高而油含量较低；西部（草原省份）生长季较短且凉爽，种子成分发生偏移，通常报道油含量较高而蛋白质较低，但区域响应因基因型和环境而异。调控种子成分的基因网络对气候条件高度敏感，这要求研究须跨越多个多样化环境以识别稳健的育种遗传靶点。然而，介导加拿大主要农业生态区之间一致种子品质差异的特定转录网络仍知之甚少。为此，研究人员采用系统生物学方法，对多环境、多年份的叶片转录组数据进行WGCNA，以期超越单基因差异表达比较，捕捉油含量这类复杂性状背后的基因协调活动，鉴定在不同环境条件下持续调控油含量的核心基因网络，并寻找能够改良固有蛋白-油权衡的分子靶点。
二、研究设计与技术方法
研究人员选取10个来自成熟度组00（maturity group 00, MG 00）的大豆基因型，这些品系代表连续宽广的种子油含量范围。田间试验于2018、2019和2021三个生长季在加拿大四个地点进行，包括一个东部试验点（安大略省渥太华）和三个西部草原试验点（马尼托巴省布兰登和莫登、萨斯喀彻温省萨斯卡通），采用随机完全区组设计，每个基因型在每个地点-年份组合至少设置三个生物学重复。成熟后测定种子油和蛋白质含量百分比。在R5发育阶段（种子充实起始期）采集幼嫩三出复叶，因叶片作为供应碳、氮及代谢前体的主要源组织，其环境响应性转录程序预期反映影响种子油含量变异的上游调控过程。RNA提取后构建100 bp双端文库，采用Illumina HiSeq 4000平台测序。数据分析流程包括：Trimmomatic和Cutadapt进行质控与接头去除，STAR比对至大豆参考基因组（Wm82.a2.v1），featureCounts生成基因水平计数，DESeq2进行rlog转换和差异表达分析。WGCNA分别构建各年份的加权基因共表达网络（符号网络，最小模块大小30个基因，软阈值幂分别为9、12和9），基于拓扑重叠矩阵（TOM）进行层次聚类和动态树切割鉴定模块，以模块特征基因（eigengene）代表模块表达谱。采用双因素方差分析（two-way ANOVA）检验每个特征基因与区域（东部vs西部）、油含量（连续变量）及区域×油交互项的关联，Benjamini-Hochberg法进行多重检验校正。模块内连通性（kWithin）最高的基因定义为枢纽基因，其功能通过SoyBase和NCBI注释，并与SoyBase QTL记录交叉比对基因组坐标以评估共定位。此外，利用PIPE4生成核心候选基因的蛋白质-蛋白质相互作用（protein-protein interaction, PPI）网络并在Cytoscape中可视化。
三、研究结果
3.1 鉴定到一个保守且受环境调控的油脂网络
研究人员分别构建2018、2019和2021年的年度特异性共表达网络。在2018和2021年网络中，与油含量相关性最强的模块标记为“cyan”（Pearson's r ≈ 0.53?0.68, p < 0.05）；2019年网络中油关联最强的模块标记为“tan”（r ≈ 0.54?0.71, p < 0.05）。油关联模块大小约为70至200个基因。尽管年度网络的模块颜色标签不同，油关联模块表现出重叠的功能谱和共享的高连通基因。2019年tan模块的特征基因分析显示显著的区域效应（p = 2.02 × 10?15）和区域×油交互效应（p = 0.00202），油含量单独并非显著预测因子（p = 0.347），模型解释46.2%的方差。在东部环境中，较高油含量与较低模块特征基因值相关（负斜率），而在西部环境中较高油含量与较高特征基因值相关（正斜率），表明环境背景强烈改变模块表达与种子油含量之间的关联。这一交互效应证明该模块捕获了环境条件依赖的调控机制：在东部，该模块高表达似乎拮抗油积累；在西部，同一模块的高表达则有助于或允许油含量增加。相比之下，2018和2021年的cyan模块在两地均与油含量呈现一致的负相关关系，无显著区域或交互效应（2018: p = 0.0011; 2021: p = 9.34 × 10?7），捕获了稳定、环境无关的油响应程序。2019年区域分离更强可能与当年西部多个地点降水较高有关。
3.2 15个核心枢纽基因构成网络骨架
跨三年共识分析鉴定出15个稳定的“核心”枢纽基因，它们在油关联模块中持续具有高模块内连通性（kWithin）。功能上这些枢纽富集调控和信号蛋白，包括转录因子（MADS-box和ARR2类蛋白）、信号组分（受体样激酶）、转运蛋白（肽转运蛋白PTR3）以及参与膜重塑和脂肪酸周转的酶（α/β-水解酶和PIGX同源物）。这支持油含量变异主要由调控和信号水平的层级控制驱动，而非核心代谢机制本身的改变。
3.3 基因位置分析
15个核心枢纽基因的基因组位置与已知QTL比对显示显著共定位。GLYMA.01G179100（Gm01）位于主要油/蛋白QTL区间内；GLYMA.03G037300（Gm03）位于油含量QTL及脂肪酸组成QTL区间内；Gm10上的串联基因GLYMA.10G092300/092400位于12–38 Mb QTL区间内，与油含量和种子大小QTL重叠；GLYMA.12G156500（Gm12）与油含量QTL及硬脂酸QTL（qSA12-1）共定位。最显著的发现是GLYMA.20G085700位于第20号染色体主要种子成分“热点”区域，紧邻经功能验证的蛋白-油-重量调控因子1（Protein, Oil, Weight Regulator 1, POWR1；GLYMA.20G085100），相距约200–250 kbp。POWR1编码CCT结构域转录因子，其天然截短等位基因以牺牲蛋白质为代价增加油含量和种子重量。GLYMA.20G085700的表达以区域依赖方式与油含量显著相关——西部正相关、东部负相关，且高蛋白亲本Saedanbaek中该基因携带提前终止密码子的推定功能缺失变异。这一证据链表明Gm20热点的表型效应可能并非仅由POWR1单独造成，而可能源于一组连锁调控基因的复合作用。
3.4 第10号染色体与蛋白和油含量相关的QTL
GLYMA.10G092300（肽转运蛋白PTR3）和GLYMA.10G092400（蛋白激酶）这对串联枢纽基因与第10号染色体多性状QTL热点共定位。Lee等人（2019）的大规模全基因组关联分析（genome-wide association study, GWAS）鉴定该位点为蛋白-油权衡的罕见例外——该区域等位基因显著提高种子蛋白含量而不伴随预期的油含量下降。PTR3参与氮再动员和向发育种子的分配，相邻激酶可能调控该转运活性，构成“转运-信号枢纽”，协调富含氮的肽和氨基酸（如天冬酰胺）向发育种子的通量，在维持三酰甘油生物合成所需碳骨架的同时促进蛋白质合成。差异表达分析显示两基因在2018年东部与西部环境间显著差异表达，且在配套多年研究中反复在高蛋白与低蛋白大豆类别间差异表达。
3.5 共表达邻域定义保守功能簇
对15个候选基因各鉴定30个最强共表达伙伴。以受体样激酶RLK2（GLYMA.01G179100）为中心的簇持续与溶血磷脂酰基转移酶LPEAT1和Ca2+信号组分共表达，提示细胞表面环境感知与膜脂动态修饰的直接联系。以ARR2类蛋白（GLYMA.03G017700）为中心的簇包括淀粉合酶4（SS4）和2019年的丙酮酸脱氢酶激酶（PDK），将该枢纽置于碳分配决策点——细胞分裂素信号调控源-库关系，直接调节碳在淀粉储存（SS4）与脂肪酸合成（PDK）之间的分配。以NB-LRR免疫受体（GLYMA.03G037300）为中心的簇连接NADP依赖苹果酸酶（NADP-ME）——产生脂肪酸生物合成所需NADPH的关键酶——以及磷脂酶D（PLD）和ABCG脂质转运蛋白，提示代谢准备与植物防御的协调系统。2019和2021年网络中出现将晚期胚胎丰富蛋白（LEA）与GDSL脂肪酶和α/β-水解酶枢纽（GLYMA.18G127700）连接的簇，体现油代谢与种子成熟末期的整合。GLYMA.20G085700与磷脂酶D和LPEAT1共表达，而拟南芥中LPEAT基因破坏会降低种子脂质含量。
四、讨论与结论
PPI网络分析显示约86个节点和87条边，呈高度模块化架构。B类MADS转录因子GLYMA.04G027200（推定APETALA3/AP3同源物）位于PPI网络中心，具有最高度中心性，与PISTILLATA（PI）同源物、AP1/CAL、AGL24/SVP、SEPALLATA1-4和SEEDSTICK（STK）等形成高序MADS复合物。在拟南芥中，过表达D类MADS基因STK增强脂质合成，而E类SEPALLATA3（SEP3）拮抗STK并降低种子油含量；油棕中EgMADS16过表达降低总脂质含量而EgMADS21激活EgDGAT2增强脂肪酸积累。GLYMA.04G027200属于油关联模块，其表达在三年中均与种子油含量负相关（2018: r = ?0.58; 2019: r = ?0.63; 2021: r = ?0.55; FDR < 0.05），表明该发育调控因子在结构发育与储存脂质合成之间微调碳分配。ARR2类蛋白位于细胞分裂素响应簇中心，与降低大豆细胞分裂素氧化酶活性（GmCKX3）可使种子油增加18%–22%的证据一致。RLK2受体样激酶与油含量负相关，其邻域包括PPR蛋白、葡萄糖-6-磷酸脱氢酶（G6PD1）和磷脂酰基转移酶，整合氧化还原状态与膜重塑。NB-LRR免疫受体与油含量正相关，连接NADP-ME、PDK和START结构域脂质转移蛋白，直接证明免疫信号与脂质代谢在发育种子中紧密交织。此外，GLYMA.14G038200编码大豆胱抑素（soyacystatin）蛋白酶抑制剂，通过限制蛋白水解维持储存蛋白完整性并调节氮再动员；GLYMA.02G260600编码葡萄球菌样核酸酶，与程序性细胞死亡和营养再动员相关。
本研究通过多年系统水平转录组分析，结合WGCNA、遗传定位和PPI数据，鉴定出一个保守且环境响应的大豆种子油含量调控基因网络。15个核心枢纽基因富集调控功能，构成网络骨架。GLYMA.20G085700位于POWR1热点区域，为调控经典蛋白-油权衡提供新候选基因；GLYMA.10G092300/092400构成的转运-信号枢纽位于可解耦该权衡的稀有QTL区域，为同时改良蛋白质和油含量提供有前景的途径。网络表达和连通性受生长环境显著调节，为加拿大东部与西部表型差异提供了分子解释。B类MADS-box转录因子GLYMA.04G027200的核心作用表明发育调控因子充当环境响应性的关键介导者。这些发现超越了单基因分析，为大豆高油品种培育提供了经过验证的高优先级基因靶点，适用于多样化环境适应性育种。
