samtools及bam文件的相关知识

samtools及bam文件的相关知识

如果flag值是0，那么说明测序为单端测序且这条read是primary line，一般是该read的最佳比对。
这里说一下secondary alignment和supplementary alignment的区别，secondary alignment代表的是该条read比对到多个位置，该条记录为次优比对，在双端测序中，代表hardclip，而supplementary alignment标记为该条比对记录为该read的补充比对记录，关于补充比对记录其实就是上文提及的嵌合比对，
这些记录是相互的补充，有可能只是这条read的一部分，我在二代测序的结果中很容易找到很多的hardclip记录，但是却没有找到supplementary alignment的记录，在三代的比对结果中就比较多了在三代测序的比对结果中常见2048的supplementary alignment。
使用samtools view -f/-F参数可以获取包含/不包含指定flag的内容的条目。
https://blog.csdn.net/qq_35696312/article/details/101760397

可以看到bam文件的header有如下几个部分组成

@HD 标题行(The header line)主要用于版本号，是否排序等标注
@SQ 参考序列字典(Reference sequence dictionary)注意它们的顺序决定了排序的优先级
@RQ read group主要记录了一文件的描述信息，例如样本信息，文库编号等
@PG program主要记录比对的软件等代码信息
@CO comments一行评论
每个字典都有一些详细的内容，具体可以参考开篇的内容，这里就不展开复述了。
原文链接：https://blog.csdn.net/qq_35696312/article/details/101760397
将sam文件转换成bam文件 $ samtools view -bS abc.sam > abc.bam $ samtools view -b -S abc.sam -o abc.bam 提取比对到参考序列上的比对结果 $ samtools view -bF 4 abc.bam > abc.F.bam 提取paired reads中两条reads都比对到参考序列上的比对结果，只需要把两个4+8的值12作为过滤参数即可 $ samtools view -bF 12 abc.bam > abc.F12.bam 提取没有比对到参考序列上的比对结果 $ samtools view -bf 4 abc.bam > abc.f.bam 提取bam文件中比对到caffold1上的比对结果，并保存到sam文件格式 $ samtools view abc.bam scaffold1 > scaffold1.sam 提取scaffold1上能比对到30k到100k区域的比对结果 $ samtools view abc.bam scaffold1:30000-100000 $gt; scaffold1_30k-100k.sam 根据fasta文件，将 header 加入到 sam 或 bam 文件中 $ samtools view -T genome.fasta -h scaffold1.sam > scaffold1.h.sam

8. 其它有用的命令
   reheader 替换bam文件的头

$ samtools reheader
cat 连接多个bam文件，适用于非sorted的bam文件

$ samtools cat [-h header.sam] [-o out.bam] [ … ]
idxstats 统计一个表格，4列，分别为”序列名，序列长度，比对上的reads数，unmapped reads number”。最后一排则显示没有比对到任何一条序列的reads number。

$ samtools idxstats

sam中的flag的详细解释：https://www.cnblogs.com/xudongliang/p/5437850.html

1 ： 代表这个序列采用的是PE双端测序

2： 代表这个序列和参考序列完全匹配，没有插入缺失

4： 代表这个序列没有mapping到参考序列上

8： 代表这个序列的另一端序列没有比对到参考序列上，比如这条序列是R1,它对应的R2端序列没有比对到参考序列上

16：代表这个序列比对到参考序列的负链上

32 ：代表这个序列对应的另一端序列比对到参考序列的负链上

64 ： 代表这个序列是R1端序列， read1;

128 : 代表这个序列是R2端序列，read2；

256： 代表这个序列不是主要的比对，一条序列可能比对到参考序列的多个位置，只有一个是首要的比对位置，其他都是次要的

512： 代表这个序列在QC时失败了，被过滤不掉了（# 这个标签不常用）

1024: 代表这个序列是PCR重复序列（#这个标签不常用）

2048: 代表这个序列是补充的比对（#这个标签具体什么意思，没搞清楚，但是不常用）
比如我们可以通过samtools view中的-f参数，找到flag值包括2的比对结果

作者：一只烟酒僧
链接：https://www.jianshu.com/p/afb51fc66658
samtools的使用技巧（持续更新）

view

samtools view [options] in.sam|in.bam|in.cram [region…]

-f 提取 ## -f 4 提取出没有mapping上的reads

-F 过滤 ## -F 4 过滤掉没有mapping上的reads，也就是说提取出mapping上的reads

-u 输出格式为未压缩的bam格式

-q 过滤掉MAPQ值低某个阈值 ## -q 1 过滤掉MAPQ值低于1的情况

-h 设定输出的SAM文件带有header

-b 输出格式设定为BAM

-S 输入格式为SAM

提取比对到参考序列的结果:

samtools view -bF 4 tmp.bam > tmp_F.bam

提取双端序列都比对到参考序列（4+8）的结果：

samtools view -bF 12 tmp.bam > tmp_F.bam

提取比对到chr1的结果

samtools view -b tmp.bam chr1 > tmp_chr1.bam

注：With no options or regions specified, prints all alignments in the specified input alignment file (in SAM, BAM, or CRAM format) to standard output in SAM format (with no header)，也就是说，没有设定输出格式的话，默认是输出SAM格式，并且是没有header的SAM

链接：https://www.jianshu.com/p/523cc26ffaa1