一起读官方文档 GATK2.1版本 篇
参考文章:
GATK介绍
GATK做什么的？
它主要用于从sequencing 数据中进行variant calling，包括SNP、INDEL。比如现在风行的exome sequencing找variant，一般通过BWA+GATK的pipeline进行数据分析。
BWA流程上一篇文章已经讲完了，这一篇主要讲一下GATK2版本的使用。
GATK参数 -- GATK2.1版本
这里只讲述了在BWA + GATK流程中的GATK2.1版本应用
STEP01: 局部重新对齐
分为两个小步骤1.RealignerTargetCreator 和2. IndelRealigner
- GATK
RealignerTargetCreator
RealignerTargetCreator 用于找到那些可能需要进行局部重新对齐的区域。
-Xmx256g分配了最多 256GB 的内存给 Java 虚拟机
-T RealignerTargetCreator使用RealignerTargetCreator工具。
-R指定参考基因组
-I指定输入的 BAM 文件
-o输出文件，其中包含了需要进行局部重新对齐的目标区域。这通常是一个.intervals文件
# 实例
java -Xmx256g -jar GenomeAnalysisTK-2.1/GenomeAnalysisTK.jar -T RealignerTargetCreator -R $ref -I bam/$file.bam -o bam/$file.realn.intervals
- GATK
IndelRealigner
IndelRealigner根据上面的RealignerTargetCreator对这些局部区域进行重新对齐
-T IndelRealigner使用IndelRealigner工具。
-targetIntervals指定由上一步RealignerTargetCreator生成的.intervals文件路径，主要表示哪些区域需要进行重新对齐
java -Xmx256g -jar GenomeAnalysisTK-2.1/GenomeAnalysisTK.jar -T IndelRealigner -R $ref -targetIntervals bam/$file.realn.intervals -I bam/$file.bam -o bam/$file.realn.bam
注意RealignerTargetCreator和IndelRealigner是一块使用的，这在GATK4版本中已经可以使用HaplotypeCaller
STEP02: picard进行MarkDuplicates
MarkDuplicates 工具对GATK的IndelRealigner后的$file.realn.bam进行 标记和处理重复的序列读取。
-Xmx4g: 为 Java 虚拟机分配最多 4GB 的内存
-XX:ConcGCThreads=2: 设置并发垃圾回收器线程的数量为 2
-XX:ParallelGCThreads=12: 设置用于并行垃圾回收的线程数量为 12
MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=8000: 设置用于存储 read ends 信息的最大文件句柄数
TMP_DIR=TMP: 设置临时目录。
INPUT=bam/$file.realn.bam: 输入 BAM 文件
OUTPUT=bam/$file.realn.marked.bam: 包含标记了重复读取信息的BAM文件
METRICS_FILE=bam/$file.realn.marked.txt: 输出的文本文件，其中包含有关重复标记的各种度量和统计信息。
VALIDATION_STRINGENCY=LENIENT: 验证严格性的设置。LENIENT表示在遇到不完全符合标准的行时发出警告，而不是终止。
java -Xmx4g -XX:ConcGCThreads=2 -XX:ParallelGCThreads=12 -jar picard.jar MarkDuplicates MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=8000 \
TMP_DIR=TMP \
INPUT=bam/$file.realn.bam \
OUTPUT=bam/$file.realn.marked.bam \
METRICS_FILE=bam/$file.realn.marked.txt \
VALIDATION_STRINGENCY=LENIENT
# samltools对bam文件建立索引后续使用
samtools index bam/$file.realn.marked.bam
STEP03: GATK UnifiedGenotyper进行基因变异检测
UnifiedGenotyper同时检测单核苷酸变异（SNPs）和短插入/缺失（indels）。最新版本GATK4已经使用更先进的工具HaplotypeCaller，之后会出一期专门的文章阅读GATK4教程。下面是对该GATK2.1中的UnifiedGenotyper参数讲解。
-Xmx256g: 为 Java 虚拟机分配最多 256GB 的内存。
-T UnifiedGenotyper: 指定使用UnifiedGenotyper工具。
--num_threads 30: 指定使用 30 个线程进行计算。
-R $ref: 指定参考基因组的路径（这里由 shell 变量$ref保存）。
-I bam/$file.realn.marked.bam -I bam/WT.realn.marked.bam: 指定输入的 BAM 文件。这里有两个输入文件：一个是样本（由$file表示），另一个是野生型（WT）。
-o VCF/$file-WT.GATK2.1.vcf: 指定输出的 VCF 文件的路径和名称。
-stand_call_conf 30.0: 设置变异检测的置信度阈值为 30。
-glm BOTH: 指定同时检测 SNPs 和 indels。
-rf BadCigar: 使用BadCigar读取过滤器来过滤具有不合规或损坏的 CIGAR 字符串的读取。
java -Xmx256g -jar GenomeAnalysisTK-2.1/GenomeAnalysisTK.jar -T UnifiedGenotyper \
--num_threads 30 \
-R $ref \
-I bam/$file.realn.marked.bam -I bam/WT.realn.marked.bam \
-o VCF/$file\-WT.GATK2.1.vcf \
-stand_call_conf 30.0 \
-glm BOTH \
-rf BadCigar
使用bcftools查看一下生成的vcf文件 : bcftools view W356-WT.GATK2.1.vcf.gz | less -SN
