【直播】我的基因组24：用GATK对SAM格式的文件进行重排
GATK教程我以前写过（:）
这个步骤分成2个小步骤：
首先用RealignerTargetCreator找到需要重新比对的区域，输出文件intervals，然后用输出的 tmp.intervals 做输入文件来进行重新比对，也就是用IndelRealigner在这些区域内进行重新比对
这个是批量运行的代码，就是用shell脚本写一个循环即可（:）。
sam/bam的染色体顺序不一样，还有可能是染色体的名字不一样，比如>1和>chr1的区别。
下面是对我的基因组bam文件用GATK进行重排，命令分成两步，请注意下面用到的参考基因组文件和软件，需要保证都是按照前面的教程安装，用的的vcf文件下载地址见文末
nohup java -Xmx60g -jar ~/biosoft/GATK/GenomeAnalysisTK.jar -R ~/reference/genome/human_g1k_v37/human_g1k_v37.fasta -T RealignerTargetCreator -I P_jmzeng.filter.rmdup.bam -o P_jmzeng.filter.rmdup.intervals -known ~/annotation/GATK/1000G_phase1.indels.b37.vcf.gz 1> P_jmzeng.filter.rmdup.RealignerTargetCreator.log
## 请保证你的bam文件比对的时候选择的参考基因组跟你在使用GATK的时候自己指定的参考基因组的一致性。
nohup java -Xmx60g -jar ~/biosoft/GATK/GenomeAnalysisTK.jar -R ~/reference/genome/human_g1k_v37/human_g1k_v37.fasta -T IndelRealigner -I P_jmzeng.filter.rmdup.bam -targetIntervals P_jmzeng.filter.rmdup.intervals -o P_jmzeng.filter.rmdup.realgn.bam 1>P_jmzeng.filter.rmdup.IndelRealigner.log
输出文件如下：
（nohup的命令在之前讲过，是后台运行的命令，其中1表示运行的日志输入log文件，如果用2的话，表示的是指输出报错的日志）
在探索这个软件的用法的时候可能会报错，根据报错日志搜索解决即可;
第一个错误：
MESSAGE: Fasta dict file /home/jianmingzeng/reference/genome/human_g1k_v37/human_g1k_v37.dict for reference /home/jianmingzeng/reference/genome/human_g1k_v37/human_g1k_v37.fasta does not exist. Please see http://gatkforums.broadinstitute.org/discussion/1601/how-can-i-prepare-a-fasta-file-to-use-as-reference for help creating it.
解决方案，在参考基因组文件所在的目录运行下面代码：
java -jar ~/biosoft/picardtools/picard-tools-1.119/CreateSequenceDictionary.jar R=human_g1k_v37.fasta O=human_g1k_v37.dict
第二个错误：
MESSAGE: An index is required, K/1000G_phase1.indels.b37.vcf.gz
因为我在broad 网站下载的文件需要gunzip解压。
文件下载地址是：
下载代码是：
## ftp://ftp.broadinstitute.org/bundle/b37/
mkdir -p ~/annotation/GATK
cd ~/annotation/variation/GATK
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/1000G_phase1.snps.high_confidence.b37.vcf.gz
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/dbsnp_138.b37.vcf.gz
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/human_g1k_v37.fasta.gz
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/NA12878.HiSeq.WGS.bwa.cleaned.raw.subset.b37.sites.vcf.gz
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/Mills_and_1000G_gold_standard.indels.b37.vcf.gz
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/hapmap_3.3.b37.vcf.gz
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/1000G_phase1.indels.b37.vcf.gz
wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b37/1000G_phase1.indels.b37.vcf.idx.gz
gunzip 1000G_phase1.indels.b37.vcf.idx.gz
gunzip 1000G_phase1.indels.b37.vcf.gz
ref:
文：Jimmy、阿尔的太阳
图文编辑：吃瓜群众