从低覆盖率的基因组序列中准确估计非模式物种的纯合子比例

　　基因组自交片段（ROH）分析在低测序深度非模式物种中的应用研究。通过优化ROHan软件的rohmu参数和窗口大小，发现3-5倍测序深度即可准确估计ROH比例和长度，但需根据种群历史调整参数。研究纠正了北极驯鹿种群Gaspésie的近交水平，揭示Plink低估ROH的问题，并提出了适用于保护遗传学的低深度分析策略。

　　
羔鹿种群遗传多样性评估与ROH分析工具优化研究

摘要：
随着测序成本的持续降低，野生动物种群的全基因组测序数据量显著增加。本研究聚焦于利用低测序深度数据（3-15倍覆盖）准确评估羔鹿种群近交水平（Runs of Homozygosity, ROH）的可行性。通过对比Plink和ROHan两种工具在不同测序深度下的表现，揭示了传统工具在低覆盖场景下的系统性偏差。研究发现，ROHan工具通过优化关键参数rohmu（容忍的杂合度阈值）和窗口大小（100kb与1Mb），可在测序深度3-5倍时仍保持较高的ROH估计精度。特别针对加拿大Gaspésie隔离亚群的研究表明，使用Plink工具可能导致近交水平被低估达87%，而经过参数优化的ROHan分析使检测准确率提升至92%。本研究为野生动物种群遗传监测提供了新的方法学框架。

1. 研究背景与科学问题
1.1 ROH在种群遗传学中的重要性
ROH指连续多代近交导致的同源纯合区，其分布特征与种群有效规模（Ne）、迁移率（M）及历史瓶颈效应密切相关。研究显示，ROH与近交衰退（inbreeding depression）存在显著正相关，且特定长度的ROH可反映近交事件的年代（Allendorf, 2017）。在濒危物种保护中，准确量化ROH对评估遗传多样性损失、制定种群恢复策略具有决定性意义。

1.2 现有方法的局限性
当前主流工具Plink在低测序深度（<10×）下存在显著偏差。其依赖的滑动窗口算法易受测序错误、低质量位点（Q<20）和SNP密度不足的干扰（Duntsch et al., 2021）。比较研究表明，Plink对ROH的误判率可达30%-50%（Silva et al., 2024），特别是在窗口长度>500kb时，容易合并短ROH形成虚假长ROH，导致近交水平高估（Kardos & Waples, 2025）。

1.3 ROHan工具的创新性
ROHan作为专为低深度数据设计的工具，采用动态权重算法和两状态隐马尔可夫模型（HMM），在处理古DNA（测序深度1-2×）时仍能保持较高准确性（Renaud et al., 2019）。本研究首次将ROHan应用于野生动物，通过22个高覆盖度（15×）羔鹿样本建立基准数据，进而模拟不同测序深度下的表现。

2. 研究方法与实验设计
2.1 数据来源与预处理
采用加拿大国家生物多样性冰冻库保存的22个高覆盖羔鹿样本（平均15×），包含从Baffin岛到George河的多个地理种群。样本涵盖不同遗传状态： Eastern Migratory（低近交）、Barren-ground（中等近交）、Boreal（强烈近交）和Neys Area（极端近交）等类型。数据预处理包括：
- BAM文件生成：通过GATK4进行基因组拼接，BamUtil进行重叠区域裁剪
- 筛选标准：保留35个主要 scaffolds（覆盖99%基因组），过滤Q<20和低深度位点（<7×）
- 深度降采样：从15×逐步降至1×，保持完整SNP数据集

2.2 工具参数优化
通过基准数据（高覆盖样本）对比Plink和ROHan的输出，建立优化流程：
- 窗口长度选择：100kb（短ROH检测）与1Mb（长ROH检测）双模式
- rohmu参数优化：测试2e-3至5e-5六个参数值
- 混合模型验证：比较基于不同rohmu值的模型似然比（LR值）

3. 关键研究发现
3.1 参数敏感性分析
- rohmu参数对ROH估计影响呈梯度变化：当rohmu=8e-4时，ROH面积与Plink 15×数据结果偏差<5%；当rohmu>1e-4时，ROH面积被高估达40%-60%（图2a,b）
- 窗口长度效应显著：100kb窗口能检测到平均长度1.2Mb的ROH，而1Mb窗口将短ROH合并导致低估约30%（表1）

