多倍体生物能源柳枝稷气候适应的基因组机制gydF4y2Ba

柳枝稷(gydF4y2BaP。gydF4y2BavirgatumgydF4y2Ba)是一种很有前途的生物燃料作物，也是北美高草草原的重要组成部分。从历史上看，高草草原是地球上最大的温带生物群落之一，它们仍然是大气碳的重要汇gydF4y2Ba^{7gydF4y2Ba，gydF4y2Ba8gydF4y2Ba}。然而，大多数现存的天然柳枝稷种群都局限于“遗迹”遗址，这些遗址代表着对未来高草草原的保护和繁殖至关重要但正在减少的遗传资源。gydF4y2Ba

生物质生产是柳枝稷作为饲料和生物能源作物的主要育种目标gydF4y2Ba^9gydF4y2Ba是种子生产和进化适应性的有力代表gydF4y2Ba^10gydF4y2Ba。自从美国能源部将柳枝稷命名为典型的草本生物燃料原料以来，生物质产量试验已经证明了柳枝稷生物能源生产的经济可行性，并且已经培育出的品种大大超过玉米和其他纤维素原料gydF4y2Ba^11gydF4y2Ba。然而，单个品种往往只在一个狭窄的气候生态位上生产。因此，为了实现收益最大化，柳枝稷育种和生物技术应重点发展气候基因型匹配gydF4y2Ba^{12gydF4y2Ba，gydF4y2Ba13gydF4y2Ba}通过鉴定种群生物量积累和气候适应的基因组基础。这将提振未来的收益率gydF4y2Ba^14gydF4y2Ba并将柳枝稷作为一种经济上和环境上可持续的生物能源产品。gydF4y2Ba

尽管大量的数量遗传变异是气候相关的胁迫耐受性和生物量生产的基础gydF4y2Ba^{15gydF4y2Ba，gydF4y2Ba16gydF4y2Ba}然而，先前柳枝稷基因组序列的碎片化和不完整性阻碍了候选基因的发现和其他分子育种工作。AP13柳枝稷基因型基因组大(单倍体基因组大小= 1129.9兆碱基(Mb))，重复(56.9%重复)(图2)。gydF4y2Ba1gydF4y2Ba扩展数据图gydF4y2Ba1gydF4y2Ba)和多倍体。与其他一些异种杂交物种，如玉米(以自交系B73参考基因组为代表)相比，AP13是近交的。其基因组在自然异交群体范围内保持了相当水平的杂合性(扩展数据图)。gydF4y2Ba1gydF4y2Ba）.尽管存在这种复杂性，我们的深度PacBio长读测序加上深度短读抛光和细菌人工染色体(BAC)克隆验证产生了高度连续的“v5”AP13基因组组装(扩展数据图)。gydF4y2Ba1gydF4y2Ba;数据可从Phytozome网站获得gydF4y2Bahttps://phytozome-next.jgi.doe.govgydF4y2Ba）.我们修剪了得到的大组合(NgydF4y2Ba_50gydF4y2Ba= 5.5 Mb)为单个代表性单倍型，然后使用两个高密度遗传图谱的共识定位和排序为染色体假分子(补充数据gydF4y2Ba1gydF4y2Ba）.染色体通过遗传距离分配到亚基因组gydF4y2Ba黍rudgeiigydF4y2Ba^17gydF4y2Ba的K亚基因组的姊妹分类群gydF4y2Bap . virgatumgydF4y2Ba)，并通过从头开始重复聚类。最终版本仅包含0.4%的漏洞，与2016年的v4版本相比减少了75倍。gydF4y2Bahttps://phytozome-next.jgi.doe.gov/info/Pvirgatum_v4_1gydF4y2Ba）.重要的是，基因组组装与三个遗传信息来源是共线性的，尽管它是独立于这三个来源组装的:一个近二倍体亲戚(gydF4y2Ba黍halliigydF4y2Ba)，伪f的标记顺序gydF4y2Ba_2gydF4y2Ba替代亚基因组的遗传图谱和基因顺序(图2)。gydF4y2Ba1gydF4y2Ba扩展数据图gydF4y2Ba1gydF4y2Ba、补充数据gydF4y2Ba2gydF4y2Ba）.这些共线性表明我们已经为每个亚基因组开发了一个单倍体组装和注释。gydF4y2Ba

至关重要的是，我们能够区分两个亚基因组之间的基因和重复序列。基因注释来自Illumina RNA测序(gydF4y2BangydF4y2Ba_库gydF4y2Ba= 88,gydF4y2BangydF4y2Ba_{条件gydF4y2Ba}= 18, bbb30亿reads)和PacBio Iso-Seq (gydF4y2BangydF4y2Ba_{条件gydF4y2Ba}= 9， > 450万次读取，补充数据gydF4y2Ba3.gydF4y2Ba) -包含80,278个初级转录本和49,664个替代转录本，并且与基因组组装一样完整(BUSCO = 99.4%)(扩展数据图。gydF4y2Ba1gydF4y2Ba）.我们利用这些注释建立了多个序列比对和时间尺度的系统发育树，将亚基因组祖先物种的分化追溯到大约670万年前(Ma)。反转录转座子亚基因组特异性增殖的长端重复序列分析确定了柳枝稷形成多倍体事件的上界≤4.6 Ma(图2)。gydF4y2Ba1 bgydF4y2Ba)，这表明四倍体柳枝稷出现于上新世或更新世早期的冰期-间冰期旋回。gydF4y2Ba

虽然有两个生殖隔离gydF4y2Ba^18gydF4y2Ba柳枝稷的细胞型有四倍体(4x)和八倍体(8x)，四倍体占多数gydF4y2Ba^19gydF4y2Ba而且比八倍体分布的地理范围更广gydF4y2Ba^20.gydF4y2Ba。为了研究气候适应、胁迫耐受性和生物量生产的遗传基础，我们建立了732个纯四倍体基因型的多样性面板(补充数据)gydF4y2Ba4gydF4y2Ba）.我们无性繁殖和移植到10个普通花园，这些花园横跨1862公里的纬度，从德克萨斯州南部到南达科他州(美国)。gydF4y2BangydF4y2Ba_{植物gydF4y2Ba}= 5,521)。gydF4y2Ba2gydF4y2Ba)，并通过深度(中位数= 59×)覆盖2 × 150-bp对端无pcr的Illumina文库对每个基因型进行重测序。重要的是，重测序覆盖范围没有偏向于任何一个亚基因组(似然比检验)gydF4y2BaχgydF4y2Ba^2gydF4y2Ba= 1.32，自由度= 1，gydF4y2BaPgydF4y2Ba= 0.25)。我们的重测序获得了3380万个单核苷酸多态性(snp)(次要等位基因频率≥0.5%)与基因组对应。我们还重新组装了这些深度重测序文库的252个基因型子集，并将其称为存在-缺失和结构变异(例如，100-1,500-bp的插入和缺失)。为了将性状和分子变异与气候联系起来，我们提取了46个气候变量gydF4y2Ba^{21gydF4y2Ba，gydF4y2Ba22gydF4y2Ba}从每个基因型的地理参考收集位置，并将这些数据聚类为七个组，这些组解释了多样性面板上的大部分气候变化(扩展数据图)。gydF4y2Ba2gydF4y2Ba）.gydF4y2Ba

柳枝稷的气候相关适应先前被假设为强调北部高地和南部低地生态型之间的差异，并以不同的叶片和整个植物形态为例gydF4y2Ba^{13gydF4y2Ba，gydF4y2Ba23gydF4y2Ba，gydF4y2Ba24gydF4y2Ba，gydF4y2Ba25gydF4y2Ba，gydF4y2Ba26gydF4y2Ba}。从形态学数据中进行的计算机分类，加上我们的多样性小组中专家的生态类型分配(补充数据)gydF4y2Ba5gydF4y2Ba)、露陆(gydF4y2BangydF4y2Ba= 268)，低地(gydF4y2BangydF4y2Ba= 99)和第三种沿海生态型(gydF4y2BangydF4y2Ba= 184)。沿海生态型与低地生态型具有广泛的同域性，但表现出高地叶片特征和低地植物构型(图2)。gydF4y2Ba2gydF4y2Ba扩展数据图gydF4y2Ba2gydF4y2Ba）.gydF4y2Ba

我们观察到强有力的证据表明，适应性进化促进了生态型分化。尽管北部高地植物的冬季死亡死亡率很少见(2.4%)，但在2018-2019年冬季，在4个最北部的花园中，近一半的沿海(42.1%)和低地(42.8%)基因型死亡。gydF4y2Ba2 bgydF4y2Ba）.在西北三个平原地区，冬季死亡尤其严重，可能是由于2019年1月下旬至3月初的一段严寒时期(扩展数据图)。gydF4y2Ba2gydF4y2Ba）.来自北部30%人群的基因型总数为218× (Fisher检验优势比= 218.17;gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−15gydF4y2Ba})在北部4个地点比最南部30%的基因型更有可能在2018-2019年冬季存活。gydF4y2Ba

在我们共同的花园中，纬度梯度也是生物量变化的主要轴。在7组相关气候变量中，生物量变化的最强预测因子始终与温度相关(扩展数据图)。gydF4y2Ba2gydF4y2Ba）.我们观察到极端30年最低温度的强烈信号，作为冬季易被杀死的低地和沿海生态型生物量的预测因子。gydF4y2Ba2摄氏度gydF4y2Ba）.对于这两种生态型，从历史极端最低温度较低的地点收集的基因型优于从北部花园中气候较温和的地点收集的基因型。然而，在耐冬杀旱地生态型中，没有观察到依赖于原产地气候的权衡。比2018-2019年更寒冷的冬季可能会在高地基因型中引入差异存活，并产生类似于在两种更南方的生态型中观察到的权衡。这些结果进一步支持了我们的观察，即对低温的敏感性既是自然选择的代理人，也是北方范围扩张的限制因素。gydF4y2Ba

