大豆种子中矿质元素含量控制位点的鉴定

背景

矿质营养素在生物系统的生化和生理功能中起着至关重要的作用。通过遗传改良提高种子矿物质含量被认为是满足饲料需求的最有希望和最具成本效益的方法。本研究的总体目标是对137个大豆品系的种子进行矿物质含量(Ca, K, P和S)的GWAS，这些品系代表了加拿大种植的早熟大豆(成熟度组000-II)的多样性。

结果

137个大豆品系在5种环境中生长，并使用便携式x射线荧光光谱仪测量种子矿物质含量。关联分析使用三种统计模型和从测序基因分型和全基因组测序组合数据集中获得的一组220万个snp进行。在使用的三种统计模型中(在两种环境中)，至少有两种模型鉴定出8个与Ca、K、P和S含量显著相关的qtl，每个qtl对表型变异的贡献率从17%到31%不等。这8个qtl中有7个在其他3个环境中具有很强的再现性。总共鉴定出三个候选基因参与这些矿物元素的运输和同化。

结论

很少有GWAS研究鉴定与大豆种子矿质元素含量相关的qtl。此外，本研究鉴定的qtl和候选基因将为大豆营养品质的遗传改良提供有益的标记辅助选择。此外，本研究还提供了在加拿大大豆种质中遇到的表型变异范围的详细信息。

大豆被广泛用于食品、饲料和工业用途，使其成为最通用的粮食作物之一。事实上，大豆是蛋白质、油脂和碳水化合物的重要来源，也是其他有益营养物质的重要来源，例如矿质元素，这些元素会影响油脂和蛋白质部分的最终用途特性，以及用于种植后续作物的种子的品质特性[1]。植物获得矿物质营养是一个非常动态和复杂的过程，受生物和非生物因素及其相互作用的影响[2]。在农业方面，据报道，缺乏必需元素可导致产量损失、疾病易感性增加、代谢受损、正常发育中断和种子质量差[3.]。为此，了解矿物元素在各种环境条件下的吸收、调节、运输和储存对于解读植物与其环境之间的复杂关系至关重要。因此，种子基因组谱是将植物的遗传特征与其对环境扰动的反应相匹配的有力工具[4]。

通过基因改良提高种子矿质营养素含量被认为是确保满足消费者饮食需求的最有希望和最具成本效益的方法。由于任何性状的育种都依赖于表型变异的存在，这就需要鉴定出具有有用的谷物矿物遗传变异的品种，并了解这些种子性状的遗传结构[5]。

自从植物营养的概念被建立以来，人们一直在努力开发定量测量生物体元素组成的方法和工具。6]。与火焰原子吸收光谱(FAAS)等湿化学方法相比，能量色散x射线荧光(ED-XRF)等光谱方法操作更容易，而且是一种非破坏性的工具，因此越来越受到人们的青睐[7]。近年来，ED-XRF已成功地用于测定大豆中Ca、K、P和S的浓度[j]。8，9]，在可可[7]和豌豆种子[10]。在这些研究中，通过光谱学获得的测量结果与以前使用更常见但昂贵的分析方法进行的研究一致。此外，这些研究表明，植物中矿物质元素的含量可以跨越很大的范围，而这种含量是由遗传和环境因素共同决定的[2，4，11，12，13]。

