759. 全基因组duplications检测¶

一句话概述：检测基因组中拷贝数变异(CNV)，找出哪些基因片段多了(duplication)或少了(deletion)——就像检查一本书有没有被多印了几页或撕掉了几页。

核心知识点速查表¶

概念	白话解释	关键工具
CNV	拷贝数变异(>50bp的缺失/重复)	多种工具联用
Duplication	基因片段多了一份或多份	读深↑
Deletion	基因片段丢失	读深↓
读深法(RD)	根据测序深度推断CNV	CNVnator, GATK gCNV
分裂读法(SR)	根据软截断reads推断断点	DELLY, Manta
读对法(RP)	根据异常配对信息推断	DELLY, Lumpy
duphold	用深度信息验证SV/CNV真假	过滤假阳性

一、原理（白话版）¶

1.1 什么是CNV/Duplication？¶

正常基因组：
  A-B-C-D-E-F（每个区域一份拷贝）

Duplication（重复）：
  A-B-C-C-D-E-F（C区域多了一份）
  → 这段区域的测序深度会是正常的2倍
  → 基因剂量增加，可能导致过表达

Deletion（缺失）：
  A-B-D-E-F（C区域丢了）
  → 这段区域的测序深度降为0
  → 基因缺失，功能丧失

CNV大小：
  小CNV: 50bp ~ 1kb → 最难检测
  中CNV: 1kb ~ 50kb → 主要检测目标
  大CNV: 50kb ~ 数Mb → 较容易检测

1.2 检测方法对比¶

方法	原理	优势	劣势	代表工具
读深(RD)	区域测序深度异常	大CNV敏感	断点不精确	CNVnator, GATK gCNV
分裂读(SR)	读段部分比对	断点精确	小事件可能漏	DELLY, Manta
读对(RP)	配对插入片段异常	中等CNV好	需配对测序	DELLY, Lumpy
组合方法	多种信号结合	最全面	计算量大	DELLY, Manta
深度学习	AI模型预测	新趋势	需训练数据	Cue

二、CNVnator分析流程¶

# ===== CNVnator：基于读深的CNV检测 =====
# 安装：conda install -c bioconda cnvnator

# Step 1: 提取mapping信息
cnvnator \
  -root sample.root \                   # 输出ROOT文件
  -tree sample_sorted.bam \            # 输入排序BAM文件
  -chrom chr1 chr2 chr3 chr4 chr5 chr6 chr7 chr8 chr9 chr10 \
         chr11 chr12 chr13 chr14 chr15 chr16 chr17 chr18 \
         chr19 chr20 chr21 chr22 chrX chrY  # 染色体列表

# Step 2: 生成直方图（计算读深分布）
cnvnator \
  -root sample.root \
  -his 100 \                            # bin大小=100bp（30x WGS推荐）
  -d /path/to/reference_dir/            # 参考基因组每条染色体的fasta目录

# Step 3: 统计学分析
cnvnator \
  -root sample.root \
  -stat 100                             # 统计bin大小=100的直方图

# Step 4: 分区分割（RD信号分段）
cnvnator \
  -root sample.root \
  -partition 100                        # 分区分析

# Step 5: 检测CNV
cnvnator \
  -root sample.root \
  -call 100 \                           # 调用CNV
  > sample_cnvnator.txt                 # 输出结果

# 输出格式解释：
# 列1: CNV类型(duplication/deletion)
# 列2: 坐标(chr:start-end)
# 列3: 大小(bp)
# 列4: 归一化读深(RD)  → duplication>1.5, deletion<0.5
# 列5-7: e-val等统计值

# ===== 过滤结果 =====
# 筛选高置信度的duplication
awk '$1=="duplication" && $4>1.5 && $5<0.01' sample_cnvnator.txt \
  > sample_duplications_filtered.txt    # 只保留读深>1.5且p<0.01的重复

三、GATK gCNV分析流程¶

# ===== GATK gCNV：基于模型的CNV检测（推荐） =====

# Step 1: 收集读深计数
gatk CollectReadCounts \
  -R reference.fa \                     # 参考基因组
  -I sample_recal.bam \                # 输入BAM
  -L intervals.interval_list \          # 区间列表
  --interval-merging-rule OVERLAPPING_ONLY \  # 区间合并规则
  -O sample.counts.hdf5                 # 输出计数文件(HDF5格式)

# Step 2: 准备区间列表（WGS用等宽bin）
gatk PreprocessIntervals \
  -R reference.fa \
  --bin-length 1000 \                   # bin大小（WGS推荐1000bp）
  --padding 0 \                         # 不加padding
  -O intervals.interval_list            # 输出区间文件