此外，每种基因型的生物量产量通常在与其采集地点最相似的气候条件下达到最大(图2)。gydF4y2Ba二维gydF4y2Ba）.因此，局部适应不仅表现在生存和抗逆性上，还表现在与每种基因型进化相似的气候条件下更高的生物量积累上。gydF4y2Ba

了解柳枝稷基因库的结构和多样性对预测分子育种的未来收益和了解气候适应的遗传基础至关重要gydF4y2Ba^{12gydF4y2Ba，gydF4y2Ba13gydF4y2Ba}。先前对柳枝稷的群体遗传研究认为，群体遗传结构与用于定义生态型的形态聚类之间存在很强的对应关系gydF4y2Ba^{20.gydF4y2Ba，gydF4y2Ba27gydF4y2Ba，gydF4y2Ba28gydF4y2Ba}。对我们的3380万个全基因组SNP数据库的分析显示，我们的多样性面板被强烈地细分为三个主要的遗传亚群，总的来说，地理上是不同的(我们称之为中西部，大西洋和海湾)gydF4y2BaFgydF4y2Ba_圣gydF4y2Ba= 0.27)。gydF4y2Ba3gydF4y2Ba）.存在-缺失和结构变异的聚类在很大程度上概括了基于snp的亚种群结构(扩展数据图)。gydF4y2Ba3.gydF4y2Ba)，提供了亚种群分化的一致证据，其中可能包括几个分子尺度上的大效应突变。gydF4y2Ba

种群遗传结构与形态生态型差异不一致，在遗传亚群内分化强烈。在大西洋(37%)和中西部(63%)的基因库中均存在旱地生态型性状的植物。同样，54%和46%的沿海生态型加入分别分配给大西洋和海湾亚种群(图2)。gydF4y2Ba3gydF4y2Ba）.所有具有低地形态的植物都聚集在海湾亚群中。然而，这些海湾低地植物在北部冬季存活和死亡的个体比例大致相等。gydF4y2Ba2摄氏度gydF4y2Ba）.因此，重要的育种遗传多样性存在于遗传亚群中，这种模式通过在几个柳枝稷育种群体中实现生物量和冬季存活率的遗传收益得到了验证gydF4y2Ba^{29gydF4y2Ba，gydF4y2Ba30.gydF4y2Ba}。gydF4y2Ba

尽管亚种群之间存在生态型趋同，但coalescent模拟将亚种群的分化时间定在更新世中期(约358,000代(0.7-1.4 Ma，假设世代时间为2-4年))(扩展数据图。gydF4y2Ba3.gydF4y2Ba）.因此，现存的柳枝稷基因库在多倍体柳枝稷近一半的进化史中一直处于分化状态。与亚种群之间的深度序列差异相反，我们观察到大西洋亚种群中陆地生态型和沿海生态型之间的分子遗传差异很小(gydF4y2BaFgydF4y2Ba_圣gydF4y2Ba= 0.03)，或在海湾亚群的低地和沿海生态型之间(gydF4y2BaFgydF4y2Ba_圣gydF4y2Ba= 0.03)(扩展数据图gydF4y2Ba3.gydF4y2Ba）.gydF4y2Ba

混合在海湾和大西洋亚群之间似乎是常见的;比较这两个亚群的沿海生态型特征的植物在分子上更为相似(gydF4y2BaFgydF4y2Ba_圣gydF4y2Ba= 0.19)比非沿岸的墨西哥湾和大西洋植物(gydF4y2BaFgydF4y2Ba_圣gydF4y2Ba= 0.24)。相比之下，中西部和大西洋亚群中具有高地形态的植物并不比这些亚群中的其他植物更相似(gydF4y2BaFgydF4y2Ba_圣gydF4y2Ba两者均= 0.30)。这种高地形态在两个高度分化的遗传亚群中的趋同可能是高地生态型独立遗传起源的结果，也可能是罕见的，但在进化上很重要gydF4y2Ba^31gydF4y2Ba掺合料的事件。我们在下面评估这些假设。gydF4y2Ba

为了检测多样性面板内气候适应和适应度的遗传基础，我们进行了多变量适应收缩gydF4y2Ba^32gydF4y2Ba遗传亚群内部和跨遗传亚群的全基因组关联图谱(GWAS)结果。多变量自适应收缩在单变量测试之间共享GWAS峰值效应大小和方向，以提高检测显著共享结果的能力。确定了适应度GWAS(映射了三个最大的常见花园(MI, MO和TX)的冬季存活率和生物量)的多变量自适应收缩结果gydF4y2Ba_2gydF4y2Ba)和气候GWAS(检测SNP变异与起源气候(7个代表性气候变量)之间的关联)。为了在亚种群(具有不同分离snp)之间进行直接比较，我们总结了12239个显著的多变量自适应收缩(log)的连锁不平衡块“峰”gydF4y2Ba_10gydF4y2Ba-变换后的贝叶斯因子> 2)gydF4y2Ba^33gydF4y2Ba分成10,090个20 kb的区域(20 kb表示连杆不平衡衰减变平的拐点)(扩展数据图。gydF4y2Ba3.gydF4y2Ba)表示气候(gydF4y2BangydF4y2Ba_{地区gydF4y2Ba}= 9,856)及体能(gydF4y2BangydF4y2Ba_{地区gydF4y2Ba}= 332) gwas (gydF4y2Ba补充数据6gydF4y2Ba）.候选基因的加权列表-包括假定的SNP效应，存在-缺失或结构变异，基因共表达和GWAS峰的物理接近性-可以在补充数据中找到gydF4y2Ba7gydF4y2Ba。gydF4y2Ba

GWAS峰解释了基因库之间和基因库内部的大部分可遗传表型和气候变异(snp遗传力)。gydF4y2Ba3 bgydF4y2Ba）.适合度的单核苷酸遗传率(gydF4y2BahgydF4y2Ba^2gydF4y2Ba= 51.5±15.4%(平均值±s.e.m.))和与气候相关的峰值(gydF4y2BahgydF4y2Ba^2gydF4y2Ba= 70.5±14.0%)共同解释了比多基因背景(适应度= 19.5±9.1%，气候= 18.2±9.5%)多三倍的变异(扩展数据表)gydF4y2Ba1gydF4y2Ba）.这些气候和生物量关联的高遗传力表明，所有变异中一小部分的亲缘关系高于预测的整体亲缘关系，为育种者提供了遗传多样性，以在当地环境中进行柳枝稷改良。gydF4y2Ba

与适应性和原始气候相关的基因座可能与当地适应有关gydF4y2Ba^34gydF4y2Ba，是培育本地适应品种的有力目标。总的来说，我们观察到与气候和适应性相关的20 kb区域的重叠几乎是偶然预期的2倍(Fisher检验优势比= 1.92，gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−6gydF4y2Ba}）.这种重叠在两个北部亚群(中西部，赔率= 11.5×gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−15gydF4y2Ba};大西洋，赔率= 17.8倍和gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−15gydF4y2Ba})(图。gydF4y2Ba3 cgydF4y2Ba)，在那里，我们预计会看到选择对寒冷冬季生存的最大影响。gydF4y2Ba

许多气候和适应度重叠的区域仅在单个遗传亚群中多态，这突出了柳枝稷气候适应的几个可能独立的遗传途径。然而，9.5%(940)的20 kb气候区间在多个遗传亚群中存在多态性。鉴于海湾和大西洋亚种群之间混合的大量证据(图2)。gydF4y2Ba3gydF4y2Ba)，我们预计当代基因流动将是共享多态性的主要贡献者。与此假设相反，所有多亚群GWAS区间的大多数(511个区域)在两个遗传差异最大的基因库(大西洋和中西部)之间共享。考虑到这些亚种群之间的深度分化时间，罕见或古老的基因流动gydF4y2Ba^35gydF4y2Ba可能创造了这些共享的适应性多态区域。gydF4y2Ba

为了明确说明基因渗入如何影响气候snp关联的分布，我们使用隐马尔可夫模型研究了基因组中物理上连续的混合区域gydF4y2Ba^36gydF4y2Ba。亚群间的基因渗入占我们重测序基因组内容的2.98%(图2)。gydF4y2Ba4gydF4y2Ba)，但在亚群中包含共享GWAS间隔的可能性比偶然预期的高1.5倍(Fisher检验比值比= 1.55，gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−8gydF4y2Ba})，这表明适应性基因渗入至少是亚种群间共享的部分遗传变异的基础。gydF4y2Ba

