文献分享33：泛基因组分析揭示了与黍稷驯化相关的基因组变异

文献

2023

Nature Genetics

Pangenome analysis reveals genomic variations associated with domestication traits in broomcorn millet

课题背景

（1）气候变化对全球粮食安全造成了严重威胁。尽管已经培育出高产的主粮作物，但孤儿作物为实现气候适应型农业和增加食品供应提供了机会。然而，尽管在低投入条件下呈现出丰富的营养多样性，孤儿作物仅由小农和边缘农民在本地种植。因此，研究这些作物可能有助于改善主粮作物的营养多样性和环境适应力。 

（2）黍稷（Panicum miliaceum L.）是一种主要在亚洲和欧洲半干旱地区种植和消费的孤儿作物。它在10000年前的中国北方被驯化，在稻米和小麦兴起之前，是当地的主粮作物。黍稷在大约3600至4000年前传播到欧洲。由于其无麸质、高蛋白含量、生长迅速和抗旱特性，黍稷成为主要谷物的替代品。然而，尽管需求和收获面积逐渐增加，但只有少数品种被推向市场。此外，黍稷的基因组多样性尚未得到广泛研究，其驯化的遗传基础仍有待探索。

亮点

本研究旨在分析群体水平黍稷的基因组，以确定其起源并探索与驯化相关的农艺性状的遗传基础。作者利用PacBio数据生成了32个代表性样本的基因组，并构建了图形泛基因组，以揭示黍稷群体中的基因组变异。

结论1 基因组组装及群体结构分析

Fig 1

为了探索黍稷的遗传多样性和群体结构，作者对516份材料（包括415份地方种、38份栽培种和63份野生种）进行了基因组测序（Fig 1a）。比对到参考基因组后，鉴定出1,890,542个高质量SNP和168,878个InDel。

为了确定群体结构，作者利用四重简并位点进行了PCA（Fig 1b-c）和STRUCTURE分析（Fig 1d）。结果表明，516份材料可以分为四簇：一个野生簇W1和三个栽培簇，即C1、C2和C3。最大的簇C3包含了来自中国西北的栽培种，这是黍稷种植的主要地区。C1簇由来自中国东北和东部的栽培种组成，而C2簇包括来自欧洲和中亚国家的栽培种（Fig 1a）。这些发现表明，黍稷的群体结构与地理位置高度相关。

进一步分析发现，黍稷的栽培种和野生种的核苷酸多样性低于水稻和大豆（Fig 1e），栽培种保留了其野生亲属中62.6%的核苷酸多样性。在系统发育树中，野生群体形成了一个与三个栽培群体不同的簇（Fig 1f）。C3簇表现出较高的复杂性，并与野生群体密切相关，暗示C3簇代表了从野生种中驯化出的基因库。此外，来自新疆和甘肃的一些材料在C3簇内构成了基部谱系或位于C2分支内，表明欧洲和中亚的标本可能起源于中国西北的材料。此外，作者还发现了栽培种与野生群体之间的基因流。这些结果表明，黍稷在中国北方被驯化，随后从中国西北向西传播。

结论2 黍稷的泛基因组分析

Fig 2a

作者选择了32份材料，包括24份栽培种和8份野生种，用来构建黍稷的图形泛基因组（Fig 2a）。作者使用PacBio HiFi数据通过hifiasm进行组装，得到的contigs的N50范围为5.16至27.25 Mb。QV、LAI和BUSCO评估表明组装的参考基因组具有高质量且高度连续。

Fig 2b-d

作者进一步发现每个基因组平均含有58.1%的重复序列，较大的基因组含有更多的转座子。在32个基因组中，共鉴定出27727个核心基因，8288个软核心基因，24494个可有可无的基因和5533个特色基因家族（Fig 2b-d）。

Fig 2e-h

结论3 黍稷的结构变异

Fig 3a-c

为了进一步探索黍稷的基因组多样性，作者采用基于组装的方法来识别基因组中的结构变异，共鉴定了207,033个SVs（Fig 3a），合并非冗余SVs后用于后续分析（Fig 3b-c）。

50515个PAVs中有59.4%仅出现在一个或两个样本中，这与水稻和大豆中低频率PAVs的情况一致，表明它们可能代表新出现的或有害突变。基于50515个PAVs，作者构建了图形泛基因组。

Fig 3d-l

作者将与TEs有90%以上重叠的PAVs分类为TE-PAVs，其余的PAVs被分类为非TE-PAVs，TE-PAVs占据了所有PAVs的大多数（68.3%）（Fig 3d）。

