生信文件格式-SAM文件

SAM文件格式

SAM（The Sequence Alignment / Map format）格式，即序列比对文件的格式，详细介绍文档

以下内容参考2019年1月SAM文件说明文档，具体细节请关注最新文档说明

SAM文件由两部分组成：头部信息和比对信息，都是以tab键分隔。

1. 头部信息介绍
   每个标题行以字符"@"开头，后面是两个字母的记录类型代码。在标题中，每一行都是由制表符分隔的。除了@CO行，每个数据字段都遵循格式"TAG:VALUE"，其中TAG是两个字母组成的字符串，定义了内容和值的格式。每个标题行应该匹配:/^ @(HD|SQ|RG|PG)(\t[A-Za-z][A-Za-z0-9]:[-~]+)+$ /或/^@CO\t.*/。小写字母是保留给用户自定义内容使用的，不会在任何版本定义中出现。以下给出了定义的记录类型和标记。当以下TAG标签出现"*"时，表示该TAG是必须的，如@SQ必须有SN和LN字段。

   • @HD 首行，输出文件的第一行
     VN* 格式版本，接受的格式：/^[0-9]+.[0-9]+$/
     SO 比对排序类型，有unknown (default), unsorted, queryname和coordinate，对于coordinate，排序的主键是RNAME，其顺序由标题中的@SQ行顺序定义，次要排序键是POS字段。如果RNAME和POS相同，顺序是任意随机的。

     # 例如
     @HD VN:1.4	SO:coordinate
     

   • @SQ 参考序列字典，@SQ行的顺序定义了比对文件排序顺序。
     SN* 参考序列名字(染色体)。每一个@SQ行必须含有一个去重的SN标签。这个字段的值是用于RNAME和PNEXT比对的记录。正则表达式[!-)+-<>-][!-]*
     LN* 参考序列长度，范围:[$1,2^{31}-1$]
     AS 基因组装配标识符。
     M5 MD5算法校验序列码
     SP 物种
     UR 参考序列的URI。这个值从一个标准的协议开始，例如，http:或ftp:，如果不是这些协议作为开始，那么就应该是文件系统路径

     # 例如
     @SQ	SN:chr1	LN:249250621	M5:1b22b98cdeb4a9304cb5d48026a85128	UR:file:/opt/biosoft/ref/sv/hg19_ref_genome.primary-assembly.fa
     

   • @RG Read Group，1个样本的测序结果为1个Read Group。允许有多个无序的@RG行
     ID* read group标识符。每一个@RG行必须含有一个唯一的ID。
     CN 测序中心提供的read名称
     DS 描述
     DT 测序的日期 (ISO8601 date or date/time)
     LB 文库名
     PL 用于产生reads的平台/技术。可用值: CAPILLARY, LS454, ILLUMINA, SOLID, HELICOS, IONTORRENT and PACBIO
     PU 平台单元（flowcell-barcode.lane for Illumina or slide for SOLiD）。唯一标识符。
     SM 样本名

     @RG	ID:1	LB:Lib	PL:Illumina	SM:sample_L001	PU:Lane
     

   • @PG 比对所使用的软件程序及版本
     ID* 程序记录标识符。每个@PG必须有一个唯一的ID
     PN 程序名字
     CL 比对操作的命令行内容
     VN 程序版本

     # 例如
     @PG	ID:bwa	PN:bwa	VN:0.7.12-r1039	CL:bwa mem -t 4 -M /opt/ref/hg19.fa read1.fq.gz read2.fq.gz
     

   • @CO 单行的text描述，是一个任意的说明信息。允许多个@CO行无序排列