特别令人感兴趣的是从中西部到大西洋亚群的一系列渗透，可追溯到距今约8700代(17 - 34000年前)，这与末次冰期极大期(约22 ka)后的北部范围扩张相吻合。具有较高中西部基因渗入水平的大西洋基因型表现出更多的高原性状(图2)。gydF4y2Ba4 bgydF4y2Ba)，并且在大西洋亚种群的其他亚热带和温带范围的北缘被过度代表(图2)。gydF4y2Ba4摄氏度gydF4y2Ba）.与其他系统中基因组渗入的适应性作用一致gydF4y2Ba^{31gydF4y2Ba，gydF4y2Ba37gydF4y2Ba}这些发现表明，假定的北方适应性等位基因从中西部进入大西洋亚群，可能促进了美国东北沿海地区较冷栖息地的柳枝稷在冰川后的殖民。为了验证这一假设，我们进行了冗余分析，将渗入块的存在与气候、地理和表型因素联系起来。总体而言，大西洋亚群中的中西部基因渗入与气候(解释方差百分比= 46.5%)的相关性是地理(11.5%)的四倍多。虽然从中西部到大西洋亚群的532和651个渗透分别与原始气候或生物量有关，但254个渗透在两个分析中都是异常值，这代表了每组之间独立性预期的近7倍富集(优势比= 6.99)。gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−15gydF4y2Ba}）.这些结果强化了中西部基因渗入塑造了北大西洋基因型的气候生态位和表型分布的假设，并支持了越来越多的证据表明适应性基因渗入可以促进范围扩大和生态型进化gydF4y2Ba^{38gydF4y2Ba，gydF4y2Ba39gydF4y2Ba}。gydF4y2Ba

多倍体在开花植物谱系中是常见的，可以增加可供选择的遗传多样性gydF4y2Ba^{40gydF4y2Ba，gydF4y2Ba41gydF4y2Ba}这可能导致适应性进化或分选，从而改变生态位特征gydF4y2Ba^42gydF4y2Ba。这一过程可以解释多倍体在极地纬度和高海拔地区普遍存在的原因，这些地区在冰川循环期间曾被冰盖覆盖gydF4y2Ba^43gydF4y2Ba。gydF4y2Ba

在多倍体形成过程中复制的基因可以亚功能化(将祖先基因的功能划分到旁系基因中)、新功能化(为旁系基因进化出新的基因功能)或简单地丢失gydF4y2Ba^44gydF4y2Ba。在多倍体物种形成之后，一个亚基因组通常保留更多的基因，并且平均表现出比另一个亚基因组更高的表达水平，这种现象被称为亚基因组优势gydF4y2Ba^45gydF4y2Ba。和其他多倍体一样gydF4y2Ba^{46gydF4y2Ba，gydF4y2Ba47gydF4y2Ba，gydF4y2Ba48gydF4y2Ba}在柳枝稷中，亚基因组显性和亚功能化表现明显。与N亚基因组相比，K亚基因组具有更高的基因密度(77.4个基因/ Mb比68.0个基因/ Mb)gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−15gydF4y2Ba})，更多的上调基因(5445对4477，二项gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−15gydF4y2Ba})和较低的突变积累率(K亚基因组中有5255个基因具有同义突变率(gydF4y2BaKgydF4y2Ba_{年代gydF4y2Ba})大于N亚基因组，而N亚基因组中有6,751个基因gydF4y2BaKgydF4y2Ba_{年代gydF4y2Ba}大于二项的K亚基因组gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−15gydF4y2Ba}）.综合起来，我们所有的11个亚基因组统计(扩展数据图)。gydF4y2Ba4gydF4y2Ba)指出了K亚基因组更强的进化约束和偏向，这表明适应进化的潜力可能在亚基因组之间被不同地划分。gydF4y2Ba

考虑到进化倾向于保留K亚基因组，我们预计会看到更强的气候适应信号gydF4y2Ba^44gydF4y2BaK亚基因组上的snp对生物量和存活率的影响。相反，在10个常见花园中，75.9%的生物量snp遗传力可归因于N亚基因组，而只有24.1%可归因于K亚基因组(扩展数据图)。gydF4y2Ba4gydF4y2Ba）.此外，在N亚基因组上发现了54.3%的中西部向大西洋亚种群的渗透，这是一个显著的富集(二项检验)gydF4y2BaPgydF4y2Ba< 1 × 10gydF4y2Ba^{−7gydF4y2Ba})，即使校正了7.5%的N亚基因组扩增(二项检验)gydF4y2BaPgydF4y2Ba= 0.0012)。归因于N亚基因组的大量渗入和可遗传的生物量变化可能与亚基因组对K亚基因组的进化倾向不一致。对于这一违反直觉的发现，一种可能的解释是，N亚基因组上放松的进化约束(减少净化选择)可能允许通过定向或多样化选择积累适应性遗传变异。因此，N亚基因组积累了可遗传的变异gydF4y2Ba^49gydF4y2Ba未来的育种制度可以以塑造天然柳枝稷种群和提高生物燃料产量为目标。gydF4y2Ba

随着气候和自然环境的变化，对驯化物种的遗传改良期望和野生种群的适应潜力进行限定变得越来越重要gydF4y2Ba^50gydF4y2Ba。事实上，植物基因组可以让我们窥见农作物和野生植物种群的过去和未来。对冰期-间冰期循环的适应为当前和未来的环境变化提供了一个有指导意义的模拟，我们在这里探讨了过去、现在和未来柳枝稷气候适应和产量提高的基因组机制。gydF4y2Ba

然而，植物基因组的复杂性也阻碍了遗传资源的开发，不利于快速有效的分子育种。我们的方法和成功测序柳枝稷复杂的基因组将促进生态和农业基因组学在几乎任何系统。例如，我们的研究结果表明，对北方气候的适应是由古代分化亚群之间的基因渗入促进的，这进一步支持了不同基因组之间的混合可以增强对新环境的适应的假设gydF4y2Ba^37gydF4y2Ba。这种适应性渗入和可遗传的亚基因组特异性遗传变异gydF4y2Ba^49gydF4y2Ba可能提供最小阻力的遗传路径，允许在环境变化期间在新栖息地定居。专性异交和多倍体——在选择基因组研究系统时经常有意识地避免的性状——是柳枝稷在自然界适应的主要驱动因素，也是遗传变异的来源，可用于选择，以通过不断变化的未来提高生物燃料产量。gydF4y2Ba

未采用统计学方法预先确定样本量。实验是完全随机的，研究人员在进行实验或测序时不知道基因型标识符。gydF4y2Ba

植物采集、繁殖、培养和表型分析gydF4y2Ba

为了形成多样性面板，从2010年到2018年收集了天然和普通园林来源的种子、根茎和无性系繁殖体。从种子长成的植物遵循标准的生长程序gydF4y2Ba^16gydF4y2Ba。简而言之，将10-15颗种子播种在9厘米的方形花盆中，其中含有ProMix BX盆栽土壤(Premier Tech Horticulture)和Turface MVP煅烧粘土(Turface Athletics)的混合物，并在4°C下春化7天。然后将花盆放置在光照的温室中，白昼长度为14小时，昼夜温度为30°C/22°C。幼苗在3叶期减薄至每罐1株，并让其生长至5分蘖期。根茎繁殖体和5分蘖苗被转移到5加仑的盆中，盆中含有细磨松皮地膜(Lone Star mulch)和缓释肥料(Osmocote 14-14-14, ScottsMiracleGro)。所有单株于2016 - 2018年通过无性系分裂在奥斯汀进行繁殖，每一株的繁殖目标为10个无性系。Cleary 3336F系统杀菌剂(Cleary Chemicals)在必要时施用于植物以控制真菌病原体。植物被放置在1加仑的花盆中进行最后的繁殖。gydF4y2Ba

2018年5月15日至7月10日在田间进行了种植，并遵循了先前公布的方法gydF4y2Ba^16gydF4y2Ba。简而言之，植物用卡车运送到每个地点，每个地点都覆盖着一层德威特除草布。植物被放置在杂草布上的洞中，这些洞被切割成蜂窝状，每个植物都有四个最近的邻居，彼此之间的距离为1.56米。为了防止边缘效应，在每个边缘位置都种植低地黑井品种。植物移植后用手浇水。所有植株的地上部分在2018-2019年冬季被保留，并在2019年春季分蘖出苗之前被移除。在2019年的季节结束时，植物被绑成一束，用镰刀割草机收割。gydF4y2Ba

我们生成了2019年生长季节的两个适应度指标:原木转化生物量(kg)和冬季存活率(kg)。gydF4y2Ba补充数据8gydF4y2Ba）.生物量数据是在2019年10月和11月收获期间从所有活个体中获得的。将估计质量<750 g的植物放在纸袋中，在60°C下干燥整株，直到没有额外的水分流失，然后称重总干生物量。在田间用±5克分辨率的悬挂标尺对估计质量为b> 750 g的植物进行湿生物量称重。为了确定这些植物的生物量，从每棵植物中提取约500克整株分蘖，称重，如上所述干燥并重新称重。然后将整个植物样品的湿生物量乘以子样品中的水分百分比，以近似总干生物量。在2018-2019年冬季，植物被认为经历了冬季死亡，因为到2019年6月1日，植物冠没有出现新的生长。2019年7月或9月，从实验中切除了死亡的植物冠，代之以布莱克威尔(Blackwell)品种的植物。gydF4y2Ba

基因组组装和抛光gydF4y2Ba

我们在美国能源部联合基因组研究所和哈德逊alpha生物技术研究所使用全基因组霰弹枪测序策略和标准测序方案对Alamo柳枝稷AP13基因型进行了测序。基因组由4,520,785 PacBio reads组装和抛光(121.66×原始序列覆盖范围来自59个P6C4 2.0和2.1化学细胞，电影时间为10小时，p-read产量为91.76 Gb)(扩展数据图)。gydF4y2Ba1gydF4y2Ba)，使用MECAT组装器gydF4y2Ba^52gydF4y2Ba和ARROW抛光机gydF4y2Ba^53gydF4y2Ba。最后用一个400 bp插入2 × 150 bp Illumina HiSeq片段文库(177.1×)进行基因组抛光和错误校正。读取>95%的简单序列重复，并在调整适配器和质量后读取<50 bp (gydF4y2Ba问gydF4y2Ba< 20,5 -bp窗口平均)被去除。最终的reads集包括1,259,053,614个reads，总共覆盖了168倍的高质量Illumina碱基。这产生了6600个支架(6600个contigs)的初始二倍体组装，其中contigs NgydF4y2Ba_50gydF4y2Ba1.1 Mb，大于100 kb的支架有3489个，总2C(二倍体)基因组大小为2,013.4 Mb。gydF4y2Ba

