sam/bam格式文件详解
1,sam文件格式介绍
sam(the sequence alignment / map format)格式,即序列比对文件的格式,详细介绍文档:http://samtools.github.io/hts-specs/samv1.pdf
sam文件由两部分组成,头部区和主体区,都以tab分列。
头部区:以’@'开始,体现了比对的一些总体信息。比如比对的sam格式版本,比对的参考序列,比对使用的软件等。
主体区:比对结果,每一个比对结果是一行,有11个主列和一个可选列。
2,头部区简要介绍
@hd vn:1.0 so:unsorted (排序类型)
头部区第一行:vn是格式版本;so表示比对排序的类型,有unknown(default),unsorted,queryname和coordinate几种。samtools软件在进行行排序后不能自动更新bam文件的so值,而picard却可以。
@sq sn:contig1 ln:9401 (序列id及长度)
参考序列名,这些参考序列决定了比对结果sort的顺序,sn是参考序列名;ln是参考序列长度;每个参考序列为一行。
例如:@sq sn:nc_000067.6 ln:195471971
@rg id:sample01 (样品基本信息)
read group。1个sample的测序结果为1个read group;该sample可以有多个library的测序结果,可以利用bwa mem -r 加上去这些信息。
例如:@rg id:zx1_id sm:zx1 lb:pe400 pu:illumina pl:miseq
id:样品的id号 sm:样品名 lb:文库名 pu:测序以 pl:测序平台
这些信息可以在形成sam文件时加入,id是必须要有的后面是否添加看分析要求
@pg id:bowtie2 pn:bowtie2 vn:2.0.0-beta7 (比对所使用的软件及版本)
例如:@pg id:bwa pn:bwa vn:0.7.12-r1039 cl:bwa sampe -a 400 -f zx1.sam -r @rg id:zx1_id sm:zx1 lb:pe400 pu:illumina pl:miseq …/0_reference/reference_sequence.fa zx_hq_clean_r1.fq.sai zx_hq_clean_r2.fq.sai …/2_hqdata/zx_hq_clean_r1.fq …/2_hqdata/zx_hq_clean_r2.fq
这里的id是bwa,pn是bwa,vn是0.7.12-r1039版本。cl可以认为是运行程序@rg是上面rg表示的内容,后面是程序内容,这里的@gr内容是可以自己在运行程序是加入的
3,主体部分介绍
主体部分有11个主列和1个可选列
qname 比对的序列名称 例如:m04650:84:000000000-b837r:1:1101:22699:1759(一条测序reads的名称)
flag bwise flag(表明比对类型:paring,strand,mate strand等) 例如:99
rename 比对上的参考序列名 例如:nc_000075.6
pos 1-based的比对上的最左边的定位 例如:124057649
mapq 比对质量 例如:60
cigar extended cigar string(操作符:midnshp)比对结果信息;匹配碱基数,可变剪接等 例如:87m
mrnm 相匹配的另外一条序列,比对上的参考序列名 例如:=
mpos 1-based leftmost mate position (相比于mrnm列来讲意思和pos差不多) 例如:124057667
isize 插入片段长度 例如:200
seq 和参考序列在同一个链上比对的序列(若比对结果在负义链上,则序列是其反向重复序列,反向互补序列) 例如:attacttggctgct
qual 比对序列的质量(ascii-33=phred base quality)reads碱基质量值 例如:-8cccgfcccf7@e-
可选的列 以tag:type:value的形式提供额外的信息
4,对于每一列内容的详细注解
(如果某一列为“0”或“*”表示这一列没有信息)
第一列:qname
进行reads比对时通常表示reads的名字,如果这条reads比对到多条序列或比对到这条序列的多个位置,相同名字会出现多次。如果是pair-end reads,相同名字会出现2次,分别表示来自于r1文件的reads和r2文件的reads,如果其matepair reads也比对2个位置,也会出现2次,则相同名字共出现4次,如果一条reads也比对2个位置,则其matepair比对1个位置,则共出现3次,如果其matepair reads没有比对上序列也会出现1次(第三列显示“*”),所以pair-end测序,r1文件和r2文件同时mapping,相同reads的id最少出现2次。
第二列:flag
数值结果如下:
1(1)该read是成对的paired reads中的一个
2(10)paired reads中每个都正确比对到参考序列上
4(100)该read没比对到参考序列上
8(1000)与该read成对的matepair read没有比对到参考序列上
16(10000)该read其反向互补序列能够比对到参考序列
32(100000)与该read成对的matepair read其反向互补序列能够比对到参考序列
64(1000000)在paired reads中,该read是与参考序列比对的第一条
128(10000000)在paired reads中,该read是与参考序列比对的第二条
256(100000000)该read是次优的比对结果
512(1000000000)该read没有通过质量控制
1024(10000000000)由于pcr或测序错误产生的重复reads
2048(100000000000)补充匹配的read
具体的flag值的解释,可以参考samtools软件提供的结果
samtools(version: 1.3.1)
其中的samtools flags用法可提供flag值的查找结果
about: convert between textual and numeric flag representation
usage: samtools flags int|str[,…]
例如:
samtools flags 10
0xa 10 proper_pair,munmap(10=2 8)
samtools flags 12
0xc 12 unmap,munmap(12=4 8)
具体的flag值的解释,也可参考如下网站:https://broadinstitute.github.io/picard/explain-flags.html
或者在必应当中搜索flag sam点击explain sam flags-github pages进入该网页,也可以输入组合flag数值会出现所存在的意思
第三列:rname
表示read比对的那条序列的序列名称(名称与头部的@sq相对应),如果这列是“*”,可以认为这条read没有比对上的序列,则这一行的第四,五,八,九 列是“0”,第六,七列与该列是相同的表示方法
第四列:pos
表示read比对到rname这条序列的最左边的位置,如果该read能够完全比对到这条序列(cigar string为m)则这个位置是read的第一个碱基比对的位置,如果该read的反向互补序列比对到这条序列,则这个位置是read的反向互补序列的第一个碱基比对的位置,所以无论该read是正向比对到该序列,或是其反向互补序列比对到该序列,比对结果均是最左端的比对位置
第五列:mapq
表示为mapping的质量值,mapping quality, it equals -10log10pr{mapping position is wrong}, rounded to the nearest integer, a value 255 indicates that the mapping quality is not available. 该值的计算方法是mapping的错误率的-10log10值,之后四舍五入得到的整数,如果值为255表示mapping值是不可用的,如果是unmapped read则mapq为0,一般在使用bwa mem或bwa aln(bwa 0.7.12-r1039版本)生成的sam文件,第五列为60表示mapping率最高,一般结果是这一列的数值是从0到60,且0和60这两个数字出现次数最多
第六列:cigar
cigar string,可以理解为reads mapping到第三列序列的mapping状态,
对于mapping状态可分为以下几类:
m:alignment match (can be a sequence match or mismatch)
表示read可mapping到第三列的序列上,则read的碱基序列与第三列的序列碱基相同,表示正常的mapping结果,m表示完全匹配,但是无论reads与序列的正确匹配或是错误匹配该位置都显示为m
i:insertion to the reference
表示read的碱基序列相对于第三列的rname序列,有碱基的插入
d:deletion from the reference
表示read的碱基序列相对于第三列的rname序列,有碱基的删除
n:skipped region from the reference
表示可变剪接位置
p:padding (silent deletion from padded reference)
s:soft clipping (clipped sequences present in seq)
h:hard clipping (clipped sequences not present in seq)
clipped均表示一条read的序列被分开,之所以被分开,是因为read的一部分序列能匹配到第三列的rname序列上,而被分开的那部分不能匹配到rname序列上。
"="表示正确匹配到序列上
"x"表示错误匹配到序列上
而h只出现在一条read的前端或末端,但不会出现在中间,s一般会和h成对出现,当有h出现时,一定会有一个与之对应的s出现
例如:
162m89s
162h89m
149m102s
149h102m
40s211m
20m1d20m211h
s可以单独出现,而h必须有与之对应的s出现时才可能出现,不可在相同第一列的情况下单独出现
n:如果是mrna-to-genome,n出现的位置代表内含子,其它比对形式出现n时则没有具体解释
m/i/s/=/x:这些数值的加和等于第10列seq的长度
第七列:mrnm
这条reads第二次比对的位置,在利用bwa mem产生sam文件时,如果该列是“”而
第3列rname不是“”则表示该reads比对到第3列显示序列名的序列上,而没有比对到其他位置,在利用bwa aln及bwa sampe比对生成的sam文件,如果和上述情况相同,则第7列为“=”,上述情况均表示该reads只比对到这一个位置
如果第3列rname和第7列mrnm都为“*”,则说明这条reads没有匹配上的序列,如果这条reads匹配两个序列,则第一个序列的名称出现在第3列,而第二个序列的名称出现在第7列
第八列:mpos
该列表示与该reads对应的mate pair reads的比对位置,如果这对pair-end reads比对到同一条reference序列上,在sam文件中reads的id出现2次,read1比对的第4列等于read2比对的第8列。同样read1比对的第8列等于read2比对的第4列。例如:
第1列(read id)····第4列(read1比对位置)····第8列(mate-pair reads比对位置)
22699:1759····124057649····124057667
22699:1759····124057667····124057649
相同的reads id一个来自read1文件,一个来自read2文件,第4列和第8列是对应的
第九列:isize
tlen:signed observed template length (可以理解为文库插入片段长度)
如果r1端的read和r2端的read能够mapping到同一条reference序列上(即第三列rname相同),则该列的值表示第8列减去第4列加上第6列的值,r1端和r2端相同id的reads其第九列值相同,但该值为一正一负,r1文件的reads和r2文件的reads,相同id的reads要相对来看。在进行该第列值的计算时,如果取第6列的数值,一定要取出现m的值,s或h的值不能取。
the unisgned observed template length equals the number of base from the leftmost mapped base to the rightmost mappedbase. theleftmost segment has a plus sign and the rightmost has a minus sign
处理bam文件的主要生信软件有
bwa,bowtie2,samtools,bedtools等
可以看mapping等多方面结果和统计,bedtools工具中genomecoveragebed的功能是:compute the coverage of a feature file among a genome