2. 比对信息介绍
   比对信息部分(alignment section)，每一行表示一个片段（segment）的比对信息，包括11个必须的字段（mandatory fields）和一个可选的字段，字段之间用tab分割。必须的字段有11个，顺序固定，不可用时，根据字段定义，可以为‘0‘或者‘*‘，这是11个字段包括：

   1. QNAME 比对片段的（template）的编号；read name，read的名字通常包括测序平台等信息
   # eg.E00528:375:HN2JTCCXY:8:1123:5457:24479
   2. FLAG  位标识，template mapping情况的数字表示，每一个数字代表一种比对情况，这里的值是符合情况的数字相加总和
   # eg.16
   3. RNAME 参考序列的编号，如果头部注释中对SQ-SN进行了定义，这里必须和其保持一致，另外对于没有mapping上的序列，这里是‘*‘；
   # eg.chr1
   4. POS  比对上的位置，注意是从1开始计数，没有比对上，此处为0;
   # eg.36576599
   5. MAPQ mapping的质量，比对的质量分数，越高说明该read比对到参考基因组上的位置越唯一;
   # eg.42
   6. CIGAR 简要比对信息表达式（Compact Idiosyncratic Gapped Alignment Report），其以参考序列为基础，使用数字加字母表示比对结果，match/mismatch、insertion、deletion 对应字母 M、I、D。比如3S6M1P1I4M，前三个碱基被剪切去除了，然后6个比对上了，然后打开了一个缺口，有一个碱基插入，最后是4个比对上了，是按照顺序的；
   # eg.36M 表示36个碱基在比对时完全匹配
   
   ###注：第七列到第九列是mate的信息，若是单末端测序这几列均无意义###
   
   7. RNEXT 配对片段（即mate）比对上的参考序列的编号，没有另外的片段，这里是‘*‘，同一个片段，用‘=‘；
   # eg.*
   8. PNEXT 配对片段（即mate）比对到参考序列上的第一个碱基位置，若无mate,则为0；
   # eg.0
   9. TLEN  Template（文库插入序列）的长度，最左边的为正，最右边的为负，中间的不用定义正负，不分区段（single-segment)的比对上，或者不可用时，此处为0(ISIZE，Inferred fragment size.详见Illumina中paired end sequencing 和 mate pair sequencing，是负数，推测应该是两条read之间的间隔(待查证)，若无mate则为0);
   # eg.0
   10. SEQ  序列片段的序列信息，如果不存储此类信息，此处为‘*‘，注意CIGAR中M/I/S/=/X对应数字的和要等于序列长度；
   # eg.CGTTTCTGTGGGTGATGGGCCTGAGGGGCGTTCTCN
   11. QUAL 序列的质量信息, read质量的ASCII编码。格式同FASTQ一样。
   # eg.F-?DD<,EEFFEEC?F5FA7FAF8FAF;F@E:DF>FEEAFFFFFCFFC4EFDFDDFEFE=EFFF
   12.第十二列及之后：Optional fields，以tab建分割。
   # eg.AS:i:-1 XN:i:0 XM:i:1 XO:i:0 XG:i:0 NM:i:1 MD:Z:35T0 YT:Z:UU
   

3. FLAG取值

1. 关于Supplementary的定义，可以看文档对chimeric alignment（嵌合比对）的描述
chimeric alignment:“嵌合比对” 的形成是由于一条测序read比对到基因组上时分别比对到两个不同的区域，而这两个区域基本没有overlap。因此它在sam文件中需要占用多行记录显示。只有第一个记录被称作"representative",其他的都是"supplementary"【Chimeric reads are also called split reads】

2. FLAG举例说明
FLAG 419 表示PAIREED(1), PROPER_PAIR(2),MREVRSE(32),READ2(128),SECONDARY(256)，其中

PAIRED 代表这条序列采用双端测序， 其值为1；

PROPER_PAIR 代表这条序列完全匹配， 其值为2；

MREVRSE 代表这条序列对应的另一端序列比对到参考序列的负链上，其值为32；

READ2 代表这条序列是R2端序列，其值为128

SECONDARY 代表这条序列不是primary alignment, 其值为256

1 + 2 +32 +128 + 256 = 419

3. bam2fq
根据FLAG搞清楚了序列究竟是R1端还是R2端之后，还有一点值得注意的地方，samtools还提供了一个bam2fq的功能，根据bam文件抽取比对时采用的fastq序列

在bam文件中，第10列代表的是序列，第11列代表的是序列的质量，根据第二列的flag还可以确定序列是R1还是R2，但是当序列本身比对到参考序列的负链时，
即flag 包含16时，bam文件中记录的序列是原始序列的反向互补序列，而且质量值也是反向的，所以根据这样的序列还原时要小心一点，以flag为339的序列为例，
因为339包含了REVERSE,所以对应的序列应该是第10列序列的反向互补序列，碱基质量值为第11列的反向序列，

4. CIGAR字符串