在近交个体中组装单倍体基因组，如AP13，通常会在杂合区域产生两个单倍体拷贝，这就需要计算步骤来将每条染色体表示为单拷贝单倍体，而不会不必要地重复重复副本。我们的初始组装大约是预期的单倍体(1C)基因组大小1.2 Gb的两倍。因此，为了检测假设的减数同源单倍型，我们鉴定并计数了恰好在组装中出现两次的共享24-mers，并相应地将其分类。共鉴定出3,152个较短和冗余的备选单倍型和2,387个重叠的contich末端，包括871.2 Mb的总序列。剩余的1,142.2 Mb序列通过比对2个可用图谱的遗传标记被排序并定向到18个染色体上(补充数据)gydF4y2Ba1gydF4y2Ba)出席MECAT大会;563个连接和57个断裂，10000个Ns代表未大小的间隙序列。总的来说，97.2%的组装序列包含在染色体中。端粒序列用(TTTAGGG)gydF4y2Ba_ngydF4y2Ba重复并正确定位。剩余的支架针对GenBank细菌蛋白和细胞器序列进行筛选，如果发现与这些序列匹配，则将其去除。为了解决小重叠区域上的contig末端，相邻的contig末端对齐到另一个使用BLATgydF4y2Ba^54gydF4y2Ba;共折叠了47个相邻的重复组对。gydF4y2Ba

我们进行了两轮纠错。首先，我们通过使用bwa mem将Illumina 2 × 150 bp文库与释放共识序列对齐来纠正纯合snp和插入和/或缺失(indels)gydF4y2Ba^55gydF4y2Ba并用GATK的统一基因型工具鉴定纯合snp和indelgydF4y2Ba^56gydF4y2Ba。共更正了690个纯合snp和80199个纯合索引。其次，我们使用ABI 3730XL毛细管测序仪组合计算完成了来自BAC克隆的11,343个组装contigs的测序gydF4y2Ba^57gydF4y2Ba和单指数Illumina克隆池，并将这组柳枝稷克隆与固定snp的基因组进行比对，以找到与其相邻的异相杂合snp。为了解决这些相位开关等位基因，将完整的原始PacBio reads与组装对齐。对于每次读取，确定每个杂合位点的期相，并纠正62,732个异期杂合位点。gydF4y2Ba

为了区分N和K亚基因组，我们使用了从头重复聚类方法，并通过与相关物种的系统发育距离验证了这一点。我们通过水母搜索“诊断”15-mersgydF4y2Ba^58gydF4y2Ba在吉普赛，Copia和Pao插入的LTR区域(由RepeatMasker识别)gydF4y2Ba^59gydF4y2Ba和LTRHarvestgydF4y2Ba^60gydF4y2Ba)，以区分每组同源染色体(一个同源染色体≤1命中，另一个同源染色体≥100命中)。具有共同15-mers的LTR序列被分组为超家族，并通过BLAST在每个超家族中进行比对。具有重大BLAST命中的超家庭成员(gydF4y2BaegydF4y2Ba< 0.01，≥90%长度)按Mafft分类并对齐gydF4y2Ba^61gydF4y2Ba。LTR家族之间的Jukes-Cantor距离由R ape包计算gydF4y2Ba^62gydF4y2Ba，并聚集成两组不同的亚基因组。LTRs和比对的聚类是相同的gydF4y2Bap . rudgeigydF4y2Ba(K.M.D.和E. Kellogg，未发表的数据)，它是K亚基因组的古代亲戚gydF4y2Ba^17gydF4y2Ba这给了我们很高的信心，我们已经有效地将所有染色体分配到正确的亚基因组中。最后，我们分配了染色体标识符，并通过与之合成来定位每条染色体伪分子gydF4y2BaSetaria italicagydF4y2Ba^63gydF4y2Ba。最终发布的单倍体5.0版本包含1125.2 Mb的序列，包含626个contigs，其中contigs N50为5.5 Mb，染色体中组装碱基的总数为97.2%。gydF4y2Ba

基因注释gydF4y2Ba

转录本组装由约20亿对2 × 150-bp的双链Illumina RNA-seq reads、约10亿对2 × 100-bp的双链Illumina RNA-seq reads和454个reads (gydF4y2Ba补充数据3gydF4y2Ba)使用PERTRAN(其详细信息先前已发布)gydF4y2Ba^64gydF4y2Ba）.简而言之，PERTRAN通过GSNAP进行基因组引导的转录组短读组装gydF4y2Ba^65gydF4y2Ba并对拼接进行对齐验证、重新对齐和校正后，构建拼接对齐图。总共约有450万个PacBio Iso-Seq循环共识序列gydF4y2Ba^66gydF4y2Ba被纠正并崩溃，产生大约677,000个假定的全长转录片段。随后，使用PASA构建了668,176个转录片段gydF4y2Ba^67gydF4y2Ba从RNA-seq reads，全长cDNA, Sanger表达的序列标签，以及校正和折叠的PacBio循环共识序列reads。基因座由EXONERATE测定gydF4y2Ba^68gydF4y2Ba柳枝稷转录物组装和蛋白质的比对gydF4y2Ba拟南芥gydF4y2Ba^69gydF4y2Ba、大豆gydF4y2Ba^70gydF4y2Ba北竹大米gydF4y2Ba^71gydF4y2Ba，gydF4y2BaSetaria冬青gydF4y2Ba^72gydF4y2Ba，gydF4y2Bap . halliigydF4y2Bavar。gydF4y2BahalliigydF4y2Ba^64gydF4y2Ba，gydF4y2Ba高粱二色的gydF4y2Ba^73gydF4y2Ba，gydF4y2BaBrachypodium distachyongydF4y2Ba^74gydF4y2Ba，葡萄和瑞士普罗特gydF4y2Ba^75gydF4y2Ba蛋白质组。使用RepeatMasker对重复软屏蔽柳枝稷基因组完成这些比对gydF4y2Ba^59gydF4y2Ba(重复库从RepeatModelergydF4y2Ba^76gydF4y2Ba和RepBasegydF4y2Ba^77gydF4y2Ba)，除非延伸到同一条链上的另一个位点，否则两端的长度可达2,000 bp。不完整的基因模型，即同源性低但转录组不完整的基因模型，或短的单外显子基因(<300 bp的编码DNA序列(CDS))，没有蛋白结构域或良好的表达。gydF4y2Ba

比较基因组学gydF4y2Ba

通过GENESPACE管道推断出两个柳枝稷亚基因组的同源同源物和旁同源物gydF4y2Ba^64gydF4y2Ba，使用默认参数和两个outgroups:gydF4y2Bap . halliigydF4y2Bavar。gydF4y2BahalliigydF4y2Ba^64gydF4y2Ba和gydF4y2Ba美国二色的gydF4y2Ba^73gydF4y2Ba。简而言之，GENESPACE将蛋白质相似度分数解析为合成块并运行正交器gydF4y2Ba^78gydF4y2Ba关于synsynconstrained爆炸结果。由此产生的块坐标和合成正交网络为进化推理提供了高置信度的锚点。gydF4y2Ba