数量有限的研究旨在确定控制种子中矿物元素积累的遗传结构，并得出结论认为这些性状很可能由许多基因控制[2，14]。Zhang等。[15[j]利用148个SSR标记和178个F / s标记，报道了4个与大豆种子钙含量相关的qtl(分别位于第7、8和20染色体上)_2:3和157 F_2:4行。King等。[13]利用916个SSR标记和92个F标记，报道了3个与磷含量相关的qtl(分别位于第7、12和17染色体上)_2:4行。最近，Ramamurthy等人。[2]利用1536个单核苷酸多态性(SNP)标记和288个大豆重组自交系(RILs)，鉴定出7个与Ca、K和S含量相关的qtl(位于4、6、15、16和18号染色体上)。使用GWAS方法，Ziegler等。[4]利用美国农业部大豆种质资源收集的1653份大豆材料和36489个snp，鉴定出大豆种子中与Ca、K、P、S含量相关的9个qtl(位于染色体1、2、5、9、10和13上)。最后，Dhanapal等。[16]报道了与大豆茎部Ca、K、P和S含量相关的20条染色体上的65个qtl，使用了31748个snp和104个大豆基因型。总的来说，在这些研究中发现的qtl几乎没有重叠。这可能是由于大豆种子中矿质元素含量的遗传决定因素在不同的种质间不同，或者一些研究受到基因组覆盖不足的影响，未能检测到共享的qtl。

在控制大豆种子矿质元素积累的qtl鉴定不完整且经常不一致的背景下[4]，特别是在早熟大豆(MG000-II)中，我们试图表征一组137个加拿大短季大豆品种的表型多样性，并鉴定该组种质中控制Ca, K, P和S含量的qtl。利用大量SNP标记和三种分析方法(CMLM、MLMM和FarmCPU)，我们共鉴定出32个控制这四个重要元素积累的qtl，其中8个qtl通过至少两种方法共同鉴定。相信本研究结果将为今后大豆种子品质和营养成分的遗传改良研究提供新的思路。

湿化学与能量色散x射线荧光法的相关性

为了验证我们选择的分析方法(能量色散x射线荧光，ED-XRF)，用ED-XRF和火焰原子吸收光谱(FAAS)对30个样品进行了Ca和K的分析，并用分光光度法对P含量进行了分析。如图3所示。1，两种方法的相关系数(r)均为正且极显著(P< 0.001)，范围为0.91 (Ca) ~ 0.94 (P)。结果表明，ED-XDF方法适用于大豆种子中Ca、P、K含量的定量测定。

表型变异及性状间的相关性

采用ED-XRF装置测定了2013年2个试验点(2个重复/试验点)137个大豆品系上Ca、K、P和S的浓度。频率分布表现为近似正态分布，并且在数量上似乎是遗传的(图2)。2)。如表所示1在干重基础上，4种元素的种子矿物质含量变化范围为:Ca 1.6 ~ 2.4 mg/g, K 17 ~ 21 mg/g, P 4.5 ~ 6.5 mg/g, S 3.5 ~ 5.5 mg/g。在所有137个品系中，Ca、K、P和S的平均值分别为1.8、18.7、5.3和4.3 mg/g。2个基因型平均值之间的最不显著差异(LSD)为Ca含量0.03 mg/g, K含量0.44 mg/g, P含量0.28 mg/g, S含量0.09 mg/g。广义遗传率在81% (K) ~ 99% (S)之间，具有较大的表型变异和较高的遗传率，表明这些性状和关联面板可以很好地揭示这些性状的遗传结构。

如表所示2，方差分析显示，基因型和环境都有极显著的影响(P≤0.001)对除Ca外的所有性状的表型变异均有显著影响，其中基因型效应是变异的唯一显著来源。在任何性状中均未观察到显著的x基因型环境相互作用。观察到的表型值显著(p< 0.001)在两个试验点之间相关，相关性范围在0.75到0.98之间。不同矿物质中的种子含量也被证明是相关的(表S1，粗体)。这些两两比较均有统计学意义(p< 0.05)， K与S的相关性最高(r²= 0.67,p< 0.001)以及P和K之间(r²= 0.65,p< 0.001)。

基因分型和SNP呼叫