3.2 测序深度与误差分析
- 误差类型分布：
- 测序深度<5×：主要误差类型为假阳性杂合（False Positive Hets）导致ROH高估（平均误差+18%）
- 5-10×深度：存在系统偏差（平均误差+12%），与Plink相比ROHan的theta值偏差降低至8%
- >10×深度：两种工具结果趋同（R2=0.92）
- 窗口大小与深度交互作用：
- 小窗口（100kb）在3×深度仍能检测到有效ROH（置信区间±5%）
- 大窗口（1Mb）在5×深度即可识别长ROH（>1Mb）

3.3 Gaspésie种群案例研究
- 基于Plink的15×数据，该种群近交水平被低估至1.1%
- ROHan优化分析（rohmu=8e-4，100kb窗口）显示：
- ROH面积达36.49% ± 5.2%
- 中等长度ROH（1-10Mb）占比达28.7%
- theta值（包含ROH）为0.0082（95% CI:0.0075-0.0090）
- 与相似地理的Boreal种群对比：
- ROH面积差异达14.2%（Gaspésie vs. Lake Superior）
- theta值差异显著（p<0.01）

4. 方法优化建议
4.1 参数选择指南
- 对于短ROH主导的种群（如Neys Area）：推荐rohmu=5e-4 + 100kb窗口
- 对于长ROH主导的种群（如Greenland）：推荐rohmu=2e-4 + 1Mb窗口
- 一般场景：rohmu=8e-4（100kb）和rohmu=2e-4（1Mb）组合使用

4.2 多窗口协同分析
建议采用双窗口策略（100kb+1Mb）：
- 100kb窗口检测短ROH（<1Mb），用于评估近期近交
- 1Mb窗口检测长ROH（>1Mb），用于识别历史瓶颈效应
- 结合使用可同时获得短/长ROH的置信区间（95% CI）

4.3 数据预处理标准流程
- 基因组组装：保持≥99%组装覆盖率
- 筛选标准：
- 深度阈值：min depth=3×（ROHan） / 7×（Plink）
- 质量阈值：min Q=20（保留10^4概率）
- SNP密度：≥1 SNP/500kb（ROHan优化条件）
- 误差校正：
- 使用高覆盖数据建立基准参数
- 通过模拟数据验证参数稳定性（图S1）
- 校正因子：1 - (实测ROH面积 - 模拟ROH面积)/模拟ROH面积

5. 管理应用与未来方向
5.1 濒危种群评估修正
- Gaspésie羔鹿案例显示：传统方法低估近交水平达87%
- 建议将ROH阈值从Plink的5%调整至8%（置信区间±2%）
- 每年监测ROH面积变化率（ΔROH/Δt）作为种群适应性的新指标

5.2 工具性能对比
| 工具 | 深度适应性 | 准确率（5-10×） | 计算效率 |
|------------|------------|----------------|----------|
| Plink | <10×不推荐 | 72% | ★★★★★ |
| ROHan | ≥3× | 89% | ★★★★☆ |
| RZooROH | ≥5× | 85% | ★★★☆☆ |

5.3 技术推广建议
- 建立"基准-优化"分析流程：
1. 高覆盖度数据（≥10×）确定基准ROH参数
2. 降采样至目标深度（3-15×）
3. 通过基准参数计算适配的rohmu值
- 开发自动化校准脚本（GitHub: https://github.com/BeckySTaylor/Phylogenomic_Analyses）
- 建议最小测序深度标准：单个体≥3×，群体≥5×

5.4 研究局限与展望
- 当前模型未考虑群体结构对ROH检测的影响（需结合EHH分析）
- 古DNA样本的突变率校正模型缺失
- 计算资源需求：1Mb窗口需20GB内存，建议采用分布式计算框架
- 计划开发多窗口并行分析模块（预计2025年Q3完成）

本研究为野生动物种群遗传监测提供了可靠的工具链和参数体系。通过建立"高覆盖基准-低覆盖优化"的标准化流程，可在测序深度≥3×时实现ROH估计的95%置信区间误差<8%。建议自然保护机构优先采用ROHan工具进行种群遗传评估，特别是在监测隔离亚群和濒危种群时，需结合双窗口分析（100kb+1Mb）和参数动态优化，以准确反映近交水平及遗传多样性现状。