为了计算CDS区域的祖先状态，我们首先使用来自Phytozome的基因组确定了具有共同祖先的序列gydF4y2Ba^79gydF4y2Ba。柳枝稷基因组的最终匹配数分别为38,960和33,772gydF4y2Bap . halliigydF4y2Ba,gydF4y2Ba美国二色的gydF4y2Ba,分别。对于任意给定的正交网络，我们在matt中建立了两个多序列比对gydF4y2Ba^61gydF4y2Ba一个不包括焦点柳枝稷序列(msa)gydF4y2Ba_0gydF4y2Ba)和一个强迫msagydF4y2Ba_0gydF4y2Ba通过——keeplength参数对准焦序列的坐标系。然后用最大似然算法在Phangorn中提取边缘特征状态gydF4y2Ba^80gydF4y2Ba。对于每次重建，只使用最接近柳枝枝条的内部节点作为祖先状态。总体而言，我们分析了40,943个柳枝稷基因模型(216,157外显子)，覆盖54.95 Mb (gydF4y2Ba补充数据9gydF4y2Ba）.gydF4y2Ba

亚基因组进化和年代测定gydF4y2Ba

为了推断柳枝稷亚基因组和四倍体柳枝稷的年龄，我们采用了一组保守的同源物，它们之间具有简单的2:1:1网络gydF4y2Bap . virgatumgydF4y2Ba，gydF4y2Bap . halliigydF4y2Ba和gydF4y2Ba美国italicagydF4y2Ba。这产生了45,045个柳枝稷蛋白，排列到24,549gydF4y2Bap . halliigydF4y2Ba通过跨物种分析，得到20,496对同源子和4,053个单子(K亚基因组2,396个，N亚基因组1,660个)。我们使用Dialign-TX对这些序列的翻译cd进行比对gydF4y2Ba^81gydF4y2Ba。将排列好的CDS序列连接到Gblocks中gydF4y2Ba^82gydF4y2Ba使用默认参数。Gblocks将18,044,244个CDS核苷酸筛选到16,321,302个位置，在50,334个区块中。然后在PhyML中使用结果对齐gydF4y2Ba^83gydF4y2Ba利用一般时间可逆模型构建最大似然树。这棵树被用作r8s的输入gydF4y2Ba^84gydF4y2Ba，计算时间树并校准gydF4y2Ba黍gydF4y2Ba- - - - - -gydF4y2BaSetariagydF4y2Ba树的节点为13.1 MagydF4y2Ba^63gydF4y2Ba。到目前为止，亚基因组分化和多倍体柳枝稷物种形成的时间，我们利用了爆发距离，它指的是LTR家族内的所有距离(而成对距离指的是相同插入的5 '和3 ' LTR之间的距离)。N-或k -亚基因组特异性反转录转座子的5 '和3 '距离用于确定这些元件的插入时间。此方法不能用于gydF4y2Bap . virgatumgydF4y2Ba特殊技能或gydF4y2Ba黍gydF4y2Ba因为这些元素最近的扩展主导了分布。相反，我们依靠比较最佳的跨物种比对来估计物种的LTR距离gydF4y2Bap . virgatumgydF4y2Ba- - - - - -gydF4y2Bap . halliigydF4y2Ba和gydF4y2Ba黍gydF4y2Ba- - - - - -gydF4y2BaSetariagydF4y2Ba节点。这样，我们就有了校准点，可以将LTR距离与物种之间更可靠的蛋白质编码基因差异进行比较。gydF4y2Ba

亚功能化和基因表达分析gydF4y2Ba

为了评估在蛋白质编码序列尺度上观察到的亚基因组进化偏差是否在表型中表现出来，我们探索了来自AP13叶片组织的生物复制同源物之间的基因表达偏差(gydF4y2BangydF4y2Ba≥5)在两个地点(TXgydF4y2Ba_2gydF4y2Ba和MI)。对Illumina配对末端RNA-seq 150-bp reads进行质量修剪(gydF4y2Ba问gydF4y2Ba≥25)，修剪后长度小于50 bp的被丢弃。对高质量序列进行比对gydF4y2Bap . virgatumgydF4y2Bav5.1参考基因组使用GSNAPgydF4y2Ba^65gydF4y2Ba使用HTSeq v.0.11.2获得与被注释基因唯一对应的reads计数gydF4y2Ba^85gydF4y2Ba。在DESeq2中通过似然比检验进行差异表达检验gydF4y2Ba^86gydF4y2Ba。在按亚基因组拆分reads之前计算文库大小;这些大小被用作分析差异表达的大小因子。从似然比检验中，亚功能化被定义为显著的亚基因组与环境的相互作用。利用条件下亚基因组间的事后wald测试对比，对田间花园和注释库的亚基因组表达偏差进行了测试。显著偏倚定义为差异表达假发现率调整gydF4y2BaPgydF4y2Ba< 0.05。用WGCNA对AP13基因注释RNA-seq文库进行加权基因共表达聚类gydF4y2Ba^87gydF4y2Ba6的幂。原始计数可在gydF4y2Ba补充数据10gydF4y2Ba。gydF4y2Ba

倍性评估gydF4y2Ba

我们使用LSRFortessa SORP流式细胞仪(BD Biosciences)来测定重测序材料的倍性水平。每株取200-300 mg幼叶组织，用刮片浸泡在培养皿中，用1 ml Cystain PI Absolute P核提取缓冲液(Sysmex流式细胞术)与1 μl 2-巯基乙醇混合处理15 min。用CellTrics 30 μm滤器(Sysmex)过滤样品，分离游离细胞核，用2 ml Cystain PI Absolute P染色缓冲液(Sysmex)、12 μl碘化丙啶和6 μl RNase a在湿冰上处理20 min，流式细胞仪测定细胞核大小，每个样品至少分析10,000个细胞核。流式细胞仪的输出用FlowJo软件(BD Biosciences)进行分析，并根据每个细胞核的平均荧光单位将样品分为三类(补充图2)。gydF4y2Ba1gydF4y2Ba）.如果细胞群有4 - 8万单位的荧光，则认为样品的倍性水平为4×， 8 - 10万单位为6×， 10 - 14万单位为8×。用流式细胞术数据建立分仓参数gydF4y2Bap . virgatumgydF4y2Ba已知倍性的材料。gydF4y2Ba

我们还通过双等位基因snp的变异等位基因频率分布来评估样本的倍性(如“变异召唤”所述)。该方法假设四倍体和八倍体遵循不同的等位基因频率分布模式，四倍体为0.5/0.5(参考深度和变异深度)，八倍体混合为0.75/0.25和0.5/0.5。如果命中的比例为0.48≤gydF4y2BaxgydF4y2Ba≤0.52 <0.035，认为文库为八倍体;≥0.035，认为文库为四倍体;870个样本中有837个(96.2%)的流式细胞术数据与这些结果相匹配。gydF4y2Ba

变量调用gydF4y2Ba

共对789份四倍体多样性样本进行了重测序，中位深度为59×(范围20×-140×)。其中，732份在过滤了缺失数据、异常杂合度升高和采集位点差异后用于进一步分析。样品采用HudsonAlpha生物技术研究所和联合基因组研究所的Illumina HiSeq X10和Illumina NovaSeq 6000对端测序(2 × 150 bp)进行测序。考虑到不同的库大小，读取被修剪到≤50倍的覆盖率，然后使用bwa-mem映射到v5程序集gydF4y2Ba^55gydF4y2Ba。gydF4y2Ba

通过将Illumina reads与AP13参考文献与BWA-mem对齐来调用snp。生成的.bam文件使用Picard (gydF4y2Bahttp://broadinstitute.github.io/picardgydF4y2Ba)，并使用GATK 3.0重新排列索引gydF4y2Ba^56gydF4y2Ba。多样本SNP调用使用SAMtools堆叠完成gydF4y2Ba^88gydF4y2Ba和Varscan V2.4.0gydF4y2Ba^89gydF4y2Ba最小覆盖率为8个，最小交替等位基因数为4个。通过二项检验来确定基因型。在24mer重复序列的25bp内的snp从进一步的分析中删除。仅保留了数据缺失率≤20%和较小等位基因频率>0.005的snp，导致75%的基因组中有33,905,042个snp，覆盖深度在8×至500×之间。使用SHAPEIT3进行分相gydF4y2Ba^90gydF4y2Ba。gydF4y2BaFgydF4y2Ba_圣gydF4y2Ba计算通过vcftools完成gydF4y2Ba^91gydF4y2Ba。我们通过生成732个文库和18条染色体每Mb的平均覆盖率来测试亚基因组读图偏差。然后，我们在lme4中对这些数据拟合了一个混合效应线性模型gydF4y2Ba^92gydF4y2Ba其中染色体数目(1-9)是随机效应，以检验亚基因组的主要效应。有和没有主效应项的模型通过似然比检验进行比较。gydF4y2Ba

使用HipMer构建了732个短读库的单独从头组装gydF4y2Ba^93gydF4y2Ba与一个gydF4y2BakgydF4y2Ba-mer的大小为101，以最大限度地在contigs之间进行单倍型分裂。由于组装体的质量和相邻性各不相同，根据组装体的总大小，将基因存在缺失和结构变异检测的样本集缩小到251个样本(泛基因组集)gydF4y2Ba_50gydF4y2Ba每个文库的长度和总基因比对。gydF4y2Ba

为了评估整个泛基因组中基因的存在-缺失变异，我们将所有AP13蛋白和一组独特的6161个来自美国的蛋白质进行了比对gydF4y2BaOropetium thomaeumgydF4y2Ba（gydF4y2BangydF4y2Ba_{蛋白质gydF4y2Ba}= 1476)gydF4y2Ba^94gydF4y2Ba，gydF4y2Ba美国italicagydF4y2Ba（gydF4y2BangydF4y2Ba= 1085)gydF4y2Ba^63gydF4y2Ba，gydF4y2BaSetaria冬青gydF4y2Ba（gydF4y2BangydF4y2Ba= 891)gydF4y2Ba^72gydF4y2Ba，gydF4y2Bap . halliigydF4y2Bavar。gydF4y2Ba菲利普gydF4y2Ba（gydF4y2BangydF4y2Ba= 1048)gydF4y2Ba^64gydF4y2Ba，gydF4y2Ba美国二色的gydF4y2Ba（gydF4y2BangydF4y2Ba= 878)gydF4y2Ba^95gydF4y2Ba和gydF4y2Bap . halliigydF4y2Bavar。gydF4y2BahalliigydF4y2Ba（gydF4y2BangydF4y2Ba= 772)gydF4y2Ba^64gydF4y2Ba。这些独特的基因是从通过正交器推断的单拷贝同源网络中提取出来的gydF4y2Ba^78gydF4y2Ba和选择由于缺乏同源的柳枝稷。使用BLAT将所有蛋白质(≥100个氨基酸)与所有从头组装对齐gydF4y2Ba^54gydF4y2Ba。如果AP13蛋白的基因序列大于或等于80%的同源性和75%的覆盖率，则认为存在，而其他草蛋白的基因序列大于70%的同源性和75%的覆盖率(允许物种之间更大的差异)。可变(泛基因组外壳)基因(被认为存在于40-60%的人群中);gydF4y2BangydF4y2Ba= 5,432)从存在-缺失变异矩阵中提取，并用于可视化来自大西洋，海湾和中西部亚群的非混合个体之间的差异。在每个亚群中检测显著过度或不足代表的基因是用agydF4y2BaχgydF4y2Ba^2gydF4y2Babenjamin - hochberg多重测试修正(gydF4y2BaPgydF4y2Ba≤0.05)。gydF4y2Ba

