基于线粒体COI基因的中国华山松大小蠹（Dendroctonus armandi）谱系地理学研究：揭示晚更新世以来的种群分化与扩张路径

《Genes》：Phylogeography of Chinese White Pine Beetle Dendroctonus armandi (Coleoptera: Curculionidae: Scolytinae) in China Hang Ning, Ruixiong Deng, Kaitong Xiao, Beibei Huang, Yu Cao and Qiang Wu

本综述通过分析中国特有华山松大小蠹（Dendroctonus armandi）的线粒体DNA（mtDNA）COI基因，揭示了该害虫在晚更新世期间的遗传分化、种群扩张与扩散路径。研究识别了与岷山（MSM）、秦岭（QLM）、米仓山（MCM）和大巴山（TPM）相对应的四个单倍型群，并推测了三条独立的迁殖路径，为理解气候变化、地理隔离与寄主分布对昆虫谱系地理格局的影响提供了关键见解。

　　

1. 1.
   引言
   在过去的300万年中，全球气候经历了剧烈波动。第四纪更新世时期冰期与间冰期的交替循环对物种分布和种群分化产生了关键性影响。在亚洲，冰川事件对物种演化的影响与欧洲、北美存在差异。在中国，复杂的地质结构与局地特征使得物种在第四纪冰期后表现出独特的进化模式。
   作为针叶林的重要害虫，华山松大小蠹（Dendroctonus armandi）专性寄生华山松（Pinus armandii），其分布局限于中国中部山区，造成严重的生态与经济损失。然而，驱动其持续发生的内部因素尚不清楚。已有研究表明，秦岭和大巴山地区（作为古北界与东洋界的分界，也是温带与亚热带的气候分界线）的物种常呈现南北（纵向）和东西（横向）的遗传分化模式。该地区在冰期可能为许多物种提供了避难所。本研究旨在利用线粒体DNA（mtDNA）COI基因序列，探究华山松大小蠹的遗传结构、估算不同谱系的分化时间，并确定其种群迁殖与扩散路线，从而构建其谱系地理学的基本框架。
2. 2.
   材料与方法
   2.1. 样品采集与DNA扩增
   研究团队在2024年至2025年间，从其整个地理分布范围内的38个地点采集了255头华山松大小蠹成虫样品。所有个体均从受侵染的树木中采集，保存于95%乙醇中，并于-80°C储存直至DNA提取。使用DNeasy组织试剂盒（Qiagen）提取DNA。mtDNA COI基因的扩增与测序使用引物C1-J-2441和T12-N-3014进行。扩增产物纯化后，在ABI Prism 3100遗传分析仪上进行测序。所有获得的核苷酸序列经比较、手工编辑，并使用Sequencher^?4.7软件进行比对。最终获得47个单倍型，用于后续分析。

2.2. 系统发育重建

使用aLRT-PHYML程序中的最大似然法（ML）重建mtDNA单倍型间的系统发育关系。以Dendroctonus simplex作为外群。利用DnaSP v6.12.03软件估算各单倍型群的分离位点、单倍型数量、单倍型多样性及平均核苷酸差异等参数。同时，采用Tajima’s D、Fu and Li’s D和Fu’s F_S统计检验来识别可能指示自然选择或种群扩张的突变-漂变平衡偏离。

2.3. 种群历史动态

采用两种方法分析有效种群大小的时间变化：中性检验（检验种群增长）和配对差异分布。其中，Tajima’s D的负值可揭示近期的人口扩张。Fu’s F_S和Ramos-Onsins and Rozas’s R₂检验也用于验证扩张信号。此外，利用不匹配分布检验突发性种群扩张。假设为阶梯扩张模型，估算参数θ₀、θ₁和τ。同时计算Harpending’s 粗糙指数（rg）。在模型有效性得到确认后，利用τ估算扩张事件发生的时间，其中昆虫世代时间设为1.5年，分歧速率采用每百万年1.1%。

