GATK 4.1.7.0 call SNP

GATK4.0之后和之前的版本存在较大差异，因此需要使用特定版本来call变异时需要有较大的改动，call变异的过程更趋向于流程化。

首先是软件安装

在github中下载gatk 4.1.7

wget https://github.com/broadinstitute/gatk/releases/download/4.1.7.0/gatk-4.1.7.0.zip
unzip gatk-4.1.7.0.zip

将zip包解压之后修改环境变量并保存。
注：这里使用的samtools版本为1.6。

GATK使用说明

## 帮助信息
gatk --help

## 列出所有的工具
gatk --list

## 工具的说明，比如以VariantAnnotator 为例
gatk HaplotypeCaller --help

GATK call变异简单流程

1. 数据准备

1. ref.fa
2. reads1.fq.gz
3. reads2.fq.gz

建立索引

bwa index ref.fa
samtools  faidx ref.fa
gatk CreateSequenceDictionary -R ref.fa -O ref.dict

2. BWA比对

#使用bwa mem比对算法
bwa mem -t 20 -M -P -R '@RG\tID:20ZS0002_clean\tSM:20ZS0002_clean\tLB:20ZS0002_clean\tPL:Illumina' reference.fa 20ZS0002_clean.out.1.fq.gz 20ZS0002_clean.out.2.fq.gz >20ZS0002_clean.sam
# -t为线程数，-R指定引号内的flag信息。ID一定要每个样本对号入座，SM号、LB号可与ID号一样，PL为测序平台，一般为Illumina

3. 用samtools将sam转换成bam

samtools view -bS 20ZS0002_clean.sam -o 20ZS0002_clean.bam
#sam转bam

samtools sort 20ZS0002_clean.bam -o 20ZS0002_clean.sort.bam
#对bam进行sort排序

4. 标记PCR重复

java -jar picard.jar MarkDuplicates REMOVE_DUPLICATES=false  \
MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=8000 \
INPUT=20ZS0002_clean.sort.bam OUTPUT=20ZS0002_clean.marked.bam \
METRICS_FILE=20ZS0002_clean.marked.bam.metrics

5. 将标记好的bam排序

samtools index 20ZS0002_clean.marked.bam

6. GATK call 变异

gatk --java-options -Xmx8G HaplotypeCaller -R ref.fa -I 20ZS0002_clean.marked.bam \
-O 20ZS0002_clean.g.vcf.gz -ERC GVCF --max-mnp-distance 0
#生成中间文件gvcf文件（gvcf格式包括所有的变异类型，包括SNP和InDel，需要进一步过滤）也就是说现在每个样本一个gvcf文件

7. 将gvcf文件进行整合

gatk --java-options -Xmx256G CombineGVCFs -R ref.fa \
--variant 20ZS0002_clean.g.vcf.gz \
--variant 20ZS0003_clean.g.vcf.gz \
--variant 20ZS0005_clean.g.vcf.gz \
-O All_samples.g.vcf.gz

8. 提取基因型

gatk --java-options -Xmx256G GenotypeGVCFs \
-R ref.fa -V All_samples.g.vcf.gz -O genotype.vcf.gz 

9. 提取SNP

gatk --java-options -Xmx256G SelectVariants \
-R reference.fa -O SNPs.vcf --variant genotype.vcf.gz \
--select-type-to-include SNP
# 同理，也可以提取indels

10. 对SNPs进行过滤

按照官方默认参数即可

gatk VariantFiltration \
    -V SNPs.vcf \
    -filter "QD < 2.0" --filter-name "QD2" \
    -filter "QUAL < 30.0" --filter-name "QUAL30" \
    -filter "SOR > 3.0" --filter-name "SOR3" \
    -filter "FS > 60.0" --filter-name "FS60" \
    -filter "MQ < 40.0" --filter-name "MQ40" \
    -filter "MQRankSum < -12.5" --filter-name "MQRankSum-12.5" \
    -filter "ReadPosRankSum < -8.0" --filter-name "ReadPosRankSum-8" \
    -O SNPs_filtered.vcf