作者进一步根据PAVs相对于蛋白编码基因的位置进行了注释，发现32.9%与基因区域重叠（Fig 3e）。在所有非TE-PAVs中，51.7%位于基因区域，而只有13.2%的TE-PAVs与基因区域相关联。此外，DNA-TE PAVs比LTR-TE PAVs更靠近基因区域（Fig 3f）。为了了解PAVs如何影响基因功能，作者比较了PAV基因和非PAV基因的表达水平，发现PAV基因的表达水平显著低于非PAV基因。此外，PAV基因中有更多的沉默基因（Fig 3g），表明PAVs与叶片和根中的基因表达降低相关。此外，在基因编码区域和基因上游位置的TE-PAVs更有可能影响基因表达，而在基因内含子和基因下游位置的TE-PAVs较少影响基因表达（Fig 3h）。综上，这些发现表明PAVs通过改变编码和顺式调控区域来影响基因表达。

此外，作者鉴定了在叶片和根中表达水平显著改变的648个PAV基因。这些差异表达的PAV基因富集了与抗性相关的结构域，如NB-ARC和Rx N-末端结构域，这也是PAV基因富集的Pfam结构域（Fig 3i）。作者发现抗性基因位于重复富集区域，周围PAVs的频率高于基因组平均水平，表明PAVs与黍稷中抗性基因的演化有关。例如，在BC027中，作者发现两个抗性基因之间存在一个4.6 kb的插入，其中longmi055791编码ENHANCED DISEASE RESISTANCE 2的同源蛋白，而longmi055792编码一个NBS-LRR基因（Fig 3j）。这个插入与两个基因的表达下降相关（Fig 3k），其等位基因频率在C2（73.3%）和C3（66.0%）群体中高于W1（17.8%）（Fig 3l）。这些结果表明这个突变可能促使了黍稷对中国西北（C3）和欧洲（C2）的适应。

结论4 黍稷驯化过程的人工选择

Fig 4a

作者整合了XP-CLR、π比值和群体分化指数来检测人工选择的信号，共识别出524个基因组区域作为人工选择的靶点，涵盖30.2 Mb的序列和3910个蛋白质编码基因（Fig 4a）。这些区域与一些已知的与驯化和适应性状相关的基因重叠，如谷物产量（GL3.1、SG1和GS1）和开花时间（Ghd2、Ehd1和Hd5）（Fig 4a）。

Fig 4b-h

黍稷是一种含有两个亚基因组的异源多倍体物种。作者发现，在A亚基因组中，受选择区域更为丰富，包含更多的蛋白编码基因。作者还观察到，B亚基因组中的蛋白编码基因比A亚基因组中的基因具有更多的PAVs（Fig 4b)。

此外，作者分析了与PAVs相关的基因丢失和拟基因化事件，以了解PAVs如何影响异源多倍体基因组中的基因丢失。作者鉴定了由PAVs删除或拟基因化的1,321个基因，发现B亚基因组经历了比A亚基因组更多的基因丢失，表明PAVs促进了基因丢失。作者还鉴定了242个基因丢失事件，这些事件在野生群体中的频率低于栽培群体。这些结果表明，黍稷基因组正在进行的二倍体化可能已经影响了基因功能，从而促成了它的驯化。

为了更好地了解基因组变异在驯化过程中如何影响基因功能，作者在225个受选择区域中鉴定了1,099个PAVs，其中包括39.9%的TE-PAVs。其中，503个PAVs与基因区域重叠。作者还发现在驯化过程中，野生和栽培群体之间有5,663个PAVs的等位基因频率发生显著变化。整合上述分析中受PAVs影响的基因，作者鉴定出4,930个可能与黍稷驯化相关的基因。在longmi031198的上游区域中发现了一个6.4 kb的TE插入（Fig 4c），它是水稻花素基因Hd3a的同源基因。这个突变（Alt）与开花表型没有显著关联（Fig 4d）；然而，栽培群体显示出不含插入的单体型（Ref）的等位基因频率增加（Fig 4c, e）。作者还鉴定到longmi040672附近的三个紧邻的缺失变异（Fig 4f）。具有缺失变异单体型与穗和主茎之间的夹角较大相关（Fig 4g），并在驯化过程中受到负选择（Fig 4h）。这些结果表明，PAVs（尤其是TE-PAVs）在黍稷驯化过程中可能起到了重要的作用。

结论5 黍稷驯化过程的基因组变异

Fig 5a-b

此外，为了将黍稷的基因组变异与表型变异联系起来，作者对两年七点的516份材料进行了43个性状的测定（Fig 5a），并基于1,890,542个SNPs和19,492个PAVs进行了GWAS。SNP-GWAS识别出与31个性状相关的139个位点，包括许多农学重要性状，如种子尺寸和植株结构，以及与驯化相关的农艺性状，如落种性（SHT）和总状花序结构（PNT）。此外，PAV-GWAS揭示了与17个性状相关的70个PAVs。作者重点分析了种子脱落性、花序和种子颜色以及总状花序结构（Fig 5b）。

5.1 种子脱落性

Fig 5c-h

种子脱落性的减弱是谷物驯化的关键步骤。作者发现栽培群体的脱落性水平低于野生群体（Fig 5b），表明在驯化过程中对非脱落性性状进行了强烈的选择。