2.4. 种群间分化时间估算

使用MDIV v1.0软件估算华山松大小蠹各单倍型群之间的分化时间和迁移率。该软件能够联合估计两个种群的θ（=2N_eμ）和T（种群间分化时间，时间单位=N_e代）。采用HKY模型，运行Metropolis-Hastings Markov chain Monte Carlo模拟500万代，其中前100万代作为burn-in。根据后验分布的众数估算θ和T，并利用公式T_pop= [(t_popθ)/2L] 1/μ计算种群分化时间，其中L为序列长度（596 bp），μ为每个位点每代的突变率。同样采用节肢动物mtDNA常用的每百万年1.1%的分歧速率。

1. 3.
   结果
   3.1. 华山松大小蠹单倍型的分子谱系
   对线粒体COI基因596个核苷酸的测序共识别出142个可变位点和47个不同的单倍型。最大似然分析构建的系统发育树（）显示，华山松大小蠹形成了四个主要的单倍型群，分别对应于岷山（MSM）、秦岭（QLM）、米仓山（MCM）和大巴山（TPM）。其中，位于秦岭中部的佛坪县拥有最多的单倍型数量（7个）。尽管部分分支节点的支持值在70至82之间，但由于aLRT检验的可靠性，这些分支仍被认为是可靠的。
   所有单倍型群均表现出较高的单倍型多样性（Hd = 0.934–0.952）。QLM单倍型群的遗传多样性最高（S = 57, Nhap = 17, π = 0.047, k = 27.8），并伴有显著为负的Fu and Li’s D值（-3.892）。结合其单峰不匹配分布，证实了该群近期发生过突发性种群扩张。相反，MSM和MCM的遗传多样性相对较低（S = 23/27, π = 0.019/0.022），Fu and Li’s D值为弱负值，表明其进化历史稳定，无明显种群扩张。MSM和MCM的Tajima’s D和Fu’s F_S值也为不显著的弱负值，反映了与突变-漂变平衡无偏离。TPM单倍型群显示出中等遗传多样性，与其温和的种群扩张信号一致。

3.2. 分化时间

各单倍型群间的种群分化时间在0.265至0.492个百万年之间。最深的分化出现在QLM与TPM之间（0.492 百万年），最浅的分化出现在MSM与MCM之间（0.265 百万年）。所有最近共同祖先（T_MRCA）的时间估计值（0.289–0.518 百万年）均略高于相应的种群分化时间（T_pop），符合祖先谱系合并先于种群分化的预期遗传模式。QLM-TPM对的θ（2.23–3.68）和有效雌性种群大小（N_e♀: 2401–4521）估计值最高，反映了这两个扩张单倍型群之间更大的遗传分化和更大的历史种群规模。

这些分化时间估计值均落在晚更新世时期，与在QLM和TPM种群中检测到的人口扩张信号（包括显著为负的Tajima’s D和Fu’s F_S值、低SSD和rg、单峰不匹配分布、θ₀≈ 0和正的τ估计值）相一致。相比之下，MSM和MCM表现出稳定的种群动态，中性检验不显著，SSD和rg值较高，不匹配分布平缓，这与其较浅的分化时间相符。

