一、sam、bam简介
1、sam
Sequence Alignment/Map format，直译就是序列比对格式；
-
序列比对是测序结果分析中最基础的一步操作，它主要基于参考序列信息回答了每条reads来自人的哪条染色体的哪个区域；
sequence mapping 目前的比对软件很多，例如bwa，star，bowtie2，算法也各有差异，但究其比对结果都是以sam纯文本格式展现；
在sam格式中，对比对情况、比对质量、PE双端比对等各种比对可能结果都有全面的记录，具体会在下面介绍。
2、bam
- 简单来说，bam格式就是sam的二进制版本；
- 因为，sam最大的缺点就是文件体积太大，而bam格式可完整保留sam信息同时，可明显降低文件大小；所以一般仅保留bam文件进行后续分析。
- samtools可实现sam与bam间的转换，以及查看bam内容，具体操作示例见笔记最后
3、序列比对小例子
-
如下图所示，准确来说是5条reads的比对结果，其中具体比对情况为--
image.png
（1）r001/1与r001/2是PE测序的pair/mate reads，现在遇到的基本都是双端测序；
（2）r002、r003、r004都是单独的一条read(在这个图片结果是这样，也许对应的pair reads比对到很远的位置，甚至其它染色体上)；
（3）r003出现两次是因为其可比对到参考序列的不同位置，但是注意 -r003表示是其reverse reads比对上的；
（4）reads序列中的大写字母表示与参考序列对应碱基相同，小写字母则表示不一致；星号或点号表示考虑insert/delete情况，reads才能比对到参考序列中； -
而sam格式有系统的规则，可简洁、全面的描述上述比对结果，如下图所示。
image.png
在PE双端测序中 pair reads与mate pairs概念基本相同，区别见https://www.biostars.org/p/77293/
2、sam格式剖析
2.1 header部分（optimal）
- 以
@开头为标志，从不同角度为比对结果补充信息； -
常见的有以下四类
sam header - （1）
@HD如果有，必须放在第一行
VN表示sam格式的版本号，SO说明比对结果是否排序了 - （2）
@SQ说明参考基因组的情况，一般会有很多行
SN表示染色体名字，对应比对结果的RNAME列；LN表示该染色体的长度 - （3）
@PG则说明比对软件情况
ID表示程序ID号(一般与PN相同)，PN表示软件名，VN表示软件版本号，CL表示产生该sam的命令语句 - （4）
@RG主要说明测序样本信息
ID，在一条lane只有一个sample情况下，可用flowcell+lane序号表示；SM表示sample名；LB表示原始测序文库名；PL表示测序平台
2.2 alignment section（required）
- alignment section是sam的主体部分，其中那个每行代表一次比对结果，以及对应详细的比对信息。
-
如下图，每条alignment record可分为9大部分，11+n列(metadata不固定)
alignment record
第一列：read name
在双端测序，一般至少有两条相同记录的read name（pair reads），但有时会有很多条(>2)记录是因为
- 首先是
Chimeric alignment情况，就是类似最一开始介绍的r003的比对情况，由于序列的特性，可能比对到不同的比较合适区域； - 然后是
Multiple mapping情况，这主要是由于reads序列与基因组的重复序列特征导致的（在人基因组中约占20%）。
read name
第二列：flag
- 简单来说就是设置该条reads的比对属性；
- 表示某条reads是否有12条特征中的1至多条(多选题)，结果用10进制数值加和表示
image.png
简单理解所有选项含义---
1表示是双端测序；
2表示pair reads均比对到了参考序列(good!)
4表示这条read没比对到参考序列
8表示这条read的mate read没比对到参考序列
16表示这条read是reverse(反向)比对到参考序列的
32表示这条read的mate read是reverse(反向)比对到参考序列的'
64表示这条read是pair read的第一条read(左边的)
128表示这条read是pair read的第二条read(右边的)
256表示是这条read的Multiple mapping reads比对情况
其余三种情况比较少见，不作记录了。关于secondary alignment与supplementary alignment的区别我认为就是分别对应Multiple mapping reads与Chimeric alignment比对情况。详见https://sourceforge.net/p/samtools/mailman/message/33235303/
- 某一加和结果数值只会是一组属性特征值的和，不会出现不同组合均能加和等于一个特定和的结果；
- 例如下图，具体flag值含义参考上面解释
113=1+16+32+64
369=1+16+32+64+256
117=1+16+32+128
417=1+32+128+256
image.png
这个网页https://broadinstitute.github.io/picard/explain-flags.html提供了根据输入的结果flag来分解原特征组合的功能，挺方便的。
第三列&第四列：read的位置
- RNAME: Reference sequence NAME of the alignment；与header部分的
@SQ SN保持一致； - POS: read比对到参考序列(染色体)的最左边位置(基于染色体，且为1-based)；
- 如果只有pairs reads的一个成功比对到参考序列，那么未比对的reads的
RNAME与POS应与比对上的reads的描述一致； - 如果pair reads均未成功比对到参考序列上，则RNAME列均为
*，POS列均为0
unmapped reads
简单来说1-based就是序列第一个碱基序号为1，例如SAM, VCF, GFF and Wiggle formats等；0-based就是序列第一个碱基序号为0，例如 BAM, BCFv2, BED, and PSL formats等。
第五列：MAPQ
- MAPping Quality. It equals −10log 10 Pr{mapping position is wrong}, rounded to the nearest
integer. - 一般来说MAPQ值越大，表示MAPping Quality越高；但255表示不可计算该read的MAPQ(我理解就是指unmapped的reads)
第六列：CIGAR
- 这一列主要来描述序列的具体比对情况，以碱基为单位，主要有如下图几种比对类型
image.png
M表示read碱基比对到了参考序列的碱基上(有趣的是允许夹杂mismatch的情况，为SNP variant calling基础)
I相比于参考序列，read里插入insert了新碱基
D相比于参考序列，read里删除了部分碱基
N比较特殊一般用于mRNA-seq结果里，表示比对到了intron内含子区域；
S/H表示read序列两端基本没能比对到参考序列，而中间部分可以比对的情况；具体区别可参考http://www.aigenetic.com/index.php/2018/03/19/soft-clip-%E4%B8%8E-hard-clip%E7%9A%84%E5%90%AB%E4%B9%89%E6%8F%8F%E8%BF%B0/
P针对是多条reads比对到同一参考序列区域，其中有一条reads存在insert的情况，具体见https://davetang.org/wiki/tiki-index.php?page=SAM#Padded_alignment - Sum of lengths of the
M/I/S/=/Xoperations shall equal the length of SEQ，其中=X不常见。并且第四列read的起始位置也是从M/I/S/=/X开始计算的
其实总结来看，
N与S/H分别表示了两种特殊比对情况spliced alignment与clipped alignment，前者是中间没比对上，而两端比对上了；后者使中间比对上，而两端没比对上。
clipped alignment -- 3S8M1D6M4S
spliced alignment -- 9M32N8M
-
举个例子，如下图表示3 bases aligned followed by 1 base deleted, 2 next ones aligned, 1 base inserted and the last one aligned.
image.png
第七、八、九列：描述mate read信息
- RNEXT是该read的mate read比对到的参考染色体，有三种情况；
（1）=：read与mate read比对到同一染色体上；
（2）若mate read比对到其它染色体上，则用相应的染色体名称即可（SN）
（3）若mate read未比对成功，用*表示 - PNEXT是该read的mate read比对到的参考染色体的起始位置；
- TLEN列主要描述pair reads的"相隔距离"，如下图。
（1）仅考虑M/I/=/X(excludes soft-clipped bases)情况；
（2）pair reads的该列绝对值相同，只是左边的reads为正值，右边的reads为负值；
（3）如果pair reads分别比对到不同染色体上，那么该值就是0
image.png
特殊的TLEN：Note: these two definitions agree in most alignments, but differ in the case of overlaps where the first segment aligns beyond the start of the last segment.
image.png
第十、十一列：序列信息
- 第十列：原始序列组成，等于原来fastq的第二行；
- 第十一列：序列的Phred质量值，等于原来fastq的第四行；详见之前fastq的笔记。
第12+n列：metadata（optimal）
-
TAG:TYPE:VALUE格式：TAG表示标签名，一般是两个大写字母；TYPE表示VALUE的数据类型；VALUE表示该read的VALUE值 - Predefined standard tags可参考：https://samtools.github.io/hts-specs/SAMtags.pdf
- 以
BWA mem比对结果产生的TAG结果为例，解释如下图。
image.png
至此，关于sam格式的基本介绍大致如上，主要参考了http://samtools.github.io/hts-specs/SAMv1.pdf教程手册，其中还有很多深入的知识，值得以后深入探索~~
三、samtools转换sam、bam
- sam转bam
samtools view -bS SRR1663608.sam > SRR1663608.bam
- bam转sam
samtools view -h -o SRR1663608.sam SRR1663608.bam
- bam结果统计
#查看bam
samtools view -h SRR1663608.bam | more
#the number of records (alignments)
samtools view -c SRR1663608.bam
#Displays basic alignment stats based on flag
samtools flagstat SRR1663608.bam
更多关于sam/bam格式的操作，可见之间生信技能树Jimmy大神布置的一些练习题，我自己也做了，详见Linux生信练习3--sam/bam
笔记图片大多来自网上，侵删~ 笔记中如有错误之处，欢迎指出！