# 批量收集所有样本的计数
for sample in sample1 sample2 sample3; do
    gatk CollectReadCounts \
      -R reference.fa \
      -I ${sample}_recal.bam \
      -L intervals.interval_list \
      -O ${sample}.counts.hdf5 &
done
wait

# Step 3: 去噪建模（Cohort模式，需要多个样本）
gatk DetermineGermlineContigPloidy \
  -L intervals.interval_list \          # 区间
  --input sample1.counts.hdf5 \         # 样本1计数
  --input sample2.counts.hdf5 \         # 样本2计数
  --input sample3.counts.hdf5 \         # 样本3计数
  --contig-ploidy-priors ploidy_priors.tsv \  # 倍性先验
  --output ploidy_model \               # 输出倍性模型
  --output-prefix ploidy                # 输出前缀

# Step 4: 运行gCNV
gatk GermlineCNVCaller \
  --run-mode COHORT \                   # Cohort模式（多样本联合）
  -L intervals.interval_list \
  --input sample1.counts.hdf5 \
  --input sample2.counts.hdf5 \
  --input sample3.counts.hdf5 \
  --contig-ploidy-calls ploidy_model/ploidy-calls \  # 倍性结果
  --output cnv_model \                  # 输出CNV模型
  --output-prefix cnv                   # 输出前缀

# Step 5: 后处理（提取每个样本的CNV调用）
gatk PostprocessGermlineCNVCalls \
  --model-shard-path cnv_model/cnv-model \  # CNV模型
  --calls-shard-path cnv_model/cnv-calls \  # CNV调用
  --allosomal-contig chrX \             # 性染色体
  --allosomal-contig chrY \
  --contig-ploidy-calls ploidy_model/ploidy-calls \
  --sample-index 0 \                    # 样本索引（0=第一个样本）
  --output-genotyped-intervals sample1_intervals.vcf.gz \  # 区间级VCF
  --output-genotyped-segments sample1_segments.vcf.gz       # 片段级VCF

四、多工具联合策略（2025最佳实践）¶

# ===== 多工具联用提高准确性 =====
# 2025年benchmark显示：单一工具的F1值通常<70%
# 多工具取交集可将精确度提升至90%+

# 工具一：DELLY（SR+RP信号）
delly call \
  -g reference.fa \                     # 参考基因组
  -o sample_delly.bcf \                # 输出BCF文件
  sample_sorted.bam                     # 输入BAM

# 提取duplications
bcftools view \
  -i 'FILTER="PASS" && INFO/SVTYPE="DUP"' \  # 只要PASS且类型为DUP
  sample_delly.bcf \
  > sample_delly_dup.vcf

# 工具二：Manta
configManta.py \
  --bam sample_sorted.bam \            # 输入BAM
  --referenceFasta reference.fa \       # 参考基因组
  --runDir manta_output                 # 输出目录

python manta_output/runWorkflow.py \
  -j 8                                  # 8线程

# 工具三：CNVnator（读深信号）
# （见上方CNVnator流程）

# ===== duphold验证（关键过滤步骤） =====
# duphold给每个SV/CNV添加深度注释，帮助区分真假
duphold \
  -t 4 \                               # 线程数
  -v sample_delly.vcf \                # 输入SV/CNV的VCF
  -b sample_sorted.bam \              # BAM文件
  -f reference.fa \                     # 参考基因组
  -o sample_delly_duphold.vcf           # 输出带深度注释的VCF

# duphold添加的关键注释：
# DHFFC: flanking fold-change（侧翼区深度变化）
# DHBFC: breakpoint fold-change（断点区域深度变化）
# 对于真实的duplication: DHBFC > 1.3（区域内深度增加）
# 对于真实的deletion: DHBFC < 0.7（区域内深度降低）

# 用duphold过滤
bcftools view \
  -i 'INFO/SVTYPE="DUP" && FMT/DHBFC > 1.3' \  # Duplication: 区域深度>1.3x
  sample_delly_duphold.vcf \
  > sample_dup_high_conf.vcf            # 高置信度duplication

bcftools view \
  -i 'INFO/SVTYPE="DEL" && FMT/DHBFC < 0.7' \  # Deletion: 区域深度<0.7x
  sample_delly_duphold.vcf \
  > sample_del_high_conf.vcf            # 高置信度deletion

五、Python可视化CNV¶

# ===== Python可视化CNV结果 =====
import pandas as pd  # 导入pandas
import matplotlib.pyplot as plt  # 导入matplotlib
import numpy as np  # 导入numpy

