《Molecular Ecology》：Patterns of Interploidy Admixture in Polyploid Complexes: Insights From Thymus Sect. Mastichina (Lamiaceae)
编辑推荐：
本研究通过整合系统基因组学（Hyb-Seq）与群体基因组学（SNPs）方法，深入揭示了伊比利亚半岛特有百里香属Mastichina组多倍体复合体的进化动态。研究证实该组由二倍体和四倍体两个姐妹群构成，其中二倍体群体存在Algarve、Cádiz、Do?ana和Hercynian四个具有网状进化特征的亚群，而四倍体群体遗传分化较低。研究首次在该复合体中检测到显著的跨倍性基因流（interploidy gene flow），其方向主要以四倍体向二倍体渗透为主。这些发现凸显了该组分类学上的隐蔽多样性（cryptic diversity），强调了结合形态学、细胞学、基因组学和生态学证据的整合分类学（integrative taxonomy）对于物种界定和濒危物种（如T. albicans）保护的重要性。
1 引言
多倍化（polyploidy）作为由全基因组复制（Whole-genome Duplication, WGD）介导的物种形成机制，在植物进化与生态中具有深远意义。尽管历史上多倍化曾被视为进化的“死胡同”，但大量近期起源的多倍体存在支持了其作为植物多样化驱动者的潜力。多倍化可通过同源多倍化（autopolyploidy）或异源多倍化（allopolyploidy）产生，后者通过杂交促进不同谱系间的基因流，从而促成新物种的出现。新建立的多倍体最初可能因基因组不稳定性或基因表达变化而经历适合度降低，但一旦克服这些问题，其快速适应的潜力可赋予选择优势。多倍化可通过增加遗传变异、表观遗传重塑以及遗传稳健性（genetic robustness）来促进短期内的成功，这些变化可能影响种间相互作用，赋予更广泛的生态耐受性，促进新栖息地的殖民，并增加入侵潜力。
多倍化可通过形成自然的生殖屏障（如三倍体阻滞，triploid block）促进植物的同域物种形成。然而，这种跨倍性杂交的屏障有时是可渗透的。因此，同域中新多倍体物种的建立可能依赖于新多倍体与其二倍体祖先之间的合子前屏障，以及从异交到自交或从有性生殖到无性生殖的转变，这可以缓解新形成多倍体的细胞型少数排斥（cytotype minority exclusion）。然而，WGD后的生殖屏障可能并不完整，特别是在中间细胞型（如三倍体）可行的情况下，这可能促进与亲本细胞型的回交（即三倍体桥梁途径，triploid bridge pathway）。此类回交事件如果形成的可育跨倍性杂种能够与其亲本回交，则可能实现双向基因流。此外，未减数配子（unreduced gametes）也可能介导从二倍体到四倍体的单向基因流。
跨倍性基因流（interploidy gene flow）被证明是许多多倍体复合体进化的特征，其中不同倍性水平的物种在不明显改变染色体数目的情况下交换遗传物质。基因渗入（introgression）和多倍化虽然在文献中通常作为独立的生物学现象进行讨论，但它们并非互斥，并且在自然界中经常共存。这些基因渗入事件允许受体物种不仅通过与近缘二倍体祖先的基因渗入，而且通过与更远缘的二倍体或多倍体供体谱系的基因渗入来积累新等位基因，从而贡献于遗传多样性并促进适应性等位基因的出现。尽管对其相互作用的认识日益增长，但关于多倍体复合体中跨倍性混合（interploidy admixture）如何塑造进化轨迹的理解仍存在显著空白。此外，由于多倍体复合体和基因渗入常常模糊系统发育关系并使表型区分变得模糊（连同其他因素如不完全谱系分选，Incomplete Lineage Sorting, ILS），传统的分类学物种区分可能具有挑战性。这对于隐蔽和受威胁物种尤其相关，因为准确的识别保证了保护策略的正确评估和实施。
多倍化、基因渗入和杂交是包括唇形科（Lamiaceae）在内的几个植物科中的常见现象，尤其是在百里香属（Thymus L.）中。该属以其高度的物种多样性、显著的经济价值、高水平的地方特有性以及生态重要性而著称，因为百里香群落是地中海灌丛最具代表性的群落之一。在伊比利亚半岛，百里香属Mastichina组为研究由多倍化和基因渗入过程驱动的进化模式提供了一个极好的模型。该组包括两个密切相关的物种（三个分类单元），其种群有时在形态上难以区分：Thymus mastichina (L.) L. 和 T. albicans Hoffmanns. and Link。Thymus mastichina是伊比利亚半岛分布最广的地方性百里香物种之一，已知会与大量分布区重叠的分类单元杂交。其内描述了两个亚种：广泛分布于伊比利亚的T. mastichina subsp. mastichina（四倍体）和主要位于西班牙西南部Do?ana国家公园和自然公园的T. mastichina subsp. donyanae（二倍体）。另一个物种是T. albicans（二倍体），发现于伊比利亚半岛西南海岸沿线，其栖息地受到沿海城市发展的高度威胁。分类学家主要基于花序直径和花萼长度等形态特征来区分这三个分类单元。然而，该组表现出高度的表型变异，包括形态和化学型水平。因此，阐明这个多倍体复合体错综复杂的关系——对于澄清其分类地位至关重要——应成为优先事项，这将有助于指导未来的保护策略。
随着高通量测序（High-Throughput Sequencing, HTS）技术变得越来越易于获取，群体基因组学方法已成为重建塑造复杂物种历史的进化模式的重要工具。在HTS技术中，基因组浅层测序（genome skimming）和靶向富集（target-enrichment）彻底改变了我们对植物系统学的理解。Hyb-Seq结合了这两种方法，由于其可承受性（合理的每样本成本）、可重复性和可复制性，同时无论组织模板质量如何都能提供高分辨率，因此在植物系统基因组学中获得了突出地位。虽然其他基因组技术，如限制性位点相关DNA测序（RADseq）或基因分型测序（GBS），已被广泛用于探索群体结构和遗传多样性，但Hyb-Seq在群体遗传学中的潜在应用在很大程度上仍未得到充分探索。然而，当前的系统基因组学和群体基因组学计算工具已被证明可以有效地连接微观和宏观进化水平，使用相同的Hyb-Seq数据，这些数据可以编码为核基因直系同源序列或单核苷酸多态性（SNPs）。通过结合这些分析方法，可以获得关于遗传连通性、基因渗入和群体微进化模式的见解，这对于识别面临近交或丧失适应性潜力风险的群体有所帮助，并有助于保护生物多样性管理。
本研究旨在理解百里香属Mastichina组多倍体复合体背后的进化动态和遗传结构，特别关注跨倍性混合的程度和模式。为实现此目标，我们采用Hyb-Seq数据结合系统基因组学和群体基因组学方法来（i）解析该组内分类单元的系统发育关系，（ii）评估复合体内二倍体和四倍体分类单元之间的遗传分化和混合，（iii）调查倍性水平内部和之间的基因流和基因渗入信号，包括导致跨倍性混合的潜在途径，以及（iv）评估这些模式对百里香属Mastichina组分类单元的分类学和保护的意义。这个整合框架旨在阐明在基因流跨越倍性水平存在的情况下，塑造地中海灌丛地方性植物物种边界的过程。
2 材料与方法
2.1 研究组和分类单元取样
百里香属Mastichina组包括两个物种和三个分类单元，均为伊比利亚半岛特有，其在该地区的分布模式各不相同。Thymus mastichina subsp. mastichina是一种四倍体分类单元（2n = 56, 58, 60），广泛分布于整个伊比利亚半岛的各种栖息地，如路边斜坡、废弃农田、松林和山麓碎石地，生长在疏松的硅质基质上，但也见于钙质土壤。Thymus mastichina subsp. donyanae是一种二倍体分类单元（2n = 30），主要位于西班牙西南部韦尔瓦省Do?ana国家和自然公园内的沿海沙洲和稳定沙丘。Thymus albicans也是一种二倍体分类单元（2n = 28, 30），栖息在伊比利亚半岛西南部加的斯省和阿尔加维海岸线上的稳定亚沿海沙丘和沿海松林中，生长在沙质基质上。Thymus albicans被列入西班牙国家和地区红色名录以及世界自然保护联盟红色名录。
我们在2023年5月至6月期间在伊比利亚半岛共取样了40个种群。其中20个种群与Girón等人研究的相同。每个种群的个体根据伊比利亚植物志进行鉴定。每个种群的标本凭证存放在SEV植物标本馆。每个种群8-12个个体的叶子在硅胶中干燥，并在DNA提取前室温保存。此外，我们还从MA植物标本馆的标本凭证中包括了另外三种百里香属物种作为外类群。
2.2 基因组大小和倍性水平评估
对来自40个研究种群的344个个体进行了核DNA含量评估。使用流式细胞术进行测量，以番茄和萝卜品种作为内标。使用LB01裂解缓冲液制备细胞核悬浮液，用碘化丙啶染色，并在流式细胞仪上测量。通过将标准品的DNA含量乘以样品和标准品在荧光强度直方图中2C峰位置的商来计算每个样品的总核DNA含量。对使用萝卜标准品获得的测量值应用了校正因子。新鲜材料和干燥叶子进行的基因组大小评估高度相关，并且两种方法分析的配对样本之间未发现显著差异。计算每个个体和每个种群的2C平均值和标准偏差。根据Girón等人的研究分配倍性水平。此外，我们还根据从基因组测序数据中调用SNP推断的等位基因比率估算了倍性水平。
2.3 分子方案
2.3.1 DNA提取和基因组文库制备
在每个种群的八个个体中进行DNA提取。使用商业试剂盒从硅胶干燥的叶组织中提取总细胞DNA。用超声破碎仪将DNA提取物打断至平均片段大小约为250 bp。使用商业试剂盒制备用于Illumina测序的双索引基因组文库。使用磁珠进行大小选择，索引PCR总共进行10个循环。使用荧光计检查所有基因组文库的浓度，使用安捷伦4200 TapeStation系统分析10%文库的片段分布。最后，将基因组文库标准化并混合。
2.3.2 Hyb-Seq实施
将文库池与被子植物353专家富集面板进行杂交。进行为期一天的杂交，然后使用链霉亲和素包被的磁珠捕获杂交的生物素标记的探针，并使用PCR试剂盒进一步富集16个循环。扩增的PCR产物用磁珠清洗，用荧光计定量，并在TapeStation系统上进行分析。然后将富集的池与普通基因组文库一起标准化和多重混合，以实现完整的Hyb-Seq实施。最后，该多重混合的Hyb-Seq池在Illumina NovaSeq X系统上测序，产生2×150 bp的双末端读长。
2.4 生物信息学分析
2.4.1 数据挖掘、靶标组装和数据矩阵优化
使用FastQC和MultiQC对解复用后的原始测序读长进行质量检查。使用fastp过滤掉低质量碱基、短读长和残留接头。使用HybPiper组装核靶标和质体编码DNA序列，使用自定义的质体CDS靶标文件。平均每个个体约有11M读长，其中约1.18M读长映射到预期的353个核靶标。靶标成功率平均约为10.5%。未检测到潜在的副同源性。使用脚本移除低覆盖度的样本和基因。最终保留了304个样本和332个超级重叠群序列用于后续分析。
2.4.2 多序列比对和过滤
使用MAFFT对平衡的数据矩阵进行比对，并使用AMAS进行质量检查。使用FastTree从生成的MSA推断探索性的基因树。应用TreeShrink自动从FastTree拓扑及其相应的MSA中移除异常分支。使用AMAS选择最大化简约信息字符比例同时最小化缺失数据百分比的阈值。排除PIC值低于中位数三分之一的靶标直系同源物。最终保留了301个个体和329个超级重叠群用于下游分析。最后，使用MAFFT重新比对选定的经过异常值过滤的MSA，并使用trimAl进一步优化。
2.4.3 基因树和物种树推断
使用IQ-TREE估算核ML基因树。使用ModelFinder Plus选择最佳拟合的序列进化模型，并使用UFBoot2和SH-aLRT算法计算支持值。使用Newick实用程序包折叠ML基因树中不支持的分支二分法。使用ASTRAL III在多重物种合并框架下估算物种树，并以局部后验概率估算分支支持。最后，使用RAxML-NG将MSC物种树缩放以呈现每位点的替换分支长度。使用FigTree可视化缩放后的物种树。由于质体组代表一个典型的合并基因，使用AMAS连接经过异常值过滤、修剪后的超级重叠群MSA。然后使用IQ-TREE推断质体ML树。
2.4.4 SNP调用
对所有包含在系统基因组分析中的百里香属Mastichina组样本进行SNP调用，采用了改进的工作流程。排除了外类群，因为此分析旨在检查仅在研究群体内分离的遗传变异。使用GATK4进行变异检测，使用每个353个靶标的最长超级重叠群序列作为伪基因组参考。我们选择了保守的SNP调用方法，消除了多等位基因变异。最初调用了总共82,137K个单核苷酸多态性。经过严格过滤后，使用BCFtools和VCFtools处理SNP矩阵。使用PLINK通过命令避免连锁不平衡。最终的过滤SNP数据矩阵包含254个个体和2223个未连锁的SNP。
2.4.5 基于基因组等位基因比率的倍性水平估算
我们按照设计的工作流程估算样品的倍性水平。使用nQuire估算等位基因频率。使用R包计算每个样品的等位基因比率和delta对数似然值。然后使用基于模型的方法推断倍性，选择具有最低对数似然值的模型对应的倍性水平。将这些估计值与流式细胞术测量值进行比较。
2.4.6 群体结构分析和汇总统计
使用上述SNP矩阵，在ADMIXTURE中进行群体结构分析，并基于交叉验证误差估计最佳聚类数。同时，使用STRUCTURE估计最佳遗传聚类数。测试了一系列K值，并通过检测ΔK统计量的最高值来选择最可能的聚类数。使用R包绘制每个K值的ADMIXTURE结果。使用R包获得饼图并将其投影到从QGIS获取的地图上。使用R包进行主成分分析和判别分析。还使用R包可视化结果。
使用PLINK计算检测到的系统发育分组之间和每个种群的固定指数、近交系数和成对身份同源值。使用R包进行分子方差分析，测试系统发育分组内部和之间的遗传变异性。使用R包进行Mantel检验，基于Pearson乘积矩相关。使用VCFtools计算核苷酸多样性和Tajima's D。为了识别系统发育分组之间私有的和共享的SNP，我们使用R包。分析使用每个分组11个个体的随机子样本进行，以标准化分组大小，从而减少由于分组大小不均导致的私有SNP数量检测偏差。
2.4.7 基因渗入和基因流分析
使用SplitsTree6重建了一个邻接网络系统发育网络，使用了用于系统基因组分析的连接基因矩阵和过滤后的SNP矩阵。该分析使我们能够可视化样本之间潜在的网状进化事件。此外，使用Dsuite在所有遗传分组之间进行了ABBA-BABA检验，该检验依赖于Patterson's D统计量来估计两个给定分组之间共享衍生等位基因的全基因组过剩，以检验基因渗入事件。
为了进一步评估跨倍性基因流的存在并量化其方向性，我们应用了fastsimcoal2，使用位点频率谱文件和easySFS软件。我们比较了四种情景：（1）无基因流，（2）从四倍体到二倍体的单向基因流，（3）从二倍体到四倍体的单向基因流，以及（4）双向基因流。我们进行了50次运行，使用ECM优化循环和合并模拟估算参数。基于最低中位AIC值选择最佳拟合模型，并使用似然分布进行验证。最后，我们使用具有100次重复的块自助法来计算参数估计的置信区间，确保我们发现的可靠性。
3 结果
3.1 基因组大小测量
流式细胞术估算的DNA含量显示每个分析个体的平均2C值范围从1.32到3.07 pg。Girón等人研究的种群的平均2C值比我们研究获得的平均值低6%，尽管这些差异在DNA量上较低。采样个体的平均2C值呈现双峰分布：2C值在1.32至1.84 pg范围内的样本被认为是二倍体，而平均2C值在2.31至3.07 pg范围内的样本被认为是四倍体。
3.2 百里香属Mastichina组的系统基因组关系
生成的物种树显示百里香属Mastichina组分为两个系统发育群，分别对应于二倍体和四倍体样本。二倍体群由被鉴定为T. albicans或T. mastichina subsp. donyanae的样本组成，以及在一些种群中被鉴定为T. mastichina subsp. mastichina的样本。四倍体群由被鉴定为T. mastichina subsp. mastichina的样本组成，分布在整个伊比利亚半岛，但花萼长度在个体内部和个体之间具有高度变异性，也符合T. mastichina subsp. donyanae和T. albicans关于花萼长度的描述。二倍体群揭示了两个姐妹亚群。一方面是Hercynian亚群，由表型上符合T. mastichina subsp. mastichina的二倍体个体组成，来自两个地理上不同的区域——埃什特雷拉山脉自然公园和科尔多瓦。另一方面是其余的二倍体种群，全部位于伊比利亚半岛西南部，它们进一步细分为三个支持度较弱但具有地理结构的亚群。Algarve亚群包括来自阿尔加维地区的个体，表型符合T. albicans。Do?ana亚群主要由来自Do?ana国家和自然公园的个体组成，表型对应于T. mastichina subsp. donyanae，以及附近的种群，其中一些呈现出T. mastichina subsp. donyanae和T. albicans之间的中间特征，而该组最西端的种群符合T. mastichina subsp. mastichina。最后，Cádiz亚群主要由来自加的斯省的个体组成。基于72个质体CDS的质体基因组树显示个体的混合分布，没有明确的分组。
3.3 基于等位基因比率的倍性水平估计
每个个体的估计等位基因比率在1到5.5之间，二倍体群的中位值从1.24到2.67，四倍体群从1.77到3.04。基于模型的方法将每个个体的中位等位基因比率分为三个水平：低、中、高。这三个等位基因比率水平分别符合nQuire管道专门测试的三倍体模型：二倍体、三倍体和四倍体。在二倍体群中发现了树类别，在后两者中发现了四倍体群。2C DNA含量在基于等
