《Plants》：Genome-Wide Association Study on Lodging Resistance-Related Traits in Oats
Lijun Zhao,
Rui Yang,
Yantian Deng,
Xiaopeng Zhang,
Lijun Shi,
Bai Du,
Mengya Liu,
Junmei Kang,
Xiao Li and
Tiejun Zhang
编辑推荐：
本研究通过全基因组关联研究 (GWAS) 对130份燕麦种质资源的7个抗倒伏相关性状进行了分析，在两种环境下共鉴定出379个显著关联的数量性状核苷酸 (QTNs)，并注释了54个候选基因，为揭示燕麦抗倒伏遗传机制及分子标记辅助育种提供了重要依据。
1. 引言
燕麦 (Avena sativaL.) 是一种重要的粮饲兼用作物，其抗倒伏性是直接影响产量和品质的关键因素。因此，选育抗倒伏燕麦新品种对于提高作物生产力和经济效益至关重要。倒伏，特别是在抽穗后频繁降雨和大风引发的倒伏，会阻碍机械化收割，降低饲草的产量和品质，已成为限制燕麦实现高产优质的主要制约因素。先前研究表明，倒伏受植物茎秆力学特性及相关性状的遗传基础影响，如株高、茎重、节间长度和茎壁厚度。茎秆穿刺强度 (SPS) 和茎秆折断强度 (SBS) 是直接衡量茎秆力学特性的关键指标。全基因组关联研究 (GWAS) 已成为解析复杂性状遗传机制的有效工具，在作物分子育种早期阶段被广泛用于识别与性状相关的标记。尽管已有研究对126份欧洲燕麦种质进行了抗倒伏性的GWAS，但燕麦中抗倒伏相关性状的更多表型评估以及持续鉴定抗倒伏位点的工作仍然必要，且抗倒伏相关基因尚未在燕麦中得到充分表征。本研究以130份燕麦种质为材料，在乳熟期对植物高度 (PH)、单茎鲜重 (FWSS)、基部第二节间长度 (LBSI)、基部第二节间直径 (DBSI)、基部第二节间壁厚 (WTBSI)、茎秆折断强度 (SBS) 和茎秆穿刺强度 (SPS) 这7个抗倒伏相关性状进行了表型评估，并利用全基因组重测序产生的高质量单核苷酸多态性 (SNPs) 进行了GWAS分析，旨在鉴定与抗倒伏性相关的遗传位点和候选基因，为理解燕麦抗倒伏的分子机制和选育抗倒伏新品种提供理论依据。
2. 结果
2.1. 表型数据分析
在两种环境下对130份燕麦种质资源的7个抗倒伏相关性状进行了评估。结果显示，这7个性状均存在显著变异，所有性状的变异系数 (CV) 均超过10%。其中，SPS的CV最低，而FWSS的CV最高。SBS在两个试验点均表现出较大的变异，CV分别为36.61%和33.53%。这些结果表明，抗倒伏相关性状在该群体中具有丰富的遗传多样性。抗倒伏相关性状的遗传力分别为0.43、0.85、0.39、0.91、0.34、0.90和0.70，表明这些性状具有较高的遗传力。此外，频率分布结果显示，所有7个性状均呈现连续变异和近似正态分布，符合数量性状的遗传模式。
2.2. 表型相关性分析
各性状间的相关性分析结果显示，燕麦茎秆抗倒伏相关性状如FWSS、PH、DBSI、SBS和SPS之间存在密切相关性。PH、DBSI、WTBSI、SBS与FWSS之间的相关系数为0.42–0.90，呈显著正相关。DBSI和WTBSI均与SBS和SPS呈正相关，相关系数为0.47–0.90。相比之下，LBSI与WTBSI (r = -0.35) 和SPS (r = -0.23) 分别呈显著负相关。这些结果强调了形态和机械性状在决定倒伏风险方面的相互关联性。
2.3. SNP分布与群体结构
为了进一步阐明性状变异的遗传基础，对130份燕麦种质资源的SNP多样性进行了分析。共检测到36,928,068个SNPs，覆盖燕麦基因组的全部21条染色体。其中，染色体4C包含的SNP标记数量最多，而染色体6D包含的SNP最少。在预先定义的模型参数范围 (K = 1 到 6) 内，交叉验证 (CV) 误差随着K值的增加先减小后增大。在K=3时观察到最小的CV误差 (0.419)，表明推断的群体结构在模型拟合和复杂性之间达到了最佳平衡。因此，130份燕麦材料被清晰地分为三个遗传上不同的亚群。故选择K=3用于后续的GWAS分析。
2.4. 抗倒伏性状的全基因组关联分析
在两个环境中，GWAS共鉴定出379个显著的QTNs，其中与WTBSI显著相关的数量最多 (94个)，与DBSI显著相关的数量最少 (18个)。与PH相关的38个显著QTNs分布在16条染色体上，其中6A号染色体上最多。与FWSS相关的34个QTNs位于12条染色体上，主要集中在4C染色体。与SBS相关的78个QTNs分布在19条染色体上，而与SPS相关的66个QTNs分布在18条染色体上，在7A (8个) 和7C (13个) 染色体上频率最高。鉴定出55个与LBSI相关的QTNs，分布在18条染色体上，主要在4A染色体 (9个)，而18个与DBSI相关的QTNs分布在12条染色体上。大多数与WTBSI相关的QTNs位于4A (16个) 和4D (55个) 染色体上。
在DT试验点，鉴定出与这7个性状相关的160个显著QTNs，分布在所有21条染色体上。QTNs数量最多的是4C染色体 (22个)，而3A染色体仅含1个。具体而言，26个QTNs与PH显著相关，20个与FWSS相关，48个与SBS相关，20个与SPS相关，18个与LBSI相关，10个与DBSI相关，18个与WTBSI相关。在ZJK试验点，GWAS鉴定出与这些性状相关的223个显著QTNs，其中在4D染色体上检测到的数量最多 (62个)。其中，12个QTNs与PH相关，16个与FWSS相关，30个与SBS相关，46个与SPS相关，37个与LBSI相关，8个与DBSI相关，76个与WTBSI相关。
进一步对两种环境进行比较分析，发现了4个在两种环境中一致检测到的QTNs。其中3个共同的位点与FWSS相关，位于5D染色体 (chr5D_90334559) 和7D染色体 (chr7D_504853535, chr7D_504853578) 上。一个与SBS相关的共享位点位于4C染色体 (chr4C_507724758) 上。这些重叠的位点可能代表了在不同环境中对抗倒伏相关性状起稳定和关键作用的关键遗传区域。
2.5. 候选基因分析与鉴定筛选
对于Manhattan图中呈现尖锐峰值的性状相关位点，在显著QTN上下游20 kb区间内对基因进行注释，最终鉴定出54个不同的候选基因。与PH相关的候选基因有5个，包括肽酶家族成员、细胞质氨肽酶家族、质体乙酰辅酶A羧化酶 (ACC2)、UDP-糖基转移酶家族和丝氨酸/苏氨酸蛋白磷酸酶。对于FWSS，在两种环境中注释了5个候选基因，包括编码MYM型蛋白1-like和核酸内切酶的基因，尽管在ZJK位点有一个注释基因功能不明确。
为SBS鉴定出7个候选基因，功能涉及磷酸核糖基转移酶、水解酶、萜类合酶、核糖核酸酶H和抗坏血酸铁依赖性氧化还原酶活性。为SPS注释了11个候选基因，与WRKY30转录因子、E3泛素连接酶、抗坏血酸铁依赖性氧化还原酶和蛋白激酶活性相关。有14个候选基因与LBSI相关，包括编码抗坏血酸铁依赖性氧化还原酶、五肽重复结构域、转运蛋白、超氧化物歧化酶、YABB蛋白和谷胱甘肽S-转移酶的基因。在两种环境中均未鉴定出与DBSI相关的注释基因。
在两个位点共有12个候选基因与WTBSI相关。功能注释显示，这些基因主要编码肽酶、水解酶、具有核糖核酸酶III活性的蛋白和3‘-核酸外切酶。这些共注释的基因可能通过调节细胞壁代谢途径，在决定WTBSI方面发挥关键的调控作用。
3. 讨论
3.1. 抗倒伏性状的表型分析
本研究中，对两种环境下呈正态分布的抗倒伏相关性状进行统计分析，结果显示表型变异率超过10%，并表明关联群体内这些性状存在相当大的多样性。相关性分析进一步揭示了PH、DBSI、WTBSI、SBS和FWSS之间存在密切的正相关关系 (r = 0.42–0.90)，强调了形态和机械性状在决定倒伏风险方面的相互关联性。例如，DBSI和WTBSI与SBS和SPS的显著正相关支持了先前的研究发现，即增加第二节间的茎秆直径和壁厚对于提高茎秆机械强度至关重要。值得注意的是，虽然PH被广泛认为是导致倒伏易感性的关键因素，但在本研究的燕麦材料中，其与SBS和直径相关性状 (DBSI, WTBSI) 的强正相关表明，较高的品系可能通过增强茎秆健壮性来补偿潜在的倒伏风险。这与小麦的研究形成对比，在小麦中株高和抗倒伏性通常呈负相关，凸显了这些性状关系的物种特异性差异。总之，这些结果表明，与茎秆粗度和机械强度相关的性状是降低燕麦倒伏敏感性的基础。
3.2. 全基因组关联分析
群体结构是影响复杂性状GWAS结果可靠性的关键因素。准确评估群体结构可以有效校正由遗传背景异质性引起的假阳性关联。在本研究中，通过10折交叉验证确定的最佳模型 (K=3) 有效地最小了此类干扰，从而确保了鉴定出的379个显著关联QTNs的可靠性。该分析也为准确识别与抗倒伏性相关的稳定位点提供了重要基础。
针对7个抗倒伏相关性状的GWAS，我们采用了BLINK模型，该模型是Farm CPU的增强版本，因为其在检测微效位点方面具有更强的统计效能，并且与常规GWAS模型相比具有更高的计算效率。该分析鉴定出379个显著QTNs，分布在全基因组21条染色体上，与目标抗倒伏性状相关。值得注意的是，我们的结果与先前的研究一致，该研究对饲用燕麦的7个农艺性状进行了基于Farm CPU的GWAS，鉴定出94个显著位点和155个候选基因。该研究进一步揭示了D亚基因组在关键生长性状中的主导作用，C亚基因组对叶片和穗部性状的显著贡献，以及A亚基因组仅在少数性状中有限的参与。类似地，我们的研究发现D亚基因组上抗倒伏相关SNP的数量 (173个) 显著高于A (106个) 和C (104个) 亚基因组。从进化角度看，D亚基因组是在燕麦多倍化过程中较晚整合的，可能保留了更多与环境适应性相关的基因。这一发现增强了我们对燕麦抗倒伏位点染色体分布的理解，并为有针对性地鉴定抗性基因，特别是在D和C亚基因组上，提供了有价值的方向。
此外，目前针对燕麦抗倒伏性状的GWAS研究仍然相对较少。先前的研究主要集中在株高、穗长和茎长等农艺性状上，而针对茎秆直径、茎壁厚度和茎秆抗弯强度等抗倒伏相关性状的GWAS主要在玉米和水稻中进行。例如，在水稻中，研究者鉴定出16个与茎壁厚度相关的数量性状位点 (QTL)。我们检测到94个与WTBSI相关的数量性状核苷酸 (QTNs)，表明燕麦的茎壁厚度受更多微效位点控制，这与该性状的数量性质一致。在玉米中，研究者利用多位点GWAS鉴定了423个与茎秆直径和抗弯强度相关的QTNs，而我们的研究分别鉴定出18个和78个与DBSI和SBS相关的QTNs。在燕麦中鉴定出的QTNs数量较少，可能反映了这种六倍体物种抗倒伏性更复杂的遗传调控。
在多种环境中一致检测到的QTNs通常被认为是稳定位点。本研究中，四个共定位的QTNs (chr5D_90334559, chr7D_504853535, chr7D_504853578 和 chr4C_507724758) 在两种环境中被重复鉴定出来。这些环境稳定的QTNs很可能代表了参与茎秆抗倒伏的关键位点，因此是进一步基因克隆和功能验证的有希望的候选者。
3.3. 候选基因功能分析
在本研究中，我们对与燕麦抗倒伏性状相关的379个QTNs两侧20 kb侧翼区域内的候选基因进行了筛选。结合基因功能注释，最终鉴定出54个具有调节燕麦抗倒伏性潜在作用的候选基因；其中一些先前被报道参与细胞壁生物合成、木质素沉积和激素信号传导等关键生物学过程，使其成为优先考虑的候选基因。PH与UDP糖基转移酶家族相关。该家族基因被证明通过介导苯丙烷生物合成途径来调节玉米中的木质素沉积，其表达水平与茎秆的维管束面积和表皮细胞厚度显著相关。此外，FWSS与一个MYM型锌指蛋白相关基因有关。先前研究证实，锌指MYM型蛋白1-like异构体X1调节节间数和穗重，从而间接提高抗倒伏性。同时，SPS与泛素相关机制和WRKY30转录因子均相关。E3泛素连接酶GmILPA1通过调节赤霉素代谢途径影响茎秆发育和强度。WRKY30直接调控关键的木质素单体生物合成基因，并参与木质素生物合成的遗传调控。LBSI与一个YABB蛋白家族相关基因有关。在水稻中，同源基因OsYABBY4负调控赤霉素信号通路，从而抑制节间伸长并优化茎秆强度。WTBSI与一个GT31家族蛋白相关基因有关，该蛋白已被证明定位于高尔基体，并介导拟南芥中初生和次生细胞壁纤维素的生物合成，这些基因的突变会导致细胞扩张缺陷和生长迟缓。这些基因被认为是燕麦抗倒伏相关性状可靠的候选基因，对其功能的进一步验证将有助于进一步阐明产量相关性状的潜在遗传和分子机制。同时，我们注意到一个新发现，即基部第二节间直径 (DBSI) 没有注释到候选基因，表明该性状可能受非编码区域调控。
尽管这些发现增强了我们对抗倒伏性遗传结构的理解，但与其他作物相比，燕麦候选基因的功能注释仍然不足。此外，本研究中的候选基因筛选依赖于参考基因组注释，燕麦基因组组装的完整性可能会限制功能预测的准确性。20 kb的窗口大小是一个经验近似值；基于群体特异性连锁不平衡 (LD) 衰减估计和实验验证 (如表达分析或功能测定) 的进一步优化，对于验证这些候选基因的作用至关重要。因此，需要包括功能验证和共表达网络分析在内的更全面的研究，以阐明茎秆倒伏调控的分子基础，并为燕麦抗倒伏性的遗传改良提供理论基础。
4. 材料与方法
4.1. 植物材料与生长条件
从国家作物种质资源库、内蒙古农业大学苜蓿研究所和中国农业科学院畜牧研究所选取130份燕麦种质资源作为材料。该群体包括115份种质材料和15个改良品种。试验点分别设在山西省大同市 (DT) 和河北省张家口市 (ZJK)。在DT点，燕麦种质于2023年4月24日播种。在ZJK点，于2023年5月16日播种。两个试验点均采用随机区组设计，三次重复。每个材料单行种植，行长1.5米，行距60厘米，株距5厘米，每个材料30株。人工播种，深度3-5厘米。在整个生长季节根据需要施用适当的农药管理病虫害和杂草。
4.2. 表型数据收集与分析
乳熟期是指燕麦开花后至成熟前的时期。此阶段植株重心最高，穗部重量最大，是评估抗倒伏性的关键时期。因此，在乳熟期测量了PH、FWSS、LBSI、DBSI、WTBSI、SBS和SPS等倒伏相关表型性状。从每个材料的每个重复中随机选择5株代表性植株，测量PH、FWSS、LBSI、DBSI、WTBSI、SBS和SPS。计算每个性状各品系的平均值。植物高度 (PH) 是在自然状态下用尺子从植株基部测量到穗尖的高度。基部第二节间 (LBSI) 是从地面向上数的第二个节间，用尺子测量其长度。基部第二节间直径 (DBSI) 在中点用数显卡尺测量。随后，小心地将节间纵向劈开，在同一位置直接测量基部第二节间壁厚 (WTBSI)。随机采样植物的主茎，在地面处切割，样品包含茎、叶、叶鞘和穗。立即使用电子天平测定单茎鲜重 (FWSS)。茎秆穿刺强度 (SPS) 和茎秆折断强度 (SBS) 均使用茎秆强度测定仪 (YYD-1；浙江托普仪器有限公司，中国杭州) 测量，使用长度均匀的节间片段。在茎秆基部，将第二和第三节间中部以恒定速度垂直插入茎秆方向，读取刺穿茎秆表皮的最大值。同样，茎秆折断强度也是在茎秆中心以均匀速度按压，记录该数值。
使用Microsoft Excel 2019整理和预处理表型数据。使用IBM SPSS Statistics 26.0软件计算每个表型指标的平均值、峰度、偏度、标准差 (SD) 和变异系数 (CV = SD/均值)，并生成频率分布图。使用Origin 2024软件进行相关性分析。根据Nyquist公式计算广义遗传力 (H2)：H2= δ2G/(δ2G+ δ2E/n)，其中δ2G和δ2E分别为遗传方差和残差方差。
4.3. SNP calling与群体结构
从田间等量采集新鲜叶片组织。样品立即放入预冷的离心管中，密封，并迅速在液氮中冷冻。使用植物基因组DNA提取试剂盒提取基因组DNA。
然后使用高质量DNA样本在Illumina NovaSeq 6000平台上进行基因组文库构建和高通量测序，每个基因型产生约59.4 GB的原始测序数据。
在原始数据预处理过程中，首先使用Trimmomatic v0.39进行严格的质量控制，去除低质量reads、测序接头和引物序列，获得高质量的clean data。随后，使用BWA-MEM v0.7.17软件将过滤后的reads比对到燕麦参考基因组，并使用GATK Haplotype Caller (v4.2.3.062) 进行SNP calling。同时，使用SAMtools (v1.13) mpileup结合VarScan (v2.4.6) 进行额外的SNP calling流程，以提高变异检测的准确性。最后，使用VCFtools (v0.1.16) 过滤缺失率 > 10%、次要等位基因频率 (MAF) < 0.05或最低测序深度 < 5的SNPs，得到36,928,068个高质量SNPs，分布在21条燕麦染色体上，用于后续GWAS分析。
使用“CM plot” R (v4.1.3) 包生成全基因组SNP密度图，以可视化SNP在染色体上的分布。使用Admixture软件 (v1.3.0) 进行群体结构分析，K值范围为1到6。最佳K值基于最小交叉验证 (CV) 误差确定。
4.4. 全基因组关联研究
使用GAPIT3 (v3.1.0) 包中的贝叶斯信息和连锁不平衡迭代嵌套关键途径 (BLINK) 模型进行全基因组关联分析。该方法增强了统计效能和计算速度。应用-log10(P) ≥ 6 (LOD ≥ 6) 的显著性阈值来识别显著关联信号。使用R包“CM plot”可视化GWAS结果，可直观展示关联位点的基因组分布和统计显著性。生成曼哈顿图和分位数-分位数 (Q-Q) 图来描述这些GWAS结果。
4.5. 候选基因鉴定
基于上述GWAS结果，筛选在两种环境中均超过显著性阈值的标记。然后对这些标记进行统计分析，提取与目标性状相关的显著QTNs信息。根据燕麦参考基因组 (OT3098) 注释文件，在每个显著QTN上下游20 kb范围内 (总共40 kb) 鉴定候选基因。该筛选区间是通过参考牧草作物全基因组关联研究中广泛采用的策略，并充分考虑燕麦基因组连锁不平衡 (LD) 衰减特征而建立的，这确保了有效覆盖与显著SNPs存在LD的功能基因，同时避免纳入过多不相关基因而增加后续功能验证的难度。随后使用EggNOG-mapper在线工具对这些基因进行功能注释，从而鉴定出可能参与燕麦抗倒伏的候选基因。
5. 结论
本研究利用36,928,068个高质量SNPs，采用BLINK模型对130份饲用燕麦种质资源进行了全基因组关联分析，在两种环境下共筛选出379个显著关联的QTNs。随后的候选基因筛选和功能注释获得了54个与燕麦抗倒伏相关的基因。功能注释表明，这些基因可能通过共同调控细胞壁代谢、茎秆发育、机械强度和生物量分配等途径来增强抗倒伏性。这些发现为燕麦抗倒伏性的遗传改良提供了宝贵的遗传资源和高效的分子标记，从而有助于选育新的抗倒伏品种。
