---------------
Nickier
2019-01-12
---------------
sam是一种序列比对后的输出格式，以tab作为分隔符，包括头部信息和比对信息。其中头部信息必须在比对信息之前。头部信息的开头是@，但是比对行不是。每一个比对行有11个重要的比对信息元素，如果比对位置和校准信息等。
关于sam/bam文件的格式，网上已经有很多介绍，但是很多没有结合igv可视化，所以不容易理解，今天，我们将sam/bam文件结合igv可视化，方便大家理解。
Example
ST-E00522:463:HT3TGCCXY:1:1210:30959:4702 163 chr1 10070 0 44S61M = 10028 61 CCCTAATCCCTAACCCCTCACCCTCACCCTCACCTTCACCCTCCCCCTAATCCTAACCCTAACCCTAACCCTAACCCTAACCCAACCCTAACCCTAACCCTAACC .==FB<+-+,CA<9*A8<-146<DFD@8==<BE+18@BCC6=E*@E5E+87D2E,AFFFECEFFFF9/DE8>BEF@FF4/A6A>/6*C14A5@EG+/>;C??DDD MC:Z:6M1I6M1I91M MD:Z:6C54 PG:Z:MarkDuplicates RG:Z:HT_11_DHE16602 NM:i:1 MQ:i:0 AS:i:56 XS:i:54
头文件
每个标题行以字符“@”开头，后面是两个字母的记录类型代码。在标题中，每一行都是由制表符分隔的，除了@CO行，每个数据字段都遵循格式“TAG:VALUE”，其中TAG是一个两个字母的字符串，定义了内容和值的格式。每个标题行应该匹配:/ ^ @[A-Za-z][A-Za-z](\ t[A-Za-z][A-Za-z0-9]:[- ~]+)+ $ /或/ ^ @CO \ t。* /。包含小写字母的标记保留给最终用户。
主要信息
每一个比对行有11个必填选项。这些字段都是以相同顺序出现，而且必须出现，但是这些值可以为0或*（取决于字段）如果无法获得相应的信息。下表概述了SAM格式的强制字段：
| Col | Field | Brief description |
|---|---|---|
| 1 | QNAME | 查询模板名称 |
| 2 | FLAG | 位标记，template mapping情况的数字表示，每一个数字代表一种比对情况，这里的值是符合情况的数字相加总和 |
| 3 | RNAME | 参考序列名称 |
| 4 | POS | 基于1的最左比对位置 |
| 5 | MAPQ | MAPping质量 |
| 6 | CIGAR | CIGAR字符串 |
| 7 | RNEXT | 比对到的参考（染色体）名字 |
| 8 | PNEXT | 配对到的第一个碱基的位置 |
| 9 | TLEN | 可以理解为文库插入片段长度 |
| 10 | SEQ | 序列片段 |
| 11 | QUAL | phred -scale基本质量+33的ASCII码 |
1.QNAME：查询模板名称。具有相同QNAME的read片段被认为来自相同的模板。QNAME ' * '表示信息不可用。
2.FLAG 位标记，下表是每一个代号代表的意义：
| Bit | Description |
|---|---|
| 1 | read是pair中的一条（read表示本条read，mate表示pair中的另一条read） |
| 2 | pair一正一负完美的比对上 |
| 4 | 片段未比对上 |
| 8 | mate没有比对上 |
| 16 | 这条read反向比对 |
| 32 | mate反向比对 |
| 64 | 这条read是read1 |
| 128 | 这条read是read2 |
| 256 | 第二次比对 |
| 512 | 没有通过质量控制 |
| 1024 | read是PCR或光学副本产生 |
| 2048 | 辅助比对结果 |
如下面这两条互相配对的reads中
ST-E00522:463:HT3TGCCXY:1:1217:29264:10662 99 chr1 10004 4 102M1I44M3S = 10195 236 CCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCTAACCCTAACCCGAACCCGAACCCGAACCCCAACCCCAACCCGAACCCCAACCCGAACCCTAACCCTAACCCTAACCCTAA >CDHDEEEEFBDEEEFBEDCEFBEEEEEBEEECFBEDEEDCFDFFDCECCE4:E-->4=BCB9-?>?>F4.EABA@>E:C9&<E(-B0<B@CF0CE*CC,11:B>*==DCC&DE9C;,=@DAB&<EBCC2>D?A?F9A?>8@AFC;>57A XA:Z:chr22,-50808000,150M,11;chr22,-50808169,3S89M1I12M8D6M1D39M,16;chr1,+180761,41M2D44M3I17M3D45M,15;chr3,+10614,28M1I43M1D6M1D28M1D29M1D15M,11;chr22,-50808294,48M5D27M1I74M,16; MC:Z:7S45M MD:Z:63T17T5T5T5T10T5T5T23 PG:Z:MarkDuplicates RG:Z:HT_11_DHE16602 NM:i:9 MQ:i:0 AS:i:99 XS:i:95
ST-E00522:463:HT3TGCCXY:1:1217:29264:10662 147 chr1 10195 0 7S45M = 10004 -236 CCCCCACCCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACC 7F=0+.<>E?;+<F?CBEBEDBBEB?/9ACBDEEAD@DC@A.0D@CCB03?< MC:Z:102M1I44M3S MD:Z:45 PG:Z:MarkDuplicates RG:Z:HT_11_DHE16602 NM:i:0 MQ:i:4 AS:i:45 XS:i:47
第一条reads的第二列是99，即
99
=64（这条reads是read1）
+32（这条reads是正向比对）
+2（pair一正一负完美的比对上）
+1（这条reads是pair中的一条）
在上图中，可以看到linux界面中第一条和第二条reads即就是前面展示相互配对的reads。可以看到两条reads都比对到了chr1上了，此时就得到了
1和2。而在igv中可以看到第一条reads（鼠标所在处）的方向是向右侧，即正向比对，就得到了数值32。至于64的由来，请看下图。
3.RNAME：比对的参考序列名称，如果@SQ头部行存在，RNAME(如果不是“”)必须出现在一个 SQ-SN标记中。没比对上此处就是“”。然而，一个未必对的片段也有一个坐标以便排序。如果RNAME 是“*”，也就没有 POS 和 CIGAR。
4.POS：于1的第一个匹配基的最左映射位置。参考序列中的第一个基的坐标是1。对于没有坐标的未映射读取，POS设置为0。如果POS为0，RNAME和CIGAR也就没有意义。
5.MAPQ：mapping质量，等于−10log 10Pr(映射位置是错的)，四舍五入到最近的整数，值255表示映射质量特别差。
6.CIGAR：CIGAR字符串。下表为CIGAR字符串的解释（‘*’表示无值）
| Op | BAM | Description |
|---|---|---|
| M | 0 | 比对匹配(可以是序列匹配或不匹配) |
| I | 1 | 插入到参考 |
| D | 2 | 从参考删除 |
| N | 3 | 参考的跳过的区域 |
| S | 4 | 软剪切(被剪切的序列存在于序列中) |
| H | 5 | 硬剪切(被剪切的序列不存在于序列中) |
| P | 6 | 填充(从填充引用中无声删除) |
| = | 7 | 序列匹配 |
| X | 8 | 序列不匹配 |
- H 值出现在最初或者最后操作中
- S 可证在他们和CIRAG末尾字符串中只有H操作
- 对于mRNA到基因组的比对，一个N操作符代表内含子。对于其他类型的比对，没有定义N的解释。
- M/I/S/=/X操作的长度之和等于SEQ的长度。
7.RNEXT:：mate的reference sequence name，实际上就是mate比对到的染色体号，若是没有mate，则是*
8.PNEXT：如果没有这个信息（没比对上）就是0
9.TLEN：如果R1端的read和R2端的read能够mapping到同一条Reference序列上（即第三列RNAME相同），则该列的值表示第8列减去第4列加上第6列的值，R1端和R2端相同id的reads其第九列值相同，但该值为一正一负，R1文件的reads和R2文件的reads，相同id的reads要相对来看。在进行该第列值的计算时，如果取第6列的数值，一定要取出现M的值，S或H的值不能取。
10.SEQ：reads片段。如果序列不存在，就是。如果不是，这个序列的长度等于CIGAR中 M/I/S/=/X的总和。=表示他的基础字段（如开始为1），和参考序列的基础字段相同。
11.QUAL：碱基质量加33的ASCII码（与Sanger FASTQ格式中的质量字符串相同）。一个碱基质量是基于错误率的phred-scaled等于−10log 10Pr(碱基是错的)。这个字段可以是“”。如果它不是，那么seq也不是“*”，它的长度与SEQ的长度得一致。