为了检测整个泛基因组的结构变异，使用nglr将每个文库中的contigs(≥2kb)与AP13参考基因组进行比对gydF4y2Ba^96gydF4y2Ba使用PacBio读取的默认设置。使用samtools对生成的.bam文件进行排序gydF4y2Ba^88gydF4y2Ba用来用鼻子召唤结构变体gydF4y2Ba^96gydF4y2Ba。使用SURVIVOR跨样本合并单个结构变体调用gydF4y2Ba^97gydF4y2Ba，最大允许距离为1 kb。使用bcftools过滤得到的.vcf文件gydF4y2Ba^88gydF4y2Ba最小等位基因频率为0.1，只考虑长度在100 - 1500 bp之间的插入和缺失。gydF4y2Ba

人口基因组学gydF4y2Ba

评价732份四倍体文库的遗传群体结构(补充资料)gydF4y2Ba4gydF4y2Ba)，我们提取了所有具有祖先状态调用的四重退化位点(相对中性)(补充数据)gydF4y2Ba9gydF4y2Ba)从祖先的状态排列。该位点列表代表我们的最高置信度中性位点，然后使用|的阈值对链接不平衡进行修剪gydF4y2BargydF4y2Ba|≤0.6，得到59,789个位点用于R包SNPRelate的下游分析gydF4y2Ba^98gydF4y2Ba。gydF4y2Ba

群体的连锁不平衡程度由snp决定gydF4y2Ba^99gydF4y2Ba在叮铃声gydF4y2Ba^{One hundred.gydF4y2Ba}。连杆不平衡(gydF4y2BargydF4y2Ba^2gydF4y2Ba)使用plink(——ld-window 500——ld-window-kb 2000)计算。的gydF4y2BargydF4y2Ba^2gydF4y2Ba每500bp取平均值。在R中使用nls函数对该数据拟合了一个非线性模型，并确定了连杆不平衡(gydF4y2BargydF4y2Ba^2gydF4y2Ba非线性曲线稳定。gydF4y2Ba

群体遗传结构分级评价。考虑到研究范围内存在高度分化的生态型，我们首先使用主成分判别分析(DAPC)分析了最广泛的遗传群体结构。gydF4y2Ba^101gydF4y2Ba在adegenet v.2.0.1中gydF4y2Ba^102gydF4y2Ba。这种方法不依赖于构成许多种群聚类方法基础的常见假设(例如，Hardy-Weinberg均衡和连锁不均衡)，因此提供了一种有价值的工具来研究广泛的结构划分。DAPC展示了一套强大的基因库，并将中西部基因型与所有其他基因型区分开来。然后，我们使用structure v.2.3.4中实现的贝叶斯聚类算法评估了剩余非中西部个体的遗传群体结构和潜在混合gydF4y2Ba^103gydF4y2Ba通过相关等位基因频率的混合模型。该分析包括20,000个燃烧步骤和30,000个1-6个基因型组的重复，每个重复运行10次。通过SNPRelate的特征向量分解，对所有亚种群的祖先系数进行事后分配。gydF4y2Ba

利用多序列马尔可夫聚结法(MSMCv.2.0)推断柳枝稷样本的人口统计学历史gydF4y2Ba^104gydF4y2Ba)，这是一种种群遗传学方法，用于从序列数据中推断人口历史和种群结构。该方法模拟了重组下聚结线虫的近似模型，并对种群大小和分化时间进行了检验。每个亚群使用4个单倍型运行MSMC，跳过模糊位点，估计rhoOverMu为0.25，时间段模式为10 × 2 + 20 × 5 + 10 × 2。我们对每个亚群的5组4个单倍型进行了100次迭代，在没有固定重组参数的情况下，估计rhoOverMu的平均值为0.25，并取其平均值。为了估计世代间的尺度分化时间，我们假设突变率为6.5 × 10gydF4y2Ba^{−8gydF4y2Ba}。为了估计初始分化时间，我们比较了相邻的相对交叉聚结率(RCCR)值(从过去到现在)(gydF4y2Ba补充数据11gydF4y2Ba）.如果在单个时间段内或在相邻时间段内或在两个交错时间段内出现下降(>0.01;观测范围为0.01 ~ 0.28)，以下相邻值几乎为零(≤0.009;观测范围:−0.1-0.009)，我们认为这是种群分离的起点。然而，如果在五个时间段内再次出现下降，我们认为后者是人口分离的开始。我们用16组不同的个体对每个亚群进行了重复分析。gydF4y2Ba

通过对距离矩阵的特征向量分解，对snp、结构变异和有无变异的种群结构进行了可视化分析。首先，计算3种变异呼叫类型的0/1/2(参考纯合子、杂合子、替代纯合子)文库×标记矩阵之间的欧氏距离矩阵。然后对欧几里得矩阵进行缩放和集中，通过在R包MDMR中实现的Gower中心相似矩阵去除库间覆盖率方差gydF4y2Ba^105gydF4y2Ba。gydF4y2Ba

生态型分类gydF4y2Ba

对成熟柳枝稷开花期或临近开花期的植株进行了16项性状调查(叶片长度、宽度、长宽比、面积、层厚和层/中脉厚度比;整株:分蘖数、分蘖高度、分蘖高度×数乘积、分蘖高度/计数比、穗高、穗高/计数比、叶冠高度、分蘖/叶高比;物候学:变绿日期和穗长出日期)，以确定2019年夏季在德克萨斯大学j.j. Pickle研究校园(PKLE;或TXgydF4y2Ba_2gydF4y2Ba(Austin, Texas, USA)和密歇根州立大学Kellogg生物站(KBSM;或MI (Hickory Corners, Michigan, USA))普通花园(见补充资料)gydF4y2Ba5gydF4y2Ba获取这些变量的详细描述)。物候测量，包括变绿(当第一个绿色营养结构从根茎冠出现时)和穗出苗(当第一个生殖结构从分蘖出现时)，每天进行分析。通过测量叶片的长度和宽度(单位:mm)、中脉和叶片厚度(单位:μm) (Mitutoyo 547-500S卡尺)和叶面积(单位:mm)，对每个植物的代表性叶片进行了详细的叶片形态评估gydF4y2Ba^2gydF4y2Ba) (Licor 3100C叶面积计)。除了这些数量性状外，我们还为叶片和整株外观生成了定性的高地-低地指数，这些指数是在2019年夏末在德克萨斯州奥斯汀收集的gydF4y2Ba_2gydF4y2Ba网站)。每种植物的特征以1-5的等级进行评估，从最低地到最高地。利用已建立的品种Alamo和Dacotah分别进行低地和高地性状的基线测定。植物性状评估包括:分蘖外观，从最厚和最像低地到最薄和最像高地;叶外观，从最宽，最长和最低地状到最短，最薄和最高地状;冠层颜色从最蓝和最典型的低地到深绿色和最典型的高地。这种视觉方法类似于柳枝稷育种者经常使用的基本选择标准。gydF4y2Ba

为了评估这些数据中的表型结构，我们使用了DAPCgydF4y2Ba^101gydF4y2Ba。首先利用主成分分析(PCA)对表型数据进行转换，确定前10个主成分gydF4y2BakgydF4y2Ba-means算法将个体分为3个可能的组，目的是最大化组之间的差异。接下来，对10个保留的主成分实施DAPC，使用两个合成变量(即判别函数)提供对生态型集群的有效描述，这两个合成变量是具有最大组间方差和最小组内方差的原始表型变量的线性组合。gydF4y2Ba

我们对调查的651个四倍体基因型中的16个性状在MI和TX进行了分类gydF4y2Ba_2gydF4y2Ba通过1层5个单元的低容量神经网络对3个生态型中的1个生态型(34个总特征，32个定量特征和2个定性有序特征)进行分析(补充数据)gydF4y2Ba5gydF4y2Ba）.该神经网络采用插入式实现gydF4y2Ba^106gydF4y2Ba并对7个已知生态型的栽培品种进行了培训(低地:Kanlow和Alamo;沿海:涨潮和斯图尔特;upland: Summer, Dacotah和Sunburst)和78个额外的基因型在同一个基于snp的遗传簇中(扩展数据图)。gydF4y2Ba3.gydF4y2Ba)，与样本品种在表型PCA空间中聚类最紧密。这些高亲和力的范例基因型在补充数据中打印gydF4y2Ba5gydF4y2Ba。对其余582个表型为生态型分类性状的基因型用插入符号进行了生态型预测gydF4y2Ba^106gydF4y2Ba。通过使用在代表北部和南部柳枝稷范围的花园收集的性状，我们希望避免当地气候对植物表型和随后的生态型分类的偏见。此外，与DAPC和专家资格相比，神经网络分类方法有一个显著的优势:因为神经网络锚定在已知和已发表的基因型上，包括这些常见品种的实验将能够更有效地概括我们的任务。gydF4y2Ba

混合和渗滤块的计算和年代测定gydF4y2Ba

