选择性清除检测 — iHS 与 EHH 统计量¶
一句话概述：iHS（integrated Haplotype Score）通过比较祖先等位基因和衍生等位基因周围单倍型的纯合度衰减速度，检测正在进行的正向选择（incomplete sweep）。
核心知识点速查表¶
| 概念 | 白话解释 |
|---|---|
| 正向选择（Positive Selection） | 有利突变在群体中频率快速上升 |
| 选择性清除（Selective Sweep） | 有利突变"拖着"周围的中性变异一起上升，降低周围多样性 |
| EHH | 扩展单倍型纯合度，衡量从一个位点向两侧纯合度衰减的速度 |
| iHS | 比较祖先/衍生等位基因的 EHH 积分比值，检测进行中的选择 |
| nSL | 类似 iHS 但不需要遗传图谱，用 SNP 数量替代物理距离 |
| XP-EHH | 跨群体 EHH，检测群体间差异选择（一个群体选完，另一个没有） |
| selscan | 计算 iHS/nSL/XP-EHH 的主流软件（2024 更新到 v2.0） |
| rehh | R 包，功能全面但速度不如 selscan |
一、选择性清除的原理（白话版）¶
比喻：一个村庄里突然出现一个超能力者（有利突变），他很快娶妻生子扩散后代。由于繁殖速度太快，他周围的邻居（连锁的中性变异）也跟着"沾光"传播了。结果就是：超能力者所在的那条染色体片段上，所有人的 DNA 都长得很像（纯合度高、多样性低）——这就是选择性清除。
iHS 的逻辑：如果一个 SNP 位点的衍生等位基因（新突变）正在被正向选择，那么携带它的单倍型还来不及被重组打断，呈现出又长又完整的特征。iHS 就是比较"衍生等位基因周围的单倍型长度"和"祖先等位基因周围的单倍型长度"。
二、分析流程¶
整体流程¶
第 1 步：数据准备¶
# 安装 selscan（主流选择性清除检测工具）
conda install -c bioconda selscan # conda安装selscan
# selscan 2.0 (2024) 新功能：支持未分型数据（--unphased）
# 但分型后的数据效果更好
# 如果VCF还没有分型（phasing），需要先分型
# 用 SHAPEIT4 或 Beagle 进行分型
conda install -c bioconda shapeit4 # 安装SHAPEIT4
# 运行 SHAPEIT4 分型
shapeit4 --input input.vcf.gz \ # 输入VCF
--map genetic_map.txt \ # 遗传图谱（cM/Mb）
--region chr1 \ # 指定染色体
--output phased_chr1.vcf.gz # 输出分型后的VCF
第 2 步：准备遗传图谱¶
# iHS 需要遗传图谱（genetic map）来衡量遗传距离
# 下载 HapMap 遗传图谱（人类数据）
wget https://ftp.ncbi.nlm.nih.gov/hapmap/recombination/2011-01_phaseII_B37/genetic_map_HapMapII_GRCh37.tar.gz
tar -xzf genetic_map_HapMapII_GRCh37.tar.gz # 解压
# 遗传图谱格式（selscan要求）：
# 染色体 位置(bp) 遗传距离(cM) 遗传位置(cM)
# 1 55550 0.000000 0.000000
# 1 568322 0.008000 0.008000
# 如果没有遗传图谱，可以用 nSL 替代 iHS（nSL不需要图谱）
第 3 步：运行 selscan 计算 iHS¶
# 计算 iHS（需要分型数据 + 遗传图谱）
selscan --ihs \ # 计算iHS统计量
--vcf phased_chr1.vcf.gz \ # 输入分型后的VCF
--map genetic_map_chr1.txt \ # 遗传图谱
--out chr1_ihs \ # 输出前缀
--threads 4 # 使用4个线程
# 输出文件：chr1_ihs.ihs.out
# 格式：位点ID 物理位置 衍生等位基因频率 iHH_ancestral iHH_derived unstd_iHS
# 如果数据未分型（selscan 2.0+）
selscan --ihs \
--unphased \ # 使用未分型模式
--vcf unphased.vcf.gz \
--out chr1_ihs_unphased \
--threads 4
第 4 步：计算其他统计量（nSL 和 XP-EHH）¶
# nSL：不需要遗传图谱，适合没有图谱的物种
selscan --nsl \ # 计算nSL统计量
--vcf phased_chr1.vcf.gz \ # 输入VCF
--out chr1_nsl \ # 输出前缀
--threads 4
# XP-EHH：跨群体比较，检测群体间差异选择
selscan --xpehh \ # 计算XP-EHH
--vcf pop1_chr1.vcf.gz \ # 群体1的VCF
--vcf-ref pop2_chr1.vcf.gz \ # 群体2的VCF（参考）
--map genetic_map_chr1.txt \ # 遗传图谱
--out chr1_xpehh \ # 输出前缀
--threads 4
第 5 步：标准化（normalization）¶
# iHS 原始值受等位基因频率影响，需要按频率分 bin 标准化
# norm 是 selscan 自带的标准化工具
# 标准化 iHS
norm --ihs \ # 指定统计量类型为iHS
--files chr1_ihs.ihs.out \ # 输入原始结果
--bp-win # 启用基于窗口的标准化
# 标准化 nSL
norm --nsl \
--files chr1_nsl.nsl.out
# 标准化 XP-EHH
norm --xpehh \
--files chr1_xpehh.xpehh.out
# 输出：*.ihs.out.100bins.norm
# 格式增加最后一列：标准化后的 iHS 值
# |iHS| > 2 表示显著（约 p < 0.05）
第 6 步：结果可视化¶
# ========== R 脚本：Manhattan plot ==========
# 读取标准化后的 iHS 结果
ihs <- read.table("chr1_ihs.ihs.out.100bins.norm", # 读取文件
header=FALSE) # 无表头
# 列名对应
colnames(ihs) <- c("SNP", "pos", "freq", # 前3列
"iHH_anc", "iHH_der", # 中间2列
"unstd_iHS", "std_iHS", # 标准化前后的iHS
"significant") # 是否显著
# 画 Manhattan plot
pdf("iHS_manhattan.pdf", width=12, height=5) # 创建PDF
plot(ihs$pos / 1e6, # X轴：位置(Mb)
abs(ihs$std_iHS), # Y轴：|iHS|绝对值
pch=20, cex=0.3, # 小点
col=ifelse(abs(ihs$std_iHS) > 2, # 显著位点用红色
"red", "gray70"),
xlab="Position (Mb)", # X轴标签
ylab="|iHS|", # Y轴标签
main="Integrated Haplotype Score (iHS)") # 标题
abline(h=2, col="blue", lty=2) # 显著性阈值线
dev.off() # 保存
# ========== 用 rehh 包可视化 EHH 衰减曲线 ==========
library(rehh) # 加载rehh包
# 从VCF读取数据
hap <- data2haplohh(hap_file="phased_chr1.vcf.gz", # 输入VCF
polarize_vcf=TRUE, # 根据AA注释极化
chr.name="1") # 染色体名
# 计算特定位点的 EHH 衰减曲线
ehh <- calc_ehh(hap, # 输入数据
mrk=which(hap@positions==12345678), # 指定位点位置
include_nhaplo=TRUE) # 包含单倍型数量
# 画 EHH 衰减图
plot(ehh, main="EHH decay at rs12345", # 画图
col=c("blue","red")) # 祖先=蓝，衍生=红
# 如果红色线（衍生等位基因）比蓝色线衰减慢，说明衍生等位基因正在被选择
三、结果解读¶
iHS 值的含义¶
| |iHS| 值 | 含义 | |----------|------| | < 1 | 无选择信号 | | 1 ~ 2 | 弱信号（可能有选择） | | 2 ~ 3 | 中等信号（显著，约 top 5%） | | > 3 | 强信号（高度显著，约 top 1%） |
注意： - iHS > 0（正值）→ 祖先等位基因的单倍型更长 → 选择作用在衍生等位基因上 - iHS < 0（负值）→ 衍生等位基因的单倍型更长 → 选择作用在祖先等位基因上 - 实际分析中通常取绝对值 |iHS|
iHS vs nSL vs XP-EHH 的选择¶
| 统计量 | 适用场景 | 是否需要图谱 |
|---|---|---|
| iHS | 检测进行中的 sweep（频率 0.2~0.8） | 需要 |
| nSL | 同 iHS，适合没有图谱的物种 | 不需要 |
| XP-EHH | 检测群体间差异选择（频率可接近 fixation） | 需要 |
四、常见报错与解决¶
| 报错信息 | 原因 | 解决方法 |
|---|---|---|
ERROR: VCF is not phased | VCF 文件没有分型信息 | 用 SHAPEIT4/Beagle 分型，或用 --unphased |
No map file provided | iHS 需要遗传图谱 | 提供图谱文件，或改用 nSL |
Segmentation fault | 染色体太大/内存不足 | 按染色体分开计算 |
All values are NaN after normalization | SNP 太少 | 增加 SNP 密度或放宽过滤条件 |
EHH does not decay to 0 | 重组率低或样本量小 | 增加 --max-extend 参数 |
| selscan 比 rehh 结果略不同 | 算法实现细节不同 | 正常现象，差异通常很小 |
五、面试高频问题¶
Q1：iHS 能检测到已经固定的选择信号吗？
不能。iHS 检测的是进行中的选择（incomplete sweep），即有利等位基因频率在 0.2~0.8 之间。已经固定（频率=1）的信号需要用 XP-EHH 或其他方法（如 Tajima's D、CLR）。
Q2：iHS 和 Fst 检测选择的区别？
iHS 是单群体内检测，基于单倍型纯合度，灵敏度高但需要分型数据。Fst 是群体间比较分化程度，不需要分型但不能区分选择方向。两者互补，通常一起用。
Q3：selscan 2.0 的主要改进是什么？
2024 年发布，最大改进是支持未分型数据（
--unphased），将基因型编码为 0/1/2 来计算 multi-locus genotype homozygosity。此外 2025 年的动态规划算法（v2.1）将速度提升 20~50 倍。
六、速查表¶
# === 选择性清除检测速查 ===
# 1. 分型（如果需要）
shapeit4 --input raw.vcf.gz --map map.txt --region chr1 --output phased.vcf.gz
# 2. 计算 iHS
selscan --ihs --vcf phased.vcf.gz --map map.txt --out ihs --threads 4
# 3. 计算 nSL（不需要图谱）
selscan --nsl --vcf phased.vcf.gz --out nsl --threads 4
# 4. 计算 XP-EHH（跨群体）
selscan --xpehh --vcf pop1.vcf.gz --vcf-ref pop2.vcf.gz --map map.txt --out xpehh
# 5. 标准化
norm --ihs --files ihs.ihs.out
norm --nsl --files nsl.nsl.out
norm --xpehh --files xpehh.xpehh.out
# 6. 筛选显著位点（|iHS| > 2）
awk 'sqrt($7^2) > 2' ihs.ihs.out.100bins.norm > significant_ihs.txt
# 显著性阈值参考：
# |iHS| > 2 → top ~5%
# |iHS| > 2.5 → top ~1%
# |iHS| > 3 → top ~0.1%
参考资料：selscan 2.0 (Szpiech 2024, Bioinformatics)、Rahman et al. 2025 (MBE)、rehh R包文档