GATK 概述与流程基础
我最早接触 GATK 是在 2015 年，那时还是 GATK3，一个巨大的 Java 单体 jar 包，跑个全基因组得先算好内存，一不留神就 OOM 了。到了 GATK4，整个架构重写了，用起来舒服很多。今天咱们就从最基础的开始，把 GATK4 长什么样、Best Practices 工作流怎么走、每一步干什么，捋一遍。
GATK4 整体架构
GATK4 的核心设计理念就是“模块化”。不是你想象的一个大工具包，而是一堆独立的小工具，每个工具干一件事。你可以像搭积木一样把它们串起来，也方便并行。
| 层级 | 组件 | 说明 |
|---|---|---|
| 工具层 | 单个 Java/Python 工具 | 比如 BwaMemIndexImageCreator、HaplotypeCaller、VariantFiltration |
| 引擎层 | Spark / plain Java 运行时 | 负责调度、内存管理、并行 |
| 编排层 | WDL / CWL / Nextflow | 把工具串成工作流 |
| 数据层 | HTSJDK / ADAM | 抽象读/写 BAM、CRAM、VCF |
说白了，GATK4 不再是“一个命令走到底”，而是给你一堆小工具，你自己或者用 WDL 脚本把它们拼起来。我个人习惯直接用 Cromwell 跑 WDL 工作流，这样集群调度、失败重试都自动搞定。
架构亮点：GATK4 原生支持 Spark 分布式引擎，对于超大 WGS 数据，可以边读 BAM 边做变异检测，不用把所有数据塞进一个 JVM。
我曾经在 30X 全基因组上测过一个 300 人的队列，用 GATK3 的 CombineGVCFs 跑了一周还没完。换到 GATK4 的 GenomicsDBImport，两天搞定。嗯，这就是架构升级带来的红利。
Best Practices 工作流概述
Broad Institute 给的两套标准流程大家应该都听过：DNA-seq 和 RNA-seq。其实 RNA-seq 的核心是从 DNA-seq 改的，主要区别在预处理阶段要加上 spliced alignment 和一些特殊过滤。
下面这个表格对比了一下两个流程的主要步骤，你想想看，差异还是蛮大的。
| 阶段 | DNA-seq | RNA-seq |
|---|---|---|
| 原始数据对齐 | BWA-MEM | STAR / HISAT2（剪接感知） |
| 预处理 | MarkDuplicates + BaseRecalibrator | MarkDuplicates + SplitNCigarReads + BaseRecalibrator |
| 变异检测 | HaplotypeCaller (GVCF mode) | HaplotypeCaller (dont-use-soft-clipped-bases) |
| 联合分析 | GenomicsDBImport + GenotypeGVCFs | 一般不做联合，按样本单独检测 |
| 过滤 | VariantFiltration + VQSR | VariantFiltration + 额外 RNA 过滤（如 lowQual） |
我记得第一次跑 RNA-seq 流程时，没加 SplitNCigarReads，结果 HaplotypeCaller 报告一堆假阳性断点。坑了一整周才发现文档里写着要切剪接 CIGAR。所以你看，流程细节真不能跳。
流程各步骤功能简介
每一个步骤都有自己的使命，我按 DNA-seq 标准流程给大家拆解一下。
- 原始数据预处理（FASTQ → BAM）：用 BWA-MEM 把 reads 贴到参考基因组上，得到 SAM/BAM。这一步决定了后续变异检测的底子。我之前调 BWA 参数时踩过一个坑：-M 标记拆分 reads 时，如果不用 -M，有些 picard 工具会报错。
- MarkDuplicates：标记 PCR 或光学重复，避免过度放大测序噪声。注意 GATK4 里 MarkDuplicates 默认不会真正删除 duplicate，只是打标记。
- BaseRecalibrator + ApplyBQSR：基础质量分数重校准。说白了，就是用已知的 SNP 位点修正测序仪的系统偏差。我建议先跑一次 BaseRecalibrator 看看报告，如果校准前后的均方残差变化很小，说明原始质量已经不错，可以跳过。
- HaplotypeCaller：核心变异检测工具。它先局部组装 haplotype，再概率性地推断基因型。最好跑 GVCF 模式，这样出来的 gVCF 可以跟其他样本一起联合分析。
- GenomicsDBImport：把多个 gVCF 导入一个数据库，省掉了传统 CombineGVCFs 的 I/O 瓶颈。注意这个工具对存储要求高，我建议用 SSD 或 Lustre 并行文件系统。
- GenotypeGVCFs：对数据库中的各样本联合做基因分型。这一步比较快，但内存消耗随样本数线性增长。
- VariantFiltration / VQSR：过滤假阳性。VQSR 基于高斯混合模型，需要大量已知变异集做训练；数据量少时直接用硬过滤条件更可靠。
下面这张图把 GATK4 DNA-seq 工作流的主要节点和流通路径画出来了。你可以对照看一下各个环节怎么衔接。
小贴士：如果你手头有 GPU 节点，可以试试 GATK4 的 HaplotypeCallerKickstart 参数，能利用 cuDNN 加速局部组装，我测过大概能快 2-3 倍。
避坑指南：我曾经在做全外显子组联合分析时，忘记设置 --merge-input-intervals，导致 GenomicsDBImport 把不同染色体的间隔合并错了，最后出来的 VCF 里很多位点被硬截断。建议每次跑之前先用 --validate 检查一下间隔列表。
为什么我要花这么多篇幅讲整体架构和工作流？因为后续的所有优化、并行策略，都是建立在你对这条流水线的理解之上的。如果你连 MarkDuplicates 和 BaseRecalibrator 谁先谁后都搞混，那后面的 Spark 分片、shuffle 优化根本无从谈起。