用broad出品的软件来处理bam文件几次遇到文件头错误

报错如下:ERROR MESSAGE: SAM/BAM file input.marked.bam is malformed: SAM file doesn't have any read groups defined in the header.  The GATK no longer supports SAM files without read groups !

有些人遇到的是bam的染色体顺序不一样,还有可能是染色体的名字不一样,比如>1和>chr1的区别,虽然很傻,但是遇到这样问题的还不少!
还有一些人是遇到基因组没有dict文件,也是用picard处理一下就好。

大部分人是在GATK遇到的,我是在RNA-SeQC遇到的,不过原理都是一样的。
都是因为做alignment的时候并未添加头信息,比如:
bwa samse ref.fa my.sai my.fastq > my.sam
samtools view -bS my.sam > my.bam
samtools sort my.bam my_sorted
java -jar ReordereSam.jar I=/path/my_sorted.bam O=/path/my_reordered.bam R=/path/ref.fa
通过这个代码可以得到排序好的bam,但是接下来用GATK就会报错
java -jar GenomeAnalysisTK.jar -T DepthOfCoverage -R /paht/ref.fa -I /path/aln_reordered.bam
就是因为没有头信息,group相关信息,解决方法有两种:
bwa samse -r @RG\tID:IDa\tSM:SM\tPL:Illumina ref.fa my.sai my.fastq > my.sam
java -jar AddOrReplaceReadGroups I=my.bam O=myGr.bam LB=whatever PL=illumina PU=whatever SM=whatever
一种是比对的时候就加入头信息,这个需要比对工具的支持。
第二种是用picard工具来修改bam,推荐用这个!虽然我其实并不懂这些头文件信息是干嘛的, 但是broad开发的软件就是需要!希望将来去读PHD能系统性的学习一些基础知识!

 

Comments are closed.