使用GATK進(jìn)行SNP和INDEL檢測(cè)
GATK是 Broad 開發(fā)的用于測(cè)序數(shù)據(jù)的變異檢測(cè)軟件镀赌，后續(xù)推廣到動(dòng)植物研究中，是目前最廣泛使用的變異檢測(cè)軟件坪创。GATK有兩種方式變異檢測(cè)承绸，一種是合并所有樣品的gVCF文件，再通過GenotypeGVCF做變異檢測(cè);另一種是每個(gè)樣品的GVCF文件先生成genomeDB文件再進(jìn)行變異檢測(cè)蜕琴。兩者的結(jié)果沒有差異，根據(jù)自己的樣品數(shù)量合理選擇變異檢測(cè)的方法宵溅。
軟件
GATK
準(zhǔn)備文件
基因組文件：genome.fasta
比對(duì)結(jié)果文件：S1.sort.rmdup.bam
參考腳本
#首先創(chuàng)建一個(gè)存放臨時(shí)文件的目錄
mkdir tmp
#每個(gè)樣品進(jìn)行HaplotypeCaller變異檢測(cè)
gatk --java-options "-Xmx10g -Djava.io.tmpdir=./tmp" \ #設(shè)置java參數(shù)
HaplotypeCaller \
-R ./genome.fasta \ #參考基因組路徑
-I ./S1.sort.rmdup.bam \ #bam文件路徑
-ERC GVCF \ # 輸出GVCF文件
-O S1.g.vcf \ # 指定輸出結(jié)果
1>S1.HC.log 2>&1 #寫到日志文件
#合并每個(gè)樣品的GVCF文件
ls ./*.g.vcf > gvcf.list \ #合并文件的列表
gatk --java-options "-Xmx10g \#限制使用內(nèi)存
-Djava.io.tmpdir=./tmp" \
CombineGVCFs \
-R ./genome.fasta \ 參考基因組
-V gvcf.list \ #GVCF文件列表
-O all.merge.g.vcf #輸出的vcf文件
SNP
SNP(單核苷酸多態(tài)性)主要是指在基因組水平上由單個(gè)核苷酸的變異所引起的DNA序列多態(tài)性凌简，包括單個(gè)堿基的轉(zhuǎn)換、顛換等恃逻，為了定位目標(biāo)性狀雏搂，每組性狀相關(guān)樣本一起call 群體SNP。
#群體變異檢測(cè)
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' \ #符合其中任意一個(gè)條件的都過濾掉
-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" \ #符合其中任意一個(gè)條件的都過濾掉
--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的檢測(cè)就完成了，接著就可以繼續(xù)做后續(xù)分析了矛市。
變異檢測(cè)結(jié)果可視化
對(duì)于比對(duì)的bam格式文件芙沥，可以使用IGV軟件對(duì)變異檢測(cè)結(jié)果進(jìn)行可視化瀏覽。使用方法參考使用說明文檔（IGV使用文檔.pdf）浊吏。
