查看原文
其他

Sam/Bam文件格式详解

生信阿拉丁 生信阿拉丁 2022-05-16

随着生物信息数据的爆发式增长,存储生物信息的文件格式也多样化起来,不同的文件格式往往有不同的目的:
为了软件之间的兼容性;
人类可读的用于数据加工、解析和处理的格式,比如.tsv,.csv等;
为了提高计算机计算效率的数据格式,一般可读性不强的二进制文件,比如今天要讲解到的bam文件。

介绍

bam文件是sam文件的二进制格式,sam 文件是Sequence Alignment/Map Format的简写,产生于比对之后的数据输出,记录了比对的具体情况。文件中以tab键分割,包括上下两部分:
  • Header section

  • Alignments section






01


Header section

  • 该部分全部以“@”开头,提供基本的软件版本,参考序列信息,排序信息等
  • @HD行:这一行中有各种不同的标识
  • 标识“VN”用以说明格式版本
  • 标识“SO”用以说明比对排序的情况,有unknown (default)、unsorted、queryname和coordinate,对于coordinate,排序的主键是Alignments section的第三列“RNAME”,其顺序由@SQ行的“SN”标识的顺序定义,次要排序键是Alignments section的第四列“POS”字段。对于RNAME和POS相等的比对,排列顺序则是任意的
  • @SQ行的“SN”标签是参考序列说明,它的值主要是用于Alignments section的第三列“RNAME”和第七列“MRNM”比对的记录
  • @PG行是使用的程序说明;该行“ID”为程序记录标识符,“PN”为程序名字,“CL”为命令行
  • @CO行是任意的说明信息


02


Alignments section

该部分包含了11列必需字段,无效或者没有的字段一般用“0”或者“*”表示。

下面,以数据详细介绍下bam文件的Alignments section中的每一列的含义。

下图是reads比对情况的示意图

  • Coor是参考基因组上的位置坐标

  • ref为参考序列

  • +r001/1和-r001/2组成1对Reads

  • +r003和-r003是反向互补关系,其实是一条read

  • +r002、+r003和+r004都是单条read比对

下图是将比对情况记录成bam文件格式的样子

图中的Alignments section有6行12列信息详细介绍了6条read的比对情况,其中前11列为必需字段,每列的含义简单汇总成下表。接下来,我们仔细解读下表1中的信息。

第1列:Qname

Read的名字

第2列:FLAG

每一个read的比对情况可以用十进制数字(或者十六进制数字)表示,如果比对情况 有多个,将多个比对情况所代表的十进制数字加和就是这一行的FLAG。比如,图1中r001的FLAG是99(1+2+32+64),则表示了“该read是pair read中的一个”,“pair read中每个都能够正确比对上”,“该read的mate read的反向互补可以比对上”,“该read是pari read中的read1”;r001的另一个FLAG是147(1+2+16+128),则表示“该read是pair read中的一个”,“pair read中每个都能够正确比对上”,“该read是原read的反向互补”,“该read是pari read中的read2”(也就是说,该read是read2的反向互补序列)。值得注意的是,r001是pair read,而且都能比对上,所以r001出现了两次,如果r001的read1比对到参考序列的2个地方,r001的名字则会出现三次;如果read1比对上一次,read2没有比对上,r001仍会出现2次,不过,其中一个r001的第三列为“*”;所以pair-end测序,read1文件和read2文件同时mapping,相同reads的id最少出现2次。另,以下网站可以通过输入FLAG值,直接找出该FLAG是那些FLAG的加和:Decoding SAM flagshttps://broadinstitute.github.io/picard/explain-flags.html

FLAG释义表如下:

第3列:RNAME

比对上的参考序列的名字,该名字出现在Header section的@SQ行的SN标识中,如果该read没有比对上,也就是说该read在参考序列上没有坐标,那么这一列则用“*”表示,那么这一行的POS和CIGAR列也会是“*”。

第4列:POS

read比对到的参考序列“RNAME”最左侧的位置坐标,也是CIGAR中第一个比对标识“M”对应的最左侧碱基在参考序列的位置,未比对上的read在参考序列中没有坐标,此列标识为“0”。

第5列:MAPQ

比对的质量值,计算方法为比对错误率的-10*log10的值,一般是四舍五入的整数值,如果是255,说明该比对值无效。

第6列:CIGAR

CIGAR标识符表示read中每个碱基的比对情况,主要有以下标识符:

  • M: alignment match (can be a sequence match or mismatch)

    read上的碱基与参考序列“RNAME”完全匹配,碱基一一对应,包括了正确匹配与错误匹配

  • I: insertion to the reference

    read上的碱基相对于参考序列“RNAME”有插入现象(如下):

REF:  CACGATCA**GACCGATACGTCCGA 
READ1:  CGATCAGAGACCGATA 
CIGAR6M2I8M
  • D: deletion from the reference

    read上的碱基相对于参考序列“RNAME”有删除现象(如下):

