如何加速训练 gcnv 模型
在训练 gcnv model 时,特别耗时,今天尝试将 interval list 拆分成 200份,
先把 1 kb 的全基因组区间做质量评估与过滤,然后把筛选后的优质区间按 200 份切片,最后对每一份并行跑 GermlineCNVCaller 建立 gCNV 模型。这样既能显著缩短单次运行时间,也能把内存消耗控制在较小范围。下面分步拆解每一句指令,告诉你为什么要这么做、关键参数是什么意思,以及它们到底改善了什么。
1️⃣ 给 1 kb 区间打“质量标签”
gatk AnnotateIntervals \
-L hg38_1kb.interval_list \
-R Homo_sapiens_assembly38.fasta \
-imr OVERLAPPING_ONLY \
-O hg38_1kb_intervals.annotated.tsv
-
AnnotateIntervals 读取 1 kb 的 bin(-L),为每个 bin 计算
- GC 含量
- mappability(比对唯一性)
- segmental duplication(段重复密度)
- -imr OVERLAPPING_ONLY 代表 “只保留原始 bin 的边界”,不合并相邻 bin——这样后续每个 bin 的注释都是一对一的,对过滤和建模精度更高。
2️⃣ 结合样本计数过滤掉“问题 bin”
inputs=$(ls /data2/pipeline/temp_data/counts/*txt | while read id; do echo " -I $id "; done)
gatk FilterIntervals $inputs \
-L hg38_1kb.interval_list \
--annotated-intervals hg38_1kb_intervals.annotated.tsv \
-imr OVERLAPPING_ONLY \
--minimum-gc-content 0.1 --maximum-gc-content 0.9 \
--minimum-mappability 0.9 --maximum-mappability 1.0 \
--minimum-segmental-duplication-content 0.0 \
--maximum-segmental-duplication-content 0.5 \
--low-count-filter-count-threshold 5 \
--low-count-filter-percentage-of-samples 90.0 \
--extreme-count-filter-minimum-percentile 1.0 \
--extreme-count-filter-maximum-percentile 99.0 \
--extreme-count-filter-percentage-of-samples 90.0 \
-O hg38_1kb_intervals_gc_filtered.interval_list
FilterIntervals 将上一步的注释 + 多个样本的 read-count TSV 一起用来打分、筛掉噪声 bin。
-
设置的阈值逻辑:
GC / mappability / 段重复:去掉极端 GC(<10% 或 >90%)、不唯一比对(<0.9) 或高度重复(>0.5) 的区间——这些 bin 通常读数离散、难以建模。
低计数过滤:≥90 % 的样本在该 bin 的 read-count ≤ 5 就剔除,避免测序深度太低的部位。
极端计数过滤:如果 ≥90 % 样本的计数落在全局最上/最下 1 %–99 % 之外,视为离群也剔除。
-
-imr OVERLAPPING_ONLY 同样保持 bin 边界一致。
结果是一个 质量不错的 1 kb interval list,后面所有训练步骤都以它为准。
3️⃣ 把高质量区间“一分为二百”
gatk IntervalListTools \
-M INTERVAL_COUNT \
--SCATTER_COUNT 200 \
-I hg38_1kb_intervals_gc_filtered.interval_list \
-O split_intervals_1kb/
IntervalListTools 带 --SCATTER_COUNT 200 会把输入 interval list 等分为 ~200 片(每片都是一个新的 .interval_list),常用来并行化 CNV 或重测序流程。
-
-M INTERVAL_COUNT 额外输出每个片段实际包含多少 bin,方便你看看均匀度。
为什么 200?官方经验是 GermlineCNVCaller 1 次处理 5k–12.5k bin(1 kb bin 即 5–12 Mb)时,内存 16–32 GB、用时数小时。拆成 200 份后,你就能把 一个需要十几小时的大任务变成 200 个可并行的“小切片”,在集群或多核机器上速度线性提速。
4️⃣ 逐片训练 gCNV 模型
small_interval=$1 # 外层 for 循环把每个切片路径传进来
output_prefix=$(basename $(dirname $small_interval))
mkdir -p /data/pipeline/ycg_gcnv_model
gatk GermlineCNVCaller \
--run-mode COHORT \
--input sample1.counts.hdf5 \
sample2.counts.hdf5 \
sample3.counts.hdf5 \
... \
--intervals $small_interval \
--annotated-intervals hg38_1kb_intervals.annotated.tsv \
--output /data/pipeline/ycg_gcnv_model/${output_prefix} \
[其他常用参数......]
- –run-mode COHORT 让 GermlineCNVCaller 同时学习 (i) 样本共享的 read-depth 背景噪声模型 和 (ii) 每个样本的 CNV 隐状态;官方推荐“把样本数控制在 ~200 以内、bin 数控制在 <12 k” ——而脚本拆片后刚好满足该条件。
- 每个切片对应输出一个像 temp_0123_of_200-model/ 的模型目录(稍后用 tar 打包)。这样你最终会得到 200 个模型 shard,完全满足 GATK-SV 流程期待的“多 tar.gz 模型文件”输入格式。
- upload 行把结果同步到对象存储(如 TOS / S3),方便后续流程拉取。
✅ 脚本带来的核心收益
步骤 | 作用 | 对性能/结果的影响 |
---|---|---|
Annotate + Filter | 去掉 GC 极端、低唯一性、重复区、计数离群的 bin | 减少噪声,提高模型质量 |
Scatter 200 份 | 把 ~3 M1 kb bin 压缩成每份几万 bin | 单次运行更快、内存更低,可并行 |
多实例 GermlineCNVCaller | 每片独立建模,再后期合并 | 整体训练瓶颈→多核/集群并行,总 wall-time 缩短近 N 倍 |
注意
- 切片训练后,要用 PostprocessGermlineCNVCalls(或 GATK-SV 管线脚本里的汇总步骤)把 200 个 shard 的 calls、model 合并/格式化,才能得到最终的 per-sample CNV calls。
- 确保所有切片使用同一批样本 & 参数,否则模型无法无缝拼接。
- 如果你还没跑 DetermineGermlineContigPloidy,记得先做,因为 gCNV 训练需要 ploidy calls 作为 priors。
这样一来,原本“跑一次全基因组 gCNV 模型要十几个小时、几十 GB 内存”的痛点就被拆解成了可以轻松排队或并发执行的小任务,既省时间也省机器。
本文由mdnice多平台发布