我们通过两步管道建立了一个外加剂信息snp数据库。首先，祖先系数在“种群基因组学”中从具有相关祖先状态呼叫的四倍退化位点计算。每个亚种群中缺失数据最少且全基因组混合比例≤0.001的30个样本用于确定亚种群特异性等位基因频率。这些文库用于寻找至少有一个成对的snpgydF4y2BaFgydF4y2Ba_圣gydF4y2Ba值>0.4，由snpRelate函数snpGdsFst中的' W&C84 '方法计算得出。其次，在每个亚群中对这些全局遗传信息位点进行分析，找出等位基因频率较小的位点(0.05)和缺失位点(0.05)。这些位点首先在亚群体内被进一步修剪到|位点gydF4y2BargydF4y2Ba| < 0.9(10个snp或1000 -bp窗口)，则为|gydF4y2BargydF4y2BasnpRelate的| < 0.95 (1000 snp或10000 bp窗口)。这一过程导致每个亚群的SNP和文库计数如下:大西洋，579,468个SNP和284个文库;GULF, 641,975个snp和215个文库;中西部地区，481563个snp和196个图书馆。gydF4y2Ba

为了测试每对亚种群之间混合块的物理位置，我们使用了Ancestry_HMMgydF4y2Ba^{36gydF4y2Ba，gydF4y2Ba107gydF4y2Ba}。这种方法利用假定亲本群体中的等位基因频率来确定测试群体中可能出现基因渗入的区域。对于三个亚种群中的每一个，我们试图确定混合块渗入的时间、程度和当前位置。在每种情况下，我们允许其他两个亚群中的每一个产生两次脉冲。Ancestry_HMM可以优化祖先脉冲发生前的代数和混合脉冲中涉及的个体比例。然而，使用>48万个站点和>150个库进行8参数优化在计算上是不可行的。因此，我们使用40个随机抽样的文库进行参数优化，这些文库的混合系数在混合比例分布的0.2-0.8分位数范围内，并且snp仅在N亚基因组的4号染色体上。我们选择这条染色体作为其他染色体的代表是因为缺乏明显的大频率渗入。所得到的祖先脉冲参数优化是建立在1万代之前的初始未混合种群上的，然后为其他两个亚种群中的每一个建立两个后续的混合脉冲;优化后的脉冲如下(来源-参考):gydF4y2BangydF4y2Ba_{一代又一代gydF4y2Ba}= 8658和gydF4y2BaPgydF4y2Ba_{混在gydF4y2Ba}= 0.001%;67和0.7%)，海湾-大西洋(85和1.1%);17和0.25%)，大西洋海湾(79和1.9%;11和0.38%)，中西部-海湾地区(79和0.86%;大西洋-中西部地区分别为66和0.27%;14和0.036%)，墨西哥湾-中西部(71和0.15%;14和0.033%)。这些脉冲提供给包含所有个体和染色体的完整模型，误差概率为0.001，出现之前的最大代数为10,000，有效种群大小为100,000。后验祖先概率被解码成单倍型块，块被分类成相似位置块的簇。gydF4y2Ba

环境基因组学gydF4y2Ba

从“自然地球”网站(gydF4y2Bahttps://www.naturalearthdata.com/gydF4y2Ba）.各种绘图例程依赖于sfgydF4y2Ba^108gydF4y2Ba和光栅gydF4y2Ba^109gydF4y2Ba包在R环境中用于统计计算gydF4y2Ba^110gydF4y2Ba。气候数据下载自WorldClimgydF4y2Ba^22gydF4y2Ba(19个生物气候变量，0.5-arcmin分辨率1960-2000)gydF4y2Ba^21gydF4y2Ba。通过动态聚类探讨了气候变量在不同收集点之间的分布gydF4y2Ba^111gydF4y2Ba然后是围绕介质聚类进行分区gydF4y2Ba^112gydF4y2Ba与gydF4y2BakgydF4y2Ba= 7。最具代表性的气候变量被定义为与每个簇内的第一个变异特征向量最相关的气候变量。七个集群中有六个包含WorldClim变量。gydF4y2Ba

从NOAA门户网站下载了距离每个花园站点最近的气象站的天气数据，这些站点拥有2018年9月1日至2019年10月31日的完整日温度(最低-最高)和降水数据。每个花园使用的NOAA气象站标识符如下:IL (USC00110338)、MI (USW00014815)、MO (USW00003945)、NE (USC00255362)、OK (USW00053926)、SD (USC00391076)、TXgydF4y2Ba_1gydF4y2Ba(USC00414810), TXgydF4y2Ba_2gydF4y2Ba(USC00410433), TXgydF4y2Ba_3.gydF4y2Ba(USC00418862)和TXgydF4y2Ba_4gydF4y2Ba(USW00003901)。gydF4y2Ba

在原始数据和估算数据的基础上，对不同花园的气候表型关联进行了研究。纬度与生存的关系(图2)gydF4y2Ba2 bgydF4y2Ba)对原始数据进行了逻辑回归，通过在R中使用二项家族的glm进行了逻辑回归，在R中使用所有可用表型的近邻完成了输入(gydF4y2BakgydF4y2Ba= 5)，专门用于园林等级顺序的检验(图2)。gydF4y2Ba2 c, dgydF4y2Ba）.在混合线性模型中，通过lmer实现了气候相似性与生物量的关联gydF4y2Ba^92gydF4y2Ba，使用似然比检验将完整模型(固定=气候距离+截距，随机=基因型标识符)与没有气候距离固定效应的简化模型进行比较。gydF4y2Ba

采用物种分布模型(SDM)模拟了青藏高原所有生态型(高地、低地和沿海)的现代潜在分布范围gydF4y2Bap . virgatumgydF4y2Ba。用于建立SDMs的最终数据集包括277个(高地)、199个(沿海)和121个(低地)的发生记录。在我们最终的SDM模型中使用了6个环境预测因子(BIO1 =年平均温度，BIO2 =平均日差，BIO4 =温度季节性，BIO5 =最温暖月份的最高温度，BIO16 =最潮湿季度的降水量，BIO17 =最干燥季度的降水量)。然后用BIOMOD2 v.3.3生成sdmgydF4y2Ba^113gydF4y2Ba采用7种建模算法:广义线性模型、增强回归树、人工神经网络、灵活判别分析、随机森林、分类树分析和多元自适应回归样条。对于每个模型，发生率数据与在建模研究区域内随机生成的500个伪缺席数据相耦合，存在和伪缺席的权重相等gydF4y2Ba^114gydF4y2Ba。模型用80%的耦合出现和伪缺席数据进行训练，并用剩下的20%进行测试。每种建模算法对总共700个模型运行100次，并通过真技能统计(TSS)对这些模型进行评估。gydF4y2Ba^115gydF4y2Ba。TSS值从0.2到0.5被认为是差的，0.6到0.8被认为是有用的，>0.8被认为是好的到优秀的gydF4y2Ba^116gydF4y2Ba。基于TSS阈值(高地TSS阈值= 0.96，低地TSS阈值= 0.93，沿海TSS阈值= 0.965)，从700个模型中约50个最佳SDMs计算出独特的整体SDMs。最终的总体sdm被投影到当前的气候层上，以显示现代的潜在范围(补充资料)gydF4y2Ba12gydF4y2Ba）.gydF4y2Ba

我们通过在素食主义者中实施冗余分析，研究了大西洋亚群中中西部基因渗入的存在如何与气候、地理和亲属关系的独立和共同影响相关联gydF4y2Ba^{117gydF4y2Ba，gydF4y2Ba118gydF4y2Ba，gydF4y2Ba119gydF4y2Ba，gydF4y2Ba120gydF4y2Ba}。为了划分由气候、亲属关系和地理因素引起的渗透存在的可解释方差，我们运行了四个模型:一个完整的模型，具有渗透存在(潜在的渗透块被编码为0，大西洋遗传为1，中西部遗传为1)，由气候(即七个代表性气候变量)、亲属关系(从一组假定中性标记计算出的前两个主要成分)和地理(纬度和经度)解释，这三个因素中的每一个都以其他两个因素为条件。比较每个模型约束矩阵的惯性(即方差)值，以确定气候、亲属关系、地理及其联合效应的相对重要性。此外，为了找到与气候和生存校正生物量密切相关的渗入区域，我们从两个额外的模型中提取了冗余分析轴的负荷:(1)一个仅由气候预测，(2)一个仅由生存校正生物量预测。根据排列检验，两个模型均显著(gydF4y2BangydF4y2Ba= 999;gydF4y2BaPgydF4y2Ba均< 0.001)，所有轴近似为正态分布。每个轴尾部的snp负载更有可能表明与预测因子(即气候或生存校正生物量)相关的选择，因此我们确定了至少2.5 s.d(双尾)的所有标记gydF4y2BaPgydF4y2Ba= 0.012)作为自然选择下的渗入gydF4y2Ba^119gydF4y2Ba。gydF4y2Ba

GWASgydF4y2Ba