为了揭示与栽培群体中非脱落性性状相关的遗传变异，作者检查了15个已知与谷物脱落性相关的基因，发现longmi009317，是水稻和高粱中控制种子脱落性的OsSh1的同源基因，以及一个相关同源基因longmi003952，在黍稷中经历了基因丢失或拟基因化。类似地，longmi003952的一个3.2 kb的缺失导致了六个外显子的丧失。作者还在受选择区域中鉴定了SSH1/OsSNB的同源基因longmi012879, 其Hap 4与野生材料中的种子脱落性显著相关。这些观察结果表明，在黍稷的驯化过程中可能已经选择了与非脱落性性状相关的多个位置的基因组变异。

为了进一步识别控制黍稷种子脱落性的基因组变异，作者利用GWAS鉴定到58个SNP与脱落性显著相关。作者发现两个基因，longmi020192（编码果胶酯酶，PE）和longmi028230（编码果胶酯酶抑制剂，PEI），位于受选择区域（Fig 5c）。PE负责中层板中果胶的降解，而PEIs可以抑制果胶甲酯酶的脱酯化。作者在longmi020192和longmi028230的编码区域中发现了两个非同义SNP。此外，携带这些非同义突变的单倍型与野生群体中的种子脱落性高度相关（Fig 5d-e）。

5.2 花序和种子颜色

Fig 5i-o

花序和种子的颜色是与植物适应性、应激响应和营养代谢相关的性状。形态学调查表明，在黍稷的栽培种群中，绿色花序和深色种子占主导（Fig 5b），表明在栽培群体中对这些特征进行了选择。作者通过GWAS鉴定了与花序颜色（IFC）相关的1211个显著的SNP（Fig 5f）。在longmi008332的调控区域发现了两个SNP，编码谷胱甘肽S-转移酶（Fig 5f），与植物中花青素积累相关。这两个SNP形成了两个主要的单倍型，大多数具有Hap2的品种具有紫色花序，而大多数Hap1的品种表现为绿色花序。此外，PAV-GWAS鉴定到与紫色花序相关的longmi008332中的一个72bp插入（Fig 5g-h）。所有没有这个插入的品种都具有Hap1，而带有这个插入的品种中有76.6%表现为Hap2，表明这段72bp的插入是引起黍稷紫色花序的变异。

在SNP-GWAS中检测到与种子颜色（SDC）相关的多个基因座，分布在chr 5-6-9-11和14上。其中，在chr 9上发现了483个相关的SNP，主要集中在色氨酸脱羧酶（TDC）基因簇周围（longmi004409、longmi004412和longmi004413）。TDC是与5-羟色胺生物合成相关的基因，其上调导致种子或叶子呈深褐色。单倍型分析表明，与深色种子相关的Hap 3、4和5主要存在于栽培品种中。此外，一个SNP导致了longmi057520中的提前终止，并与深色种子高度相关。Longmi057520与水稻中合成硫酸雪松醇和丁香醇前体的GH2同源。GH2突变体的种子呈金黄色，而野生型的种子呈浅黄色。因此，longmi057520中的提前终止可能导致黍稷栽培群体的深色种子。此外，longmi057520位于一个选择清除区域。这些观察结果表明，SDC是一个由多个遗传因子控制的复杂性状，在黍稷驯化过程中经历了选择。

5.3 花序构造

花序构造是影响谷物产量的关键因素，也是作物驯化和改良的焦点。在黍稷中，野生品种具有敞开的花序，而栽培品种具有闭合的花序（Fig 5b, i）。作者对八个与花序形状相关的性状进行了SNP-GWAS分析（Fig 5b），在chr 14上鉴定出与四个与花序相关的性状相关的55个基因（Fig 5j）。其中，longmi038370编码一个SBP结构域蛋白，是玉米中控制叶片角度、穗支数和穗支角度的LG1的同源基因。单倍型分析表明，Hap 1和5与较低的花序密度（DOI）以及谷穗和主茎的较大分支强烈相关。此外，longmi038370位于一个受选择区域中，表明它在黍稷的驯化过程中受到了强烈的选择。为了验证longmi038370的功能，作者构造了三个Seita.3G022100.1（SiLG1）的CRISPR–Cas9突变体，该基因是黍稷在早熟禾中的同源基因。Silg1-1突变体显示失去了箨舌（Fig 5k, l），花序支与主茎之间的角度较小（Fig 5m, n），花序较为紧凑（Fig 5o）。因此，作者认为longmi038370控制了黍稷的花序形状。

此外，作者还鉴定了与谷穗主轴方向（MED）相关的chr 1和chr2上的基因座（Fig 5j），其中包含一个候选基因，即longmi054645，编码一个无顶端分生组织蛋白。单倍型分析表明，Hap 3与野生品种中的低MED强烈相关。另一个候选基因 (longmi011403)编码一个依赖钙的磷三酯酶蛋白，也在chr1上被鉴定出来。Longmi011403的Hap 2仅存在于野生品种中，它与较低的DOI和敞开的花序强烈相关，表明它是多功能的，控制了黍稷中多个花序性状。