Op	BAM	描述	                                消耗待比对序列	消耗参考序列
M	0	位置能比对上                                 yes	            yes
I	1	相对参考序列有插入	                    yes	            no
D	2	相对参考序列有缺失	                    no	            yes
N	3	从参考序列上跳过一段	                    no	            yes
S	4	软切割（被切割的序列保留在SEQ中）	            yes	            no
H	5	硬切割（被切割的序列不出现在SEQ中）	    no	            no
P	6	补丁（打了补丁的参考序列中的沉默缺失）	    no	            no
=	7	read碱基与参考序列相同	                    yes	            yes
X	8	read碱基与参考序列不同	                    yes	            yes

1. "消耗查询序列"与"消耗参考序列"分别指CIGAR是否引起比对沿着查询序列和参考序列的方向向前前进一个或几个碱基；
2. H 只能出现在CIGAR的开始或最后；
3. S的两边必为H，否则必须位于CIGAR的两端；
4. 对于mRNA到基因组的比对，N表示内含子。对于其他类型的比对，N的解释未被定义；
5. MIS=X的长度和应该等于SEQ长度。

5. 比对信息第12列及之后：可选字段
   所有的可选字段服从TAG:TYPE:VALUE格式，其中TAG是两个字符组成的字符串，符合正则表达式/[A-Za-z][A-Za-z0-9]/。每个TAG在一个比对行中只能出现一次。小写字母组成的TAG保留给终端用户自定义。在可选字段中，TYPE是大小写敏感的单个字母，用来定义VALUE的格式

   TYPE	相应VALUE的正则表达式	描述
   A	[!-~]	可打印字符
   i	[-+]?[0-9]+	有符号的整数（SAM格式对这里整数的范围没有要求，但BAM要求范围是[-2^31, 2^32 ]，所以SAM也得在这个范围内）
   f	[-+]?[0-9]*.?[0-9]+([eE][-+]?[0-9]+)?	单精度浮点数
   Z	[ !-~]*	可打印字符串，包括空格
   H	([0-9A-F][0-9A-F])*	十六进制格式的二进制数列表（例如，二进制列表[0x1a, 0xe3, 0x1]对应十六进制串 ‘1AE301’）
   B	[cCsSiIf](,[-+]?[0-9]*.?[0-9]+([eE][-+]?[0-9]+)?)*	整数或数字数组

   对于一个整数或数字列表（类型B），第一位字母指明后面逗号分隔的列表的数字类型，可以是‘cCsSiIf’中的一个，分别对应int8 t (有符号的8位整数), uint8 t (无符号的8位整数), int16 t, uint16 t, int32 t, uint32 t and float（Explicit typing eases format parsing and helps to reduce the file size when SAM is converted to BAM）。输入输出过程中，如果其他类型也能与列表兼容，列表元素的类型可能会改变。

6. 可选字段：预定义的标签（tags）
   在SAM可选字段说明中有描述，其中描述了已有的标准TAG字段和自定义字段的细节。以X、Y、Z开头的TAG和包含小写字母的TAG保留为终端用户自定义使用，不会在本说明未来版本中正式定义。

   ## 例如
   AS:i 匹配的得分
   XS:i 第二好的匹配的得分
   YS:i mate 序列匹配的得分
   XN:i 在参考序列上模糊碱基的个数
   XM:i 错配的个数
   XO:i gap open的个数，针对于比对中的插入和缺失
   XG:i gap 延伸的个数，针对于比对中的插入和缺失
   NM:i 编辑距离。但是不包含头尾被剪切的序列。一般来说等于序列中error base的个数
   YF:i 该reads被过滤掉的原因。可能为LN(错配数太多，待查证)、NS(read中包含N或者．)、SC(match bonus低于设定的阈值)、QC(failing quality control，待证)
   YT:Z 值为UU表示不是pair中一部分(单末端？)、CP(是pair且可以完美匹配)、DP(是pair但不能很好的匹配)、UP(是pair但是无法比对到参考序列上)
   MD:Z 比对上的错配碱基的字符串表示
   
   ## 比对软件会定义一些自定义可选字段，比如BWA，详细情况请查看bwa说明文档
   ## bwa定义，如
   XT:A: 比对Type: Unique/Repeat/N/Mate-sw 如 XT:A:U 表示唯一比对
   XM:i: 比对中mismatch的数目
   ...
   
   

生信文件格式-SAM文件

原文：https://www.cnblogs.com/dataanaly/p/13177364.html