# 读取CNVnator结果
cnv_data = []  # 存储CNV数据
with open("sample_cnvnator.txt") as f:
    for line in f:
        parts = line.strip().split()  # 按空格分割
        cnv_type = parts[0]           # duplication或deletion
        coords = parts[1]            # chr:start-end
        size = float(parts[2])       # 大小
        rd = float(parts[3])         # 归一化读深
        cnv_data.append({
            "type": cnv_type,
            "coords": coords,
            "size": size,
            "rd": rd
        })

df = pd.DataFrame(cnv_data)  # 转为DataFrame

# 可视化：CNV大小分布
fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# 左图：大小分布
for cnv_type, color in [("duplication", "red"), ("deletion", "blue")]:
    subset = df[df["type"] == cnv_type]  # 按类型筛选
    axes[0].hist(np.log10(subset["size"]), bins=50, alpha=0.6,
                 color=color, label=cnv_type)  # 对数尺度直方图
axes[0].set_xlabel("log10(CNV size / bp)")  # x轴标签
axes[0].set_ylabel("Count")  # y轴标签
axes[0].set_title("CNV Size Distribution")  # 标题
axes[0].legend()  # 图例

# 右图：读深分布
axes[1].scatter(
    np.log10(df["size"]),  # x=大小
    df["rd"],              # y=归一化读深
    c=df["type"].map({"duplication": "red", "deletion": "blue"}),  # 颜色
    alpha=0.5, s=10
)
axes[1].axhline(y=1.0, color="gray", linestyle="--", label="Normal RD=1")  # 正常线
axes[1].axhline(y=1.5, color="red", linestyle="--", alpha=0.5, label="DUP threshold")  # DUP阈值
axes[1].axhline(y=0.5, color="blue", linestyle="--", alpha=0.5, label="DEL threshold")  # DEL阈值
axes[1].set_xlabel("log10(CNV size / bp)")
axes[1].set_ylabel("Normalized Read Depth")
axes[1].set_title("Read Depth vs CNV Size")
axes[1].legend()

plt.tight_layout()
plt.savefig("cnv_overview.png", dpi=300)  # 保存图片
plt.show()

六、常见报错与解决¶

报错信息	原因	解决方案
`CNVnator: ROOT error`	ROOT未安装	conda安装cnvnator自带ROOT
`GATK gCNV: too few samples`	Cohort模式需多样本	至少3个样本，推荐>10个
`duphold: no reads`	BAM路径错误	检查BAM文件和索引
`Duplication false positives`	重复区域/低复杂度	用duphold过滤DHBFC<1.3的
`bin size选择`	bin太大漏小CNV	30x WGS: bin=100-1000
`小CNV检测率低`	<1kb CNV难以检测	多工具联用，增加测序深度

七、面试高频问题¶

Q1: 如何提高CNV检测的准确性？¶

A: 多工具联用是关键。2025年benchmark显示，组合GATK gCNV + DELLY + CNVnator，取≥2个工具一致的结果，可以显著提高精确度。另外duphold工具利用深度信息过滤假阳性非常有效：真duplication区域内读深应增加(DHBFC>1.3)。

Q2: 检测duplication和deletion哪个更难？¶

A: Duplication更难检测。2025年临床基准测试显示，所有工具对duplication的敏感性都明显低于deletion。原因：①deletion导致读深为0，信号清晰；②duplication只是读深从1x变成2x，变化较小且容易受GC含量影响；③串联重复的比对困难。

Q3: WGS和WES检测CNV有什么区别？¶

A: WGS覆盖全基因组，能检测基因间区的CNV，分辨率更高。WES只覆盖外显子（~2%基因组），CNV检测受限于目标区域，断点精度差，但成本低。对于临床外显子CNV检测，WES配合专门工具（如XHMM、ExomeDepth）可以满足需求。

八、速查表¶

# ===== CNV/Duplication检测速查 =====

# CNVnator（读深法）
cnvnator -root s.root -tree s.bam
cnvnator -root s.root -his 100 -d ref_dir/
cnvnator -root s.root -stat 100
cnvnator -root s.root -partition 100
cnvnator -root s.root -call 100 > cnv.txt

# GATK gCNV
gatk CollectReadCounts -I s.bam -L intervals -O s.hdf5
gatk DetermineGermlineContigPloidy ...
gatk GermlineCNVCaller --run-mode COHORT ...
gatk PostprocessGermlineCNVCalls ...

# duphold验证（关键！）
duphold -v sv.vcf -b s.bam -f ref.fa -o sv_dh.vcf
# DUP: DHBFC > 1.3 → 真
# DEL: DHBFC < 0.7 → 真

# 多工具联用推荐组合：
# DELLY(SR+RP) + CNVnator(RD) + GATK gCNV(模型)
# 取≥2工具一致 → 高置信度CNV