使用GATK进行SNP和INDEL检测
GATK是 Broad 开发的用于测序数据的变异检测软件，后续推广到动植物研究中，是目前最广泛使用的变异检测软件。GATK有两种方式变异检测，一种是合并所有样品的gVCF文件，再通过GenotypeGVCF做变异检测;另一种是每个样品的GVCF文件先生成genomeDB文件再进行变异检测。两者的结果没有差异，根据自己的样品数量合理选择变异检测的方法。

软件

GATK

准备文件

基因组文件：genome.fasta
比对结果文件：S1.sort.rmdup.bam

参考脚本

#首先创建一个存放临时文件的目录
mkdir tmp

#每个样品进行HaplotypeCaller变异检测
gatk --java-options "-Xmx10g -Djava.io.tmpdir=./tmp" \  #设置java参数
HaplotypeCaller \  
-R ./genome.fasta \  #参考基因组路径
-I ./S1.sort.rmdup.bam  \ #bam文件路径
-ERC GVCF \ # 输出GVCF文件
-O S1.g.vcf  \ # 指定输出结果
1>S1.HC.log   2>&1  #写到日志文件


#合并每个样品的GVCF文件
ls ./*.g.vcf > gvcf.list  \ #合并文件的列表
gatk  --java-options "-Xmx10g \#限制使用内存
-Djava.io.tmpdir=./tmp"    \
CombineGVCFs \
-R ./genome.fasta \ 参考基因组
-V gvcf.list  \ #GVCF文件列表
-O all.merge.g.vcf  #输出的vcf文件

SNP

SNP(单核苷酸多态性)主要是指在基因组水平上由单个核苷酸的变异所引起的DNA序列多态性，包括单个碱基的转换、颠换等，为了定位目标性状，每组性状相关样本一起call 群体SNP。

#群体变异检测
gatk  --java-options "-Xmx10g -Djava.io.tmpdir=./tmp"   \
GenotypeGVCFs -R ./genome.fasta \
-v all.merge.g.vcf \ #群体的GVCF文件
-O all.merge_raw.vcf  #输出文件

#提取SNP
gatk  --java-options "-Xmx10g -Djava.io.tmpdir=./tmp"  \
SelectVariants  \
-R ./genome.fasta \
-V all.merge_raw.vcf \#合并后的VCF文件
--select-type SNP \#指定提取的类型SNP
-O all.raw.snp.vcf 

#过滤SNP
gatk  --java-options "-Xmx10g -Djava.io.tmpdir=./tmp"  \
VariantFiltration \
-R ./genome.fasta \
-V all.raw.snp.vcf \
--filter-expression "QD < 2.0 || MQ < 40.0 || FS > 60.0 || SOR > 3.0 || MQRankSum < -12.5 || ReadPosRankSum < -8.0" --filter-name 'SNP_filter' \ #符合其中任意一个条件的都过滤掉
-O all.filter.snp.vcf

#提取过滤好的SNP
gatk  --java-options "-Xmx4g -Djava.io.tmpdir=./tmp"  \
SelectVariants  -R ./genome.fasta \
-V all.filter.snp.vcf --exclude-filtered \
 -O all.filtered.snp.vcf

INDEL

InDel是小型的Insertion和Deletion的总称，利用SelectVariants命令来提取 InDel 突变信息。

#提取InDel文件
gatk  --java-options "-Xmx10g -Djava.io.tmpdir=./tmp"  \
SelectVariants  -R ./genome.fasta \#参考基因组
-V all.merge_raw.vcf \#合并后的VCF文件
--select-type INDEL \#筛选提取类型
-O all.raw.indel.vcf  #输出文件

#过滤InDel
gatk  --java-options "-Xmx10g -Djava.io.tmpdir=./tmp"  \
VariantFiltration -R ./genome.fasta  \
-V all.raw.indel.vcf  \ #提取后的InDel文件
--filter-expression "QD < 2.0 || FS > 200.0 || SOR > 10.0 || MQRankSum < -12.5 ||  ReadPosRankSum < -8.0"  \ #符合其中任意一个条件的都过滤掉
--filter-name 'INDEL_filter'  \ 
-O all.filter.indel.vcf

#提取过滤后的InDel
gatk  --java-options "-Xmx10g -Djava.io.tmpdir=./tmp"  \
SelectVariants  -R ./genome.fasta   \
-V all.filter.indel.vcf   \#过滤后的InDel文件
--exclude-filtered   \
-O all.filtered.indel.vcf

SNP和InDel的检测就完成了，接着就可以继续做后续分析了。

变异检测结果可视化

对于比对的bam格式文件，可以使用IGV软件对变异检测结果进行可视化浏览。使用方法参考使用说明文档（IGV使用文档.pdf）。

欢迎关注Bioinfor 生信云微信公众号！