关联小组的细胞系最初通过GBS方法进行基因分型，共产生56 K高质量snp。在第二步中，使用包含4.3万个snp的参考面板对原始gbs衍生snp集进行缺失位点插入。去除InDels后，MAF < 0.05，杂合度> 0.1的标记，共保留2.18 M个SNP，平均标记密度为每435个碱基1个SNP。这些2.18 M snp在大豆20条染色体上的物理分布如图所示。S1。由此获得的基因型数据随后用于表征该面板内的群体结构并寻找标记-性状关联。

人口结构

这137个加拿大大豆核心系的群体结构最初是通过fastSTRUCTURE和亚群数(k)为7(图。3)。此外，从图中可以看出。3b和c系统发育树和基于pca的种群结构分析结果与fastSTRUCTURE分析结果一致。总之，这些结果表明k= 7可以很好地评价种群结构，采用相应的Q矩阵进行GWAS。

大豆种子矿质元素含量的全基因组关联扫描

为了发现导致表型变异的染色体区域，我们使用了三种分析工具来测量标记-性状关联:FarmCPU、CMLM和MLMM。如分位数-分位数(QQ)图所示(图2)。S2)，三个模型都成功地限制了观察到的混淆效应p -值只偏离对角线(预期)p-值)在最极端的值(超过3E-03的几乎所有性状)。

这些关联分析的结果显示为图中FarmCPU、CMLM和MLMM的曼哈顿图。4。根据错误发现率阈值(蓝色水平线，FDR≤0.05)，我们检测到32个qtl，其中7个与Ca含量相关，10个与K含量相关，5个与P含量相关，10个与S含量相关(表S2)。有趣的是，有一个共同的QTL对K (K_#1)和P (P_#1)都有贡献。的不正确p -Ca含量为1.35E-06 ~ 2.84E-21, K含量为1.89E-05 ~ 8.05E-19, P含量为1.17E-06 ~ 3.61E-12, S含量为1.75E-05 ~ 6.63E-15。

在这32个qtl中，共有8个qtl被至少两个模型共同识别(图2)。5)，这8个稳健qtl的特征总结于表中3.。表型变异部分解释(R²)， Ca为20% ~ 21%，K为17% ~ 31%，P为22%，S为18% ~ 23%。Ca、K、P和S的等位基因效应大小分别为0.06 ~ 0.07 mg/g、0.30 ~ 0.57 mg/g、0.30 mg/g和0.15 ~ 0.46 mg/g。遗传变异(加性)解释为狭义遗传力(h²Ca为41%，K为82%，P为78%，S为93%。

八个共鉴定的QTL在三种环境中的验证

为了验证至少两个模型检测到的8个qtl中的每个qtl的稳定性，从另外三个试验中获得数据。总的来说，在三个新环境中，七个qtl在至少两个环境中得到了验证(图2)。S3和表S3)。只有Ca的QTL#4无法在三种新环境中的任何一种中验证。I_18环境的验证率最低，在该环境中成功检测到5个qtl(图2)。6)。在24个可能的qtl -环境组合(8个qtl × 3个环境组合)中，有18个在对照SNP峰值的株系平均表型之间存在显著差异。这些结果表明，所鉴定的qtl在广泛的环境中具有鲁棒性。

改进共鉴定QTL的GWA扫描

为了更深入地探索这些健壮qtl的变异，我们从2.2万个snp的大目录中提取了七个最健壮qtl周围的单倍型区域内的所有snp。将这些数据与修剪后的数据集(243 K)合并，再次使用三个模型进行GWAS。在这7个例子中的6个，观察到更强的关联信号，前一个和新的峰值SNP之间的物理距离在1到311 kb之间(表S4)，但总是位于同一单倍型块内。

稳健QTL区域内候选基因的预测