REF: AGCTAGCATCGTGTCGCCCGTCTAGCATACGC
READ:             GTGTAACC-GACTAG
CIGAR:8M1D6M
  • N: skipped region from the reference

    read上的碱基相对于参考序列“RNAME”存在连续没有比对上的空缺,这些空缺用N来表示,跟“D”相似但远比“D”缺失的更多,这种比对类型也叫“Spliced alignment”,常见cDNA与参考序列比对(如下):"…"表示intron

REF:AGCATCGTGTCGCCCGTCTAGCATACGCATGATCGACTGTCAGCTAGTCAGACTA
READ:     GTGTAACCC................................TCAGAATA
CIGAR:9M32N8M
  • S: soft clipping (clipped sequences present in SEQ)

  • H: hard clipping (clipped sequences NOT present in SEQ)

    read的开头或者结尾部分没有比对到参考序列"RNAME”上,但这部分未比对上的连续序列仍保留在sam文件的该read序列中,用“S”来表示;如果未保留,则用“H”表示,也即“hard cliping”(如下所示,也可同图2中r003的比对CIGAR中看出)

REF: AGCTAGCATCGTGTCGCCCGTCTAGCATACGCAT
READ:          gggGTGTAACC-GACTAGgggg
CIGAR:3S8M1D6M4S (在sam中存储为GGGGTGTAACCGACTAGGGGG)
CIGAR:3H8M1D6M4H (在sam中存储为GTGTAACCGACTAG)
  • P: padding (silent deletion from padded reference)

    多条read比对到参考序列的同一位置时,如果不同read单独同该参考序列比对时,参考序列的情况也不同,比如下方READ1同参考序列比对时,“GA”属于插入(6M2I8M),READ2同参考序列比对时,“A”属于插入(4M1I9M ),READ3同参考序列完全匹配(10M),没有插入,但是三条read之前却没有可比性。因此,当参考序列“比对情况包含完整”且序列唯一时,所有read同时进行比对,read3这种原本没有插入却默认插入的比对称之为Padded alignment,这种情况用“P”表示。

REF:  CACGATCA**GACCGATACGTCCGA 
READ1:  CGATCAGAGACCGATA 
READ2:    ATCA*AGACCGATAC 
READ3:   GATCA**GACCG  
The padded CIGAR are different: 
READ1: 6M2I8M 
READ2: 4M1P1I9M 
READ3: 5M2P5M


  • =:sequence match 正确匹配

  • X:sequence mismatch 错误匹配

第7列:MRNM

该read的mate read比对上的参考序列的名字,该名字出现在Header section的@SQ行的SN标识中,
  • 如果和该read所在行的第三列“RNAME”一样,则用“=”表示,说明这对read比对到了同一条参考序列上;
  • 如果mate read没有比对上,第七列则用“*”表示;
  • 如果这对read没有比对到同一条参考序列,那么这一列则是mate read所在行第三列的“RNAME”。

第8列:MPOS

该read的mate read比对到的参考序列“RNAME”最左侧的位置坐标,也是mate read CIGAR中第一个比对标识“M”对应的最左侧碱基在参考序列的位置,未比对上的read在参考序列中没有坐标,此列标识为“0”。

第9列:ISIZE

表示pair read完全匹配到同一条参考序列时,两个read之间的长度,可简单理解为测序文库的长度。这个定义有两种情况(虚线表示未比对上的序列,即soft-clipped bases):
  1. 第一种情况如图4左所示,Read1和Read2比对到同一条参考序列,此时ISIZE即为Read2的最右侧比对坐标减去Read1最左侧比对坐标;
  2. 第二种情况如图4右所示,由于至今没有明确的定义和共识,因此ISIZE可以是TLEN#1,也可以是TLEN#2,视情况而定。



第10列:SEQ

存储的序列,没有存储,此列则用“*”标识。该序列的长度一定等于CIGAR标识中“M”,“I”,“S”,“=”,“X”标识的碱基长度之和。

第11列:QUAL

序列的每个碱基对应一个碱基质量字符,每个碱基质量字符对应的ASCII码值减去33(Sanger Phred-33 质量值体系),即为该碱基的测序质量得分(Phred Quality Score)。不同Phred Quality Score代表不同的碱基测序错误率,如Phred Quality Score值为20和30分别表示碱基测序错误率为1%和0.1%。

【参考资料】

  • Sam格式详解手册:http://samtools.github.io/hts-specs/SAMv1.pdf

  • Sam格式相关文献:https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2723002/

  • Sam格式的wiki介绍(着重介绍了比对类型):https://davetang.org/wiki/tiki-index.php?page=SAM

  • Sam格式相关的博客:https://blog.csdn.net/genome_denovo/article/details/78712972

  • 快速查询FLAG的网站:https://broadinstitute.github.io/picard/explain-flags.html

  • ASCII码对应表:https://snipcademy.com/binf/img/lessons/dna-file-formats/ascii.svg

  • 关于测序质量体系的博客:https://blog.csdn.net/xxxxy314/article/details/48948891/


作者:麦茬道

审稿:童蒙

编辑:amethyst

往期回顾


您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存