对于sam文件格式，已经有很多人写过相应的中文介绍文章，但是基本上都是简单粗暴的讲了一下每个字符代表什么，但是如果真去研究的话有时候还是感觉有点这些写的还是让人懵逼，比如不知道sam文件格式内容和比对的关系到底是怎么样的？
sam/bam格式：源于别人的最全的 - 简书 (jianshu.com)
生信文件格式-SAM文件 - WarningMessage - 博客园 (cnblogs.com)
sam/bam格式：源于别人的最全的 - 简书 (jianshu.com)
但实际上李恒大佬已经写的很清楚了
不过值得注意的是由于遵循古老的编程习惯，位置是要从0开始计算的
https://samtools.github.io/hts-specs/SAMv1.pdf
如果开始比对后的序列结果长这样
image.png
那么对应的sam格式其实是这样的，是将上面的比对结果按照一定规则进行编码
image.png
有了这个认知之后再去看这些列就很容易明白了
@HD File-level metadata. Optional. If present, there must be only one @HD line and it must be the first line of the file
@SQ Reference sequence dictionary. The order of @SQ lines defines the alignment sorting order
header内容不多，也不会太复杂，每一行都用‘@’ 符号开头，里面主要包含了版本信息，序列比对的参考序列信息，如果是标准工具（bwa，bowtie，picard）生成的BAM，一般还会包含生成该份文件的参数信息（如上图），@PG标签开头的那些这里需要重点提一下的是header中的@RG也就是Read group信息，这是在做后续数据分析时专门用于区分不同样本的重要信息。它的重要性还体现在，如果原来样本的测序深度比较深，一般会按照不同的lane分开比对，最后再合并在一起，那么这个时候你会在这个BAM文件中看到有多个RG，里面记录了不同的lane，甚至测序文库的信息，唯一不变的一定是SM的sample信息，这样合并后才能正确处理
第一列：read name，read的名字通常包括测序平台等信息；
第二列：sum of flags，比对flag数字之和，比对flag用数字表示，分别为：
- 0 read是Single-read且正向比对
- 1 read是pair中的一条（read表示本条read，mate表示pair中的另一条read）
- 2 pair一正一负完美的比对上
- 4 这条read没有比对上
- 8 mate没有比对上
- 16 这条read反向比对
- 32 mate反向比对
- 64 这条read是read1
- 128 这条read是read2
- 256 第二次比对
- 512 比对质量不合格
- 1024 read是PCR或光学副本产生
- 2048 辅助比对结果
通过这个和可以直接推断出匹配的情况。假如说标记不是以上列举出的数字，比如说83=（64+16+2+1），就是这几种情况值和。
77=64+8+4+1 read1没比对上
141=128+8+4+1 read2没比对上
当然要是自己算一下还是比较麻烦，有专门的网站可以查这个计算
http://www.samformat.info/sam-format-flag
- 第三列：RNAM，reference sequence name，实际上就是比对到参考序列上的染色体号。若是无法比对，则是*；
- 第四列：position，read比对到参考序列上，第一个碱基所在的位置。若是无法比对，则是0；
- 第五列：Mapping quality，比对的质量分数，越高说明该read比对到参考基因组上的位置越唯一；(数字越高越好）
- 第六列：CIGAR值，read比对的具体情况，（这个很重要）
“M”表示 match或 mismatch； “I”表示 insert； “D”表示 deletion； “N”表示 skipped（跳过这段区域）； “S”表示 soft clipping（被剪切的序列存在于序列中）； “H”表示 hard clipping（被剪切的序列不存在于序列中）； “P”表示 padding； “=”表示 match； “X”表示 mismatch（错配，位置是一一对应的）；
- 第七列：MRNM(chr)，mate的reference sequence name，实际上就是mate比对到的染色体号，若是没有mate，则是*；
- 第八列：mate position，mate比对到参考序列上的第一个碱基位置，若无mate,则为0；
- 第九列：ISIZE，Inferred fragment size.详见Illumina中paired end sequencing 和 mate pair sequencing，是负数，推测应该是两条read之间的间隔(待查证)，若无mate则为0；
- 第十列：Sequence，就是read的碱基序列，如果是比对到互补链上则是reverse completed eg. CGTTTCTGTGGGTGATGGGCCTGAGGGGCGTTCTCN
- 第十一列：ASCII，read质量的ASCII编码。
- 第十二列之后：Optional fields，可选的自定义区域
AS:i 匹配的得分
XS:i 第二好的匹配的得分
YS:i mate 序列匹配的得分
XN:i 在参考序列上模糊碱基的个数
XM:i 错配的个数【Cui add：如果有错配，这里的数字非0】
XO:i gap open的个数
XG:i gap 延伸的个数
NM:i 经过编辑的序列
YF:i 说明为什么这个序列被过滤的字符串
MD:Z 代表序列和参考序列错配的字符串(例如MD:Z:45A2C3 失配碱基的位点，45,45+2两个位点失配)
XT:A:U read只有一个完整比对，U unique【Cui add：只有一个我完整比对】
XT:A:R read有一个以上位置完整比对，R repeat
NM:i:2 read有2个碱基mismatch
X0:i:2 read有2个位置完整比对(与XT:A有对应关系)
X1:i:2 read有2个位置以1个碱基失配比对
XA:Z:Chr3,+1530, 50M,0;Chr4,-7568,50M,1 有0/1个碱基失配的详细比对情况(与XT:A、X0:i、X1:i有对应关系)
BAM文件的一些操作
- 选取任一区域，提取BAM文件
samtools view
- 计算基因组每一个位置的覆盖度
samtools depth
- 计算基因组每一个位置的详细覆盖情况（包据突变信息等）
它的主要功能主要是生成BCF、VCF文件或者pileup一个或多个bam文件。比对记录以在@RG中的样本名作为区分标识符。如果样本标识符缺失，那么每一个输入文件则视为一个样本。
samtools mpileup
#第一步：把sam文件转换成bam文件,我们得到map.bam文件
samtools view -bS map.sam > map.bam;
#第二步：sort 一下 BAM 文件,得到map.sorted.bam
samtools sort map.bam map.sorted;
#第三步：创建一个关于bam的索引文件,我们得到一个map.sorted.bam.bai的文件
samtools index map.sorted.bam;
#第四步：找snp，这里用的是sort以后的bam文件，如果不是，就会不断的报错
samtools mpileup -t DP -t SP -uvf ref.fa map.sorted.bam --output map.sorted.bam.vcf
samtools mpileup -ugf ref.fa map.sorted.bam | bcftools view -vcg -D100 ->snp.vcf
如果我们要获取全部的位点的信息，而不是仅仅snp位点，那么我们只需要把最后一行的-v (bcftools) 去掉就可以了，如下：
samtools mpileup -ugf ref.fa map.sorted.bam | bcftools view -cg -D100 ->snp.vcf
- 快速获取各条染色体的比对情况
samtools idxstats