基于GWAS结果，我们研究了大豆基因组中标注的基因，以便从与每个性状显著相关的位点中鉴定出可能的候选基因。为了建立候选基因列表，我们只关注那些位于由最左和最右侧标记所划分的区域内的基因，这些标记处于完美LD (D ' = 1)，具有上述七个qtl的峰值SNP。这些基因组区域(大小在32到360 kb之间)是从Wm82.a2中提取的。我们检测了v1和这些区域内基因的GO注释(表1)4)。这种方法的一个例子如图2所示。7。每个区域中(全部或部分)存在的基因数量在4到43个之间，这些基因的完整列表及其注释见表S5。

为了确定候选基因，我们寻找符合以下两个标准之一的基因:1)被注释为参与给定矿物元素的运输并在根、芽或叶中表达的基因;2)被注释为参与感兴趣元素的摄取、易位和/或平衡并主要在种子中表达的基因。总共确定了三个有希望的候选基因，这些基因参与了这些矿物元素的运输或同化。我们首先发现Glyma.06G046000(位于Ca_#3 SNP峰上游132 kb)， Glyma.10G020000(位于K_#3 SNP峰下游222 kb)。这两个基因都被注释为参与运输，并在根尖和根毛中表达。此外，Glyma.06G046000在嫩叶、花、主根、豆荚和种子中均有表达(图2)。S4B和d)。最后，Glyma.20G151500 (S_#10下游32 kb)被注释为参与硫酸盐同化，并在花、根、根瘤和种子中表达(图2)。S4f). qtl K_#1、P_#1、S_#4和S_#7未找到符合标准的候选基因。

候选基因的结构和核苷酸变异及其预测的功能影响

为了确定候选基因内或重叠的遗传(结构或核苷酸)变异是否可能构成因果变异，我们检查了56个品系的全基因组测序数据中建立的此类变异目录。没有发现与这三个候选基因完全或部分重叠的结构变异(> 51 bp)。在核苷酸变异方面，在两个基因的编码区共发现18个snp(1个在Glyma.06G046000中，17个在Glyma.20G151500中)。所有这些变异都被预测为对蛋白质功能具有“修饰”或“低”影响。因此，观察到的表型变异不太可能是由于这些候选基因的功能丧失。

为了更深入地了解这些候选基因在观察到的表型变异中的作用，进行了单倍型分析。如图1所示。S5显著的表型差异(P各候选基因的单倍型间差异均< 0.05)。例如，在Glyma.06G046000中鉴定的4个单倍型(A、B、C和D)中，携带B单倍型的7个品系与携带A、C或D单倍型的130个品系相比，Ca含量显著不同。同样，Glyma.10G020000中携带C单倍型的5个品系，其K含量高于携带其他4个单倍型中的一个的品系。在Glyma.20G151500的例子中，它是更常见的单倍型A (n= 99)， S含量显著高于其他4个单倍型。这些结果支持每个候选基因都很有希望，因为它们都对表型变异有贡献。

表型变异及性状间的相关性

在最初发现标记-性状关联的两种环境中，所有元素(Ca、K、P和S)的种子含量均呈正态分布，表明它们是定量遗传的。在干重基础上，本研究的表型变异范围为Ca含量1.6 ~ 2.4 mg/g, K含量17 ~ 21 mg/g, P含量4.5 ~ 6.5 mg/g, S含量3.5 ~ 5.5 mg/g。在以前的研究中，不同的范围被报道过。Otaka等。[8]和Homura等人。[17Ca(分别为1.5 ~ 3.5 mg/g和1.5 ~ 3.2 mg/g)和K(分别为19.8 ~ 22.1 mg/g和16.7 ~ 21.2 mg/g)的种子含量相似。同样，King等人。[13]报告的P含量值范围相当(3.9至5.6 mg/g)，而Dhanapal等人[16]报告的数值范围较低(1.3至4.9 mg/g)。然而，本研究中S含量的范围略高于以往报道。我们的S含量值普遍高于Fageria [18]和Bellaloui等人。[19(分别为2.5至4.0毫克/克和1.6至3.1毫克/克)。总体而言，本文报告的值与以往研究报告的值相当一致，范围的细微差异可能只是反映了在不同环境中生长的特定植物组以及测量方法的选择所造成的差异。尽管如此，在关联研究的背景下，表型值的准确性不如充分捕捉值在整个小组中的变化那么重要。

在我们的研究中，两种环境的广义遗传力估计都很高(H²> 0.80)。这种相对较高的广义遗传力表明，表型在很大程度上是由基因型效应决定的[20.]。在以前的研究中已经报道了这些性状的类似遗传力(H²= 0.48 ~ 0.93) [14，16]。此外，K和P浓度之间存在显著的正相关关系(0.65，表S1)。这与先前报道的0.80至0.94的结果一致[2，19]。

大豆种子矿质元素含量的全基因组关联扫描

共鉴定出32个与Ca、K、P和S含量相关的qtl(表S2)，其中8个qtl被至少两个模型检测到。更重要的是，这8个qtl中有7个被证明是非常健壮的，因为它们可以在另外3个试验中成功地确认为影响矿物元素含量。在之前的GWAS研究中，Ziegler等。[4]和Ning等。[14]分别报道了22个和9个qtl。使用链接映射方法，Ramamurthy等。[2]报道7和Bellaloui等人。[19]报道了11个与大豆种子钙、钾、磷、硫含量相关的qtl。在这项工作中检测到的大量初始qtl(32个)可能归因于更详尽的基因组覆盖(2.2 M snp)以及使用多种模型来检测标记-性状关联。这不太可能是由于我们的关联小组中表型值的范围特别大，因为这个范围与以前研究报告的范围相当。

本研究中描述的qtl通常解释了相当大一部分表型变异(18-31%)。先前报道的qtl解释的表型变异从2%到18%不等。有趣的是，我们发现其中一个qtl对钾和磷含量都有影响。考虑到钾和磷元素之间的高度正相关，这种钾和磷含量的共享QTL并不奇怪。同样，Dhanapal等人。[16]也报道了一个与这两个元素的内容相关的QTL。这种相关性可能是由于共同的生理机制和代谢途径[4，21]。换句话说，它可能是通过控制这些矿物质浓度的同一基因的多效性发生的，如共转运蛋白[20.或者仅仅是因为在同一区域存在独立的基因。

为了将我们的结果与之前确定的QTL区域进行比较，我们在SoyBase基因组浏览器上查询了之前的QTL，并定义了它们的物理位置。本研究中鉴定的7个强效QTL均与先前报道的基于家族的定位或GWAS中鉴定的QTL间隔不一致。因此，当前的qtl可以被认为是新颖的。通过这项工作确定的qtl与之前报道的qtl之间没有重叠，这可能反映了这些性状是由不同作图实验中使用的实验材料中的不同基因决定的。

矿物元素积累的候选基因及其功能

如上所述，我们关注的是在根、芽或叶中也表达的运输相关基因，以及主要在种子中表达的参与营养吸收、转运和/或平衡的基因。对于钾含量，在Gm 10上的QTL K_#3上发现了一个与转运相关的基因Glyma.10G020000。该基因被标注为K⁺钾转运蛋白及其同源物答:芥（AT4G13420.1)编码高亲和力的K⁺转运体5 (HAK5)。在水稻方面，Yang等。[22研究表明，HAK5在控制K的流入和K的迁移中都起着重要作用⁺进入根部并将其输送到植物的空中部分。我们的候选基因的两个类似基因(Glyma.02 g154100和Glyma.07 g042500)被功能表征为参与根对K的吸收⁺在大豆中[23，24]。有趣的是，在Dhanapal等人的研究中。[16]， K含量的qtl确实与这些HAK5相似物的基因组位置重叠。在我们的研究中，我们没有发现这些平行体与钾含量之间的任何关联，这表明是10号染色体(Glyma.10G020000)上HAK5平行体的变异导致了加拿大早熟大豆品系钾积累的差异。

对于S含量，我们发现了一个编码ATP硫化酶1的基因(Glyma.20G151500)，这是已知的第一个参与硫酸盐同化途径的酶拟南芥(AT3G22890.1;ASA1) [25]。一个一个芥互补脱氧核糖核酸编码ASA1成功地补充了酿酒酵母ATP硫化酶突变体(met3)，从而恢复蛋氨酸异养和硫酸盐运输[26]。有趣的是，我们的候选基因在Gm10上被[16(Glyma.10 g242600)。这再次表明，相同的酶活性有助于S积累，但不同的基因拷贝控制不同的种质资源S含量。

最后，对于Ca含量，候选基因Glyma.06G046000被注释为参与钙运输的钙转运atp酶。这个基因与an是同源的答:芥轨迹(AT1G27770.1;ACA1)编码叶绿体包膜Ca2⁺- atp酶，已知与钙调蛋白结合，导致Ca2激活⁺泵(27]。研究表明Ca2⁺- atp酶是真核细胞中主动运输Ca2+的酶[28]并参与植物生命周期的所有阶段，包括生长和发育[29]。

与以往的研究相比，本研究中使用的高密度标记有助于重复检测与大豆种子中矿质元素含量相关的几个新位点。除了提供有关加拿大大豆种质中种子矿质元素含量的表型变异范围的详细信息外，本研究还提供了有关其积累的遗传结构的更多信息。本研究鉴定的标记和基因将为大豆通过标记辅助选择进行遗传改良提供参考。

植物材料和实验设计

根据Sonah等人的群体结构分析，从304份更大的群体中选择了137份早熟大豆品系(成熟度组000-II, MG000-II) [30.来代表加拿大短季大豆的遗传多样性。大豆品系来源于博士。Louise O 'Donoughue (CÉROM, St-Mathieu-de-Beloeil, QC)， Elroy Cober(加拿大农业和农业食品部，渥太华，安大略省)，Istvan Rajcan(圭尔夫大学，圭尔夫，安大略省)和Éric Gagnon先生(Semences program Inc.， st - c辽阔，QC)。在第一阶段(发现)，品系在两种环境中表型化，即2013年在加拿大的Woodstock (ON)和St-Mathieu-de-Beloeil (QC)。试验设计为广义格型，所有品系种植在单行地块上，每个位置有两个重复。在第二(验证)阶段，2017年(17)和2018年(18)在渥太华(ON)的中央实验农场的三种环境中种植了相同的品系。这些品系以改良的增强型设计种植，作为四行地块，有一个重复。在每一年内，采用两种不同的处理:不灌溉(N)或滴灌(I)。由于N_17试验没有完整的品系/种子，因此仅使用I_17、I_18和N_18试验的数据对发现的qtl进行鲁棒性检验。

校正及验证

能量色散x射线荧光(ED-XRF)光谱仪的校准采用经验校准方法[31]，其中使用了与感兴趣的样品具有相似成分和形态的标准物。在这里，使用火焰原子吸收光谱法(FAAS)测量了一组样品中Ca和K的元素浓度。32]和分光光度法测定磷含量，按[33]。对于S的含量，我们连续加入Na₂所以₄参考物质(WEPAL, IPE 885 (Maize))由美国国家标准与技术研究所(NIST)提供，其S含量已知。由此获得的值可作为校准ED-XRF装置(Niton XL3t955 GOLDD)标准曲线的基线。为了确定我们的ED-XRF测量的准确性和可靠性，用ED-XRF和湿化学方法分析了30个种子样品中Ca、K和p的浓度，而不需要湿化学验证S的含量。

表型和统计分析

对于每个样品，使用研磨机(Foss a /S: Cyclotec™1093 sample Mill)研磨10g整粒种子。从每条生产线上取0.3 g的均匀细粉样品，在液压制粒机(Carver 4350.L)中使用不锈钢制粒模具压制，以生产紧凑的13毫米(~ 0.2毫米厚)的颗粒。待用ED-XRF测定钙、钾、磷、硫含量时，将微球保存。

利用R软件包«lmer»进行描述性统计、基因型变异、环境和环境效应基因型以及这些矿物质含量之间的相关性分析。为了结合来自不同环境的信息，使用META-R中的限制最大似然来计算最佳线性无偏预测(blup) [34]。广义的遗传性H²跨环境计算如下:

在哪里\ ({\ upsigma} _ {\ mathrm {g}} ^ 2, {\ upsigma} _ {\ mathrm {g} \ mathrm {e}} ^ 2 \)和\({\upsigma}_{\ mathm {e}}^2 \)分别为基因型、基因型与环境交互作用和误差方差分量。nEnv是环境的数量，nrep是复制的数量。

基因分型和snp归算

来自测序192 plex GBS文库的总计约2.03亿个100 bp Illumina HiSeq2000单端reads可用于137个系(如先前详细介绍]30.])。简而言之，限制性内切酶猿K1用于生成单个192 plex GBS文库(包含额外的不相关样本)，该文库在Illumina HiSeq2000测序仪的单线上测序。137个种群共获得2.03亿个100-bp单端reads。Fast-GBS管道[35]和Wm82.a2。V1参考基因组[36]用于最小读取深度为2次的SNP调用，并删除丢失数据超过80%的位点。使用BEAGLE v5对这组gbs衍生的snp进行缺失基因型的第一步插入[37]。其次，使用从102个部分重叠(56个共享)系的全基因组重测序中获得的4.3万个snp的参考面板，对gbs衍生的snp集进行缺失位点插入[38]。以前曾评估过这种未分型基因座的输入准确性[38的比例为96.4%。

人口结构与亲属关系分析

对于种群结构分析，我们使用了修剪后的(r)²> 0.5)使用Plink v1.9获得的14 K snp目录[39]。程序fastSTRUCTURE中实现的贝叶斯模型[40]用于人口分层分析。亚种群的数量(k)设为1 ~ 12，独立迭代3次。亚种群的数量(k)是使用fastSTRUCTURE中实现的python脚本(" choosek ")确定的。除了fastSTRUCTURE之外，还使用了两种不同的方法来推断种群结构:(i)基于2000次迭代的Tamura-Nei模型计算的共识系统发育树，在MEGA7中实现[41]及(ii)在GAPIT计划中实施的主成分分析(PCA) [42]。

全基因组关联分析

在GAPIT中，利用243个K snp的修剪目录评估了标记与表型之间的全基因组关联²> 0.9)和各性状的BLUP值。单位点GWAS采用CMLM，多位点GWAS采用MLMM和FarmCPU。通过亲属矩阵(K)和通过fastSTRUCTURE估计的群体结构矩阵(Q)传递的系间遗传亲缘关系用于控制假阳性关联。三种模型(CMLM、MLMM和FarmCPU)的标记-性状关联显著性阈值进行了调整p-值，错误发现率(FDR)设为≤0.05，根据Benjamini和Hochberg程序[43]。我们假设所有显著的标记-性状关联标记在同一单倍型区(由最左和最右的标记所划分的区域，这些标记处于完美LD (D ' = 1)， SNP峰值)。当不同模型检测到的峰值不同时，选择两种模型检测到的峰值。此外，当至少两个模型检测到相同的峰值SNP时，最低的未校正p-value报告。我们选择仅报告和研究共检测的qtl，即每个性状至少有两个模型检测到的qtl。为了鉴定候选基因，将SNP峰值周围的感兴趣基因组区域定义为在与SNP峰值处于完美LD (D ' = 1)的最左和最右侧标记之间延伸。

三种环境下共检测qtl的等位基因效应验证

为了评估在发现阶段鉴定的qtl的可重复性，我们使用同一组137个品系，验证了qtl在三种不同环境(I_17、I_18和N_18)中的等位基因效应。标记-性状相关性采用t检验。根据SNP峰值的等位基因类别将群体分为两组。我们对两组的平均表型值进行t检验。标记-性状关联的显著性阈值根据多个检验进行调整(α = 0.05/n，其中n为一个性状的qtl总数)。显著性检验采用R 3.5.1版本中实现的t检验函数，根据[44]: Y = μ + f(标记)+误差，其中Y等于性状值，μ等于总体均值，f(标记)是显著标记的函数。

候选基因及其功能分析

通过使用数据挖掘算法[45]，从SoyBase浏览器中提取了所有存在于感兴趣的单倍型块中的基因，并检查了它们的GO注释。在确定候选基因后，使用电子荧光象形文字(eFP)浏览器(www.bar.utoronto.ca)。

候选基因列表中由潜在功能缺失(LOF)等位基因引起的转录本改变通过检查[46]。LOF分析仅使用全基因组测序数据集(56个品系)。对于基因区域内核苷酸变异的研究，SnpEff [47]与全套核苷酸变异(修剪前，2.18 M snp)一起使用。

为了进一步支持候选基因参与观察到的表型变异，采用基因中心单倍型(GCH)方法，使用HaplotypeMiner识别候选基因内定义单倍型的多态性标记[48]和2.18个snp的完整目录。然后对每个单倍型的平均表型值进行t检验。

本研究中使用和分析的数据集可应通讯作者的合理要求向其提供。

BLUP:

最好的线性无偏预测器

CMLM:

压缩混合线性模型

FarmCPU:

固定与随机模型循环概率统一

GWAS:

全基因组关联研究

H²：

广义遗传

h²：

狭义遗传

LD:

连锁不平衡

LOF:

功能丧失

加:

次要等位基因频率

MLMM:

多轨迹混合线性模型

QTL:

数量性状位点

R²：

表型变异解释

SNP:

单核苷酸多态性

作者要感谢来自加拿大魁省拉瓦尔大学Richard Belanger实验室和franois Belzile实验室的Jean Martin和Martin Lacroix。

本研究由马里粮食安全农业培训(FASAM in French)通过合作研究与发展赠款资助。资助者在研究设计、分析或手稿准备中没有任何作用。

从属关系

1. 加拿大魁北克省魁北克市拉瓦尔大学农业与食品科学学院和综合与系统生物学研究所(IBIS)植物学硕士

   Sidiki Malle和franois Belzile

2. 加拿大农业和农业食品渥太华研究与发展中心，加拿大安大略省渥太华

   马尔科姆·莫里森

贡献

S.M.和F.B.构思了这个项目，解释了数据并撰写了手稿。S.M.进行了大部分表型分析和所有数据分析。M.M.进行了用于验证qtl的试验(三个地点)。所有作者已阅读并同意稿件。

相应的作者

对应到弗朗索瓦Belzile。

伦理批准并同意参与

不适用。

发表同意书

不适用。

相互竞争的利益

作者声明不存在利益冲突。

出版商的注意

施普林格·自然对已出版的地图和机构关系中的管辖权要求保持中立。

开放获取本文遵循知识共享署名4.0国际许可协议，该协议允许以任何媒介或格式使用、共享、改编、分发和复制，只要您适当地注明原作者和来源，提供知识共享许可协议的链接，并注明是否进行了更改。本文中的图像或其他第三方材料包含在文章的知识共享许可协议中，除非在材料的署名中另有说明。如果材料未包含在文章的知识共享许可中，并且您的预期用途不被法律法规允许或超过允许的用途，您将需要直接获得版权所有者的许可。如欲查阅本许可证副本，请浏览http://creativecommons.org/licenses/by/4.0/。创作共用公共领域免责声明(http://creativecommons.org/publicdomain/zero/1.0/)适用于本文中提供的数据，除非在数据的信用额度中另有说明。

转载及权限