由于我们常用的花园数据集的规模很大，我们开发了一个管道- switchgrassGWAS R包(gydF4y2Bahttps://github.com/Alice-MacQueen/switchgrassGWASgydF4y2Ba)——允许在多样性面板上进行快速、低内存密集型的GWAS，并分析在不同位点测量的表型上SNP效应相似或不同的程度。这个包利用了bigsnprgydF4y2Ba^121gydF4y2Ba对编码为矩阵的大量SNP阵列进行快速(比TASSEL快300倍)统计分析。它还结合了当前人类遗传学文献中SNP质量控制、修剪和插入以及GWAS群体结构校正的金标准。为了测试在许多条件下(例如，多个地点、气候变量等)许多影响的重要性，我们使用了mashrgydF4y2Ba^32gydF4y2Ba这是一种灵活的数据驱动方法，可以共享任何数据集中效应大小和符号模式的信息，这些数据集可以在许多条件和snp的逐个条件基础上估计其影响。我们使用局部错误标记率确定哪些snp具有显著表型效应的证据，这类似于错误发现率，但更保守(因为它们也反映了效应标记估计的不确定性)。gydF4y2Ba^122gydF4y2Ba。我们用这些值用log找到snpgydF4y2Ba_10gydF4y2Ba-变换后的贝叶斯因子>2。这里，贝叶斯因子是一个SNP出现一个或多个显著表型效应的可能性与该SNP只产生零效应的可能性之比。根据之前的工作gydF4y2Ba^33gydF4y2Ba，贝叶斯系数为bb10gydF4y2Ba^2gydF4y2Ba被认为是支持SNP具有一个或多个显着表型效应的假设的决定性证据。gydF4y2Ba

为了计算气候和健康相关snp的区域遗传力，我们采用了先前描述的两步法gydF4y2Ba^123gydF4y2Ba。方差成分分析使用ASReml (VSN International)软件完成，使用van Raden方法计算的基因组关系矩阵gydF4y2Ba^124gydF4y2Ba。计算每个亚群内和整个多样性面板的基因组关系矩阵。基于单变量GWAS中使用的所有snp计算亲属关系矩阵(gydF4y2BaGgydF4y2Ba)，以及基于该亚群中与气候显著相关的snp的亲属矩阵(loggydF4y2Ba_10gydF4y2Ba-变换后的贝叶斯因子>2;gydF4y2Ba问gydF4y2Ba_{气候gydF4y2Ba})和基于与该亚种群生物量或冬季存活率显著相关的snp的亲缘关系矩阵(loggydF4y2Ba_10gydF4y2Ba-转换贝叶斯因子> 2，或>1.385海湾亚群;gydF4y2Ba问gydF4y2Ba_{健身gydF4y2Ba}）.这些亲缘关系矩阵用于区域遗传力作图gydF4y2Ba^123gydF4y2Ba就像之前的出版物一样gydF4y2Ba^125gydF4y2Ba，使用混合模型的形式:gydF4y2Ba

其中向量gydF4y2BaygydF4y2Ba表示生物量值，gydF4y2BaZgydF4y2Ba是随机效应的设计矩阵，gydF4y2BaugydF4y2Ba是整个基因组的加性遗传效应，gydF4y2BavgydF4y2Ba区域基因组的加性遗传效应和gydF4y2BaegydF4y2Ba是残差。矩阵gydF4y2BaGgydF4y2Ba为全基因组加性效应中使用所有snp的全基因组关系矩阵。矩阵gydF4y2Ba问gydF4y2Ba区域基因组关系是否如上所述得到:之一gydF4y2Ba问gydF4y2Ba_{气候gydF4y2Ba}或gydF4y2Ba问gydF4y2Ba_{健身gydF4y2Ba}。gydF4y2Ba我gydF4y2Ba等级是——gydF4y2BaygydF4y2Ba单位矩阵，其中gydF4y2BaygydF4y2Ba等于生物量值的个数。全基因组、区域基因组和残差方差gydF4y2Ba\({\σ}_{你}^ {2}\)gydF4y2Ba，gydF4y2Ba\({\σ}_ {v} ^ {2} \)gydF4y2Ba和gydF4y2Ba\({\σ}_ {e} ^ {2} \)gydF4y2Ba,分别。表型变异(gydF4y2Ba\({\σ}_ {{rm \ p {}}} ^ {2} \)gydF4y2Ba)是gydF4y2Ba\({\σ}_{你}^ {2}\)gydF4y2Ba+gydF4y2Ba\({\σ}_ {v} ^ {2} \)gydF4y2Ba+gydF4y2Ba\({\σ}_ {e} ^ {2} \)gydF4y2Ba。全基因组遗传力、区域遗传力和总遗传力gydF4y2Ba\ ({h} _{你}^ {2}\)gydF4y2Ba= (gydF4y2Ba\({\σ}_{你}^ {2}\)gydF4y2Ba/gydF4y2Ba\({\σ}_ {{rm \ p {}}} ^ {2} \)gydF4y2Ba)，gydF4y2Ba\ ({h} _ {v} ^ {2} \)gydF4y2Ba= (gydF4y2Ba\({\σ}_ {v} ^ {2} \)gydF4y2Ba/gydF4y2Ba\({\σ}_ {{rm \ p {}}} ^ {2} \)gydF4y2Ba),gydF4y2Ba\ ({h} _ {u + v} ^ {2} \)gydF4y2Ba= (gydF4y2Ba\({\σ}_{你}^ {2}\)gydF4y2Ba+gydF4y2Ba\({\σ}_ {v} ^ {2} \)gydF4y2Ba/gydF4y2Ba\({\σ}_ {{rm \ p {}}} ^ {2} \)gydF4y2Ba),分别。gydF4y2Ba

这些模型是在进行亚种群GWAS的三个地点运行的:哥伦比亚，密苏里州;密歇根州的山胡桃角;以及德克萨斯州的奥斯汀。这产生了80个模型:4组种群(全多样性面板和3个亚种群)，2种模型类型(一个模型与gydF4y2BaGgydF4y2Ba只有和一个gydF4y2BaGgydF4y2Ba+gydF4y2Ba问gydF4y2Ba模型)，10种表型(3个地点的生物量和7个环境变量)。gydF4y2Ba

方差成分分析也用于划分K亚基因组和n亚基因组之间的方差。只有具有祖先状态的snp调用(gydF4y2Ba补充数据9gydF4y2Ba)在该分析中使用，结果每个种群子集使用了460,429个snp。根据特定染色体上的所有snp计算亲属矩阵(gydF4y2Ba问gydF4y2Ba_{Chr01KgydF4y2Ba}来gydF4y2Ba问gydF4y2Ba_{Chr09KgydF4y2Ba},gydF4y2Ba问gydF4y2Ba_{Chr01NgydF4y2Ba}来gydF4y2Ba问gydF4y2Ba_{Chr09NgydF4y2Ba})，得到18个亲属矩阵。这些亲属关系矩阵被用于区域遗传力映射，使用混合模型的形式:gydF4y2Ba

其中向量gydF4y2BaygydF4y2Ba表示生物量值，gydF4y2BaZgydF4y2Ba是随机效应的设计矩阵，gydF4y2BavgydF4y2Ba_{1 kgydF4y2Ba}(gydF4y2BavgydF4y2Ba_{9 kgydF4y2Ba})或gydF4y2BavgydF4y2Ba_{1 ngydF4y2Ba}(gydF4y2BavgydF4y2Ba_{9 ngydF4y2Ba})(统称为gydF4y2BavgydF4y2Ba_我gydF4y2Ba)是染色体特异性基因组加性遗传效应和gydF4y2BaegydF4y2Ba是残差。矩阵gydF4y2Ba问gydF4y2Ba_我gydF4y2Ba是N和K亚基因组的9条染色体的染色体特异性基因组关系矩阵。染色体特异性和残差是gydF4y2Ba\({\σ}_ {{v} _{我}}^ {2}\)gydF4y2Ba和gydF4y2Ba\({\σ}_ {e} ^ {2} \)gydF4y2Ba,分别。染色体特异性遗传力是gydF4y2Ba\ ({h} _ {{v} _{我}}^ {2}\)gydF4y2Ba= (gydF4y2Ba\({\σ}_ {{v} _{我}}^ {2}\)gydF4y2Ba/gydF4y2Ba\({\σ}_ {{rm \ p {}}} ^ {2} \)gydF4y2Ba)，而亚基因组特异性遗传力是每个亚基因组内9条染色体上这些差异的总和。gydF4y2Ba

候选基因探索gydF4y2Ba

我们整合了多个数据结构，对渗入间隔内的候选基因进行排序，并提供有意义的筛选标准，以及与数量性状位点峰值的物理接近度。在GWAS峰的情况下，候选基因被定义为围绕mashr峰的20 kb间隔内的位点。基因组渗入的候选基因必须至少部分重叠于渗入间隔。由于GWAS和基因渗入的推断是在遗传亚群内进行的，所有统计数据均报道于gydF4y2Ba补充数据7gydF4y2Ba(候选基因列表)也具有亚种群特异性，但基因共表达分析除外(该分析仅对用于注释目的的AP13 rna测序文库进行)(gydF4y2Ba补充数据3gydF4y2Ba）).对于给定的区间，我们给出一组统计量。首先，将峰值位置的物理接近度计算为基因到间隔(渗入)或GWAS峰值位置的中点。其次，由于一个亚群体中GWAS峰值的因果位点在该亚群体中必然是可变的，我们提取了候选基因模型内和邻近的所有snp。这些变体用SNPeff进行了注释gydF4y2Ba^126gydF4y2Ba以及三大类变量的加权和(高、中、低;这些的描述可以在gydF4y2Bahttps://pcingola.github.io/SnpEff/se_inputoutput/#effect-prediction-detailsgydF4y2Ba)计算为SNPeff_score =高× 20 +中× 5 +低× 1。第三，对于每个基因，我们计算了结构变异和存在缺失变异的次要等位基因频率。第四，我们包括每个基因的WGCNA集群的身份载体。最后，如果候选基因是先前发表的开花时GWAS候选基因的同源物gydF4y2Ba^127gydF4y2Ba，包括重叠区间或基因的身份。gydF4y2Ba

报告总结gydF4y2Ba

有关研究设计的进一步资料，请参阅gydF4y2Ba自然研究报告摘要gydF4y2Ba链接到本文。gydF4y2Ba