GATK4 —— 获取短变异 （call SNP+indel）
GATK是一款用于基因组数据分析的软件，其强大的处理引擎和高性能计算功能使其能够承担任何规模的项目。
GATK的功能非常强大，这里不详细介绍，大家可以根据自己的要求，从首页进入对应的模块，说明书还是很详细的。
官网：
https://gatk.broadinstitute.org/hc/en-us
GATK 4.2.1.0 命令说明：
https://gatk.broadinstitute.org/hc/en-us/sections/4404575821339-4-2-1-0?page=10#articles
1. 下载安装
1.1下载
下载链接：
https://github.com/broadinstitute/gatk/releases
1.2 安装
解压：
$ unzip gatk-4.2.4.1.zip
查看帮助：
$ cd /your/path/gatk-4.2.4.1
$ ./gatk --help
查看工具列表：
$ ./gatk --list
查看指定工具帮助：
$ ./gatk ToolName --help
2.获取短变异（SNP+indel）
说明书：
https://gatk.broadinstitute.org/hc/en-us/articles/360035535932-Germline-short-variant-discovery-SNPs-Indels-
主要流程如下：
2.1 准备数据
-
这里使用在之前已经比对好并排序去重后的bam文件（.bam），详情参见以下链接：
测序质量分析 —— FastQC - 简书 (jianshu.com)
gff格式与gtf格式转换——NBISweden / AGAT - 简书 (jianshu.com)
序列比对 —— Hisat2 - 简书 (jianshu.com)
sam文件转换为bam文件——SAMtools - 简书 (jianshu.com)
SAMtools——bam文件排序 - 简书 (jianshu.com)
SAMtools——bam文件去重 - 简书 (jianshu.com) -
bam文件索引（.bai）
$ samtools index LPF1_R1_MP.sort.dup.bam
- 参考基因组（.fa）
- 参考基因组索引（.fai）
$ samtools faidx Mparg_v2.0.fa
$ less Mparg_v2.0.fa.fai
Mpar_chr1 41449654 11 100 101
Mpar_chr2 50537509 41864173 100 101
Mpar_chr3 46570116 92907069 100 101
Mpar_chr4 46692298 139942898 100 101
Mpar_chr5 50691827 187102130 100 101
Mpar_chr6 65498417 238300887 100 101
Mpar_chr7 34026340 304454300 100 101
Mpar_chr8 48943667 338820915 100 101
五列分别对应，序列名、序列长度、第一个碱基的偏移量、行的碱基数、行宽。
- 参考基因组索引（.dict）
$ gatk CreateSequenceDictionary -R Mparg_v2.0.fa -O Mparg_v2.0.dict
2.2 HaplotypeCaller
HaplotypeCaller是GATK4中的核心工具，可以利用单倍型区域的局部从头组装同时调用种系snv和小indel。详细原理参见说明书。
https://gatk.broadinstitute.org/hc/en-us/articles/360050814612-HaplotypeCaller#--sample-name
$ ./gatk HaplotypeCaller --help
Required Arguments:
--input,-I <GATKPath> BAM/SAM/CRAM file containing reads This argument must be specified at least once.
Required.
--output,-O <GATKPath> File to which variants should be written Required.
--reference,-R <GATKPath> Reference sequence file Required.
示例：
$ gatk --java-options "-Xmx100g -XX:ParallelGCThreads=4" HaplotypeCaller \
-I ~/LPF1_R1_MP.sort.dup.bam \
-R ~/Mparg_v2.0.fa \
-ERC GVCF \
-O ~/vcf/LPF1_R1_MP.g.vcf.gz
报错：A USER ERROR has occurred: Argument emit-ref-confidence has a bad value: Can only be used in single sample mode currently. Use the --sample-name argument to run on a single sample out of a multi-sample BAM file.
解决办法：
picard——修改BAM文件的Read Group - 简书 (jianshu.com)
更改命令：
$ samtools index LPF1_R1_MP.rg.sort.bam
$ gatk --java-options "-Xmx100g -XX:ParallelGCThreads=4" HaplotypeCaller \
-I ~/LPF1_R1_MP.rg.sort.bam \
-R ~/Mparg_v2.0.fa \
-ERC GVCF \
-O ~/LPF1_R1_MP.g.vcf.gz
gatk好像最多只支持4个线程，如果我我理解有错误的话，还请大佬指正。
这里注意HaplotypeCaller只能处理单样本文件，当有多样本时，官方建议使用HaplotypeCaller对单bam文件分别进行变异检测，生成GVCF文件，GVCF会记录每一个位点到情况，包括有无突变，VCF只记录突变位点情况，之后在下一步对GVCF文件进行合并。
VCF文件和GVCF文件的格式的详细说明参见下链接：
VCF和GVCF格式说明 - 青萍，你好 - 博客园 (cnblogs.com)
3.GVCF文件合并
官方说明书：https://gatk.broadinstitute.org/hc/en-us/articles/4404604801947-CombineGVCFs
gatk --java-options "-Xmx100g -XX:ParallelGCThreads=4" CombineGVCFs \
-R ~/ref/Mparg_v2.0.fa \
-V LPF1_R1_MP.g.vcf.gz \
-V LPF1_R2_MP.g.vcf.gz \
-V LPF1_R3_MP.g.vcf.gz \
-O LPF1_MP.g.vcf.gz