1. 4.
   讨论
   分子数据分析表明，华山松大小蠹包含四个主要谱系，分别分布于岷山、秦岭、米仓山和大巴山地区，并拥有几乎独一无二的单倍型集合和明显的地理分布。先前对中国潜在第四纪冰川避难所的研究已识别出多个区域，包括横断山、秦岭和大巴山。本研究中所有华山松大小蠹种群也局限在这些地理区域，并沿岷山和秦巴山脉自西向东分布。受青藏高原隆升的影响，该地区的地形地貌变得极为复杂。地理隔离阻碍了基因流，并产生了遗传分化，这是限制华山松大小蠹独特地理分布的关键因素。
   单倍型和核苷酸多样性结果表明，MSM、QLM、MCM和TPM是华山松大小蠹种群的古老避难所。位于青藏高原东北部的岷山（属于横断山脉）可能是在冰河时期幸存下来，随后向周边地区扩张的来源地。QLM的单倍型和核苷酸多样性高于MCM，因此推测华山松大小蠹首先自西向东传播并定殖于QLM。此外，该害虫可能通过另一条自西向东的路线从MSM扩散到MCM。然而，TPM（直至神农架）的单倍型和核苷酸多样性高于MCM。类似的结果在其他北美和欧洲的树皮甲虫中也有发现，许多分布区扩张事件都与冰川期后从避难所的迁殖有关。
   系统发育拓扑结构显示，支系I和II处于基础位置且分支长度相似，表明华山松大小蠹各单倍型群的初始分化可能在时间上接近。来自广东茂名盆地上更新世的Pinus L.木化石表明，耐寒植物华山松在末次冰期曾扩张到华南低纬度地区，随后在间冰期再次收缩。同时，分布于秦岭、大巴山以及中国中西部云贵高原的华山松受到冰期降温效应的影响，使得喜凉的华山松能够扩散到较低海拔地区。冰期结束后气候变暖，其分布逐渐退缩到合适的高海拔栖息地。
   本研究通过中性检验和不匹配分布两种方法，专门检验了假设的人口扩张。中性检验表明，来自QLM和TPM的种群经历了人口扩张事件。结合所有指标，结果强有力地支持QLM和TPM单倍型群发生人口扩张的假设。相反，MSM和MCM表现出稳定的种群规模，没有近期人口扩张的证据。所有山地单倍型群推定扩张事件的估计时间均在晚更新世，可能是在过去冰川作用导致的严重瓶颈期之后。
   青藏高原的持续隆升是建立东亚气候模式的关键。自第四纪以来，周边山脉受到冰川、流水侵蚀等外力的影响，最终形成了现今的地形。无疑，气候波动和地理隔离可能是华山松大小蠹分化的原始驱动力。QLM和TPM单倍型群之间的分化时间估计为晚更新世的0.492 百万年。一种可能的解释是，这两个单倍型群内的谱系经历了显著分化，但随着时间的推移存在一些间歇性接触，并且在末次冰盛期之前就持续向东扩张至神农架。TPM内的谱系因气候变化直接导致了寄主可获得性的差异。在这种情况下，TPM地区华山松的东部范围可能提供了多个避难所，而冰川覆盖了其大部分西部范围。
   总之，许多因素影响着物种的分布格局（遗传结构），包括采样个体的分布、采样原则、所研究基因片段的进化速率、以及自然选择或种群的历史动态变化。此外，还需考虑环境和人为干扰引起的种群隔离或分布碎片化。除了这些，长期的气候振荡、地质事件、地形异质性和栖息地变化会促进长期的种群隔离并限制基因流。同时，加剧的人为干扰（如森林砍伐、农业扩张、城市化、道路建设、过度开发和栖息地丧失与破碎化）进一步破坏了种群连通性，降低了有效种群规模，并加速了遗传分化。这些自然和人为过程共同塑造了物种当代的遗传结构和分布模式。
2. 5.
   结论
   对来自38个种群的255个个体进行的mtDNA分析，识别出45个具有地理结构的单倍型，它们形成了对应于MSM、QLM、MCM和TPM的四个单倍型群。在这四个单倍型群中，QLM和TPM单倍型群经历了人口扩张，所有单倍型群间的分化时间均可追溯到晚更新世。华山松大小蠹展现出三条独立的迁殖路线：两条横向的（MSM→QLM；MSM→MCM→TPM）和一条纵向的（QLM→TPM）。除了气候波动和地理隔离外，其寄主华山松的斑块状分布进一步强化了这种谱系地理格局。这些结果增进了对华山松大小蠹在山地系统中进化与适应的理解。