一、首先需要知道以下几个知识点:

1.1-based coordinate system

A coordinate system where the rst base of a sequence is one. In this coordinate

system, a region is specied by a closed interval. For example, the region between the 3rd

and the 7th bases inclusive is [3; 7]. The SAM, VCF, GFF and Wiggle formats are using the 1-based coordinate system.

2.0-based coordinate system

A coordinate system where the rst base of a sequence is zero. In this

coordinate system, a region is specied by a half-closed-half-open interval. For example, the region

between the 3rd and the 7th bases inclusive is [2; 7). The BAM, BCFv2, BED, and PSL formats are using the 0-based coordinate system.

3.模板(Template):

由测序仪测序所得或由原始序列组装所得的DNA/RNA序列

4.片段(Segment)

一段连续的序列或者子序列

5.Read

一段由测序仪测序所得的原始序列。一条Read可能由多个片段组成,在测序数据中,reads是根据它们被测的顺序来建立索引的。

6.Linear alignment(线性比对)

一个Read单向地比对到参考基因组上,这个比对结果中可以有插入、缺失、跳跃等,但是不能存在“双向”的比对结果,即Read的一段比对到正链参考基因组、一段匹配到负链,这种方向切换是不允许的,在SAM文件中,线性比对的特性就是:只用一行来记录。

7.Chimeric alignment(嵌合比对)

就是当一条Read对比时,比对到了多个区域,但是这些区域并没有重叠的部分,也即由多个“线性比对”结果组成了一个集合,这个集合就组成了一个嵌合比对,嵌合比对中只有一个“线性比对”结果是具有代表性的,其余的都以补充的身份出现,嵌合比对的特征就是多个“线性比对”记录中的Read对应的Qname(Read的名字,每个Read只有一个Qname)都是相同的,且这些“线性比对”集合中的每个记录的flag值都是一样的。

8.Read alignment(Read 比对)

无论是上面提到的线性比对还是嵌合比对,只要能够完整的表现出一条Read的对比情况,就是一个Read 比对。

9.Multiple mapping(多次比对)

由于序列的重复性,导致一个Read在比对时会被比对到多个区域上,其中只有一个比对质量最好的会被当做比对结果的代表性结果,目前来看,这种决定方式还不是很严谨。

多次比对和嵌合比对是有根本性区别的,多次比对是因为序列本身具有重复性引起的,属于正常比对结果,而嵌合比对更多是由于实验、测序、结构突变、融合、以及其他因素引起的。

10.Phred scale:

一个可能性区间,如果一个碱基被检测正确的正确率是99%,那么错误的概率就是1%=0.01=10(^-2),此时phredscore=-10log10(0.01)=20,

Phred score

二、SAM文件内容解析

SAM文件由两部分组成,头部和主体,都以tab分列。

头部内容主要以各种说明为主,比如说明所用软件啦,参考基因组信息啦,排序信息啦等等,下面表格是SAM文件中涉及的一些专有名词的解释。

名词定义

2.1头部内容:

头部内容说明信息均是以@符合开头的,在头部是以“@说明类型码TAB键TAG:Value” 开始的,每个标签和说明类型码都是由两个字母组成的,下面将列举说明类型码和TAG

2.1.1说明类型码有HD、SQ、RG、PG、CO五种;

2.1.1.1 HD(此文件的数据信息)

HD:代表意思为:SAM文件的开头标志、一般只要出现就会在第一行

HD中的标签(TAG)有:

VN*:注释版本信息

SO:比对结果的排序类型:有unknown、unsorted、queryname、coordinate四种排序类型

GO:比对结果的分组信息:相似的比对结果会被分到一组,这里的分组结果中并不是需要全部进行过排序的,排序的类型有:none (default)、query (alignments are grouped by QNAME)、reference (alignments are grouped by RNAME/POS  三种类型

SS:子排序类型,比如在某次算法中需要根据coordinate排序,而在每个coordinate排序的结果中又根据QNAME进行了排序,则在SAM文件中就应该表示为:@HD SO:coordinate SS:coordinate:queryname.如果在SO的基础排序中,排序的类型不在之前定义的四种排序类型中时,则SO对应的就应该时unsorted,此时SS就会起主要作用,比如,如果基础排序是根据一个辅助标签MI排序的,之后又根据coordinate排序的,则在SAM头部中就应该表现为@HD SO:unsorted SS:unsorted:MI:coordinate.

2.1.1.2 SQ(参考序列信息)

SQ:说明类型码代表的含义为参考序列字典,SQ的排序决定了比对结果的排序顺序

SQ中又以下标签:

1.SN*:参考序列名称

2.LN*:参考序列长度

3.AH、AN、AS、DS、M5、SP(种族)、TP、UR,这些不常用

例如:@SQSN:JF-PLAC8_CT_convertedLN:88

2.1.1.3 RG:样本信息

Read Group,每个Sample都有一个RG ID,一个Sample可以在多个库中进行测序。

RG对应的有以下TAG:

ID*:Read group的唯一ID

BC:辨别样本或文库的标签序列

SM:样品名

LB:文库名

PU:测序仪

PL:测序平台

CN:产生read的测序中心的名称

2.1.1.4 PG (运行程序信息)

PG:程序 说明类型码,PG中有以下几个标签:

ID*:程序记录标识,每个程序记录都只有一个ID

PN:程序名称

CL:命令行内容(utf-8编码)

PP:好像是指前一个 @PG-ID:不怎么用这个,有了解的大哥大姐可以帮忙做个备注

DS:描述

VN:所用程序版本

2.1.1.5 CO (备注或评论信息)

以上就是五个说明类型码代表的意思以及常用的标签意义。下面截图作为例子,可以对照上面的内容看一下:

SAM文件头部信息举例

2.2 主体内容

SAM文件主体内容中主要有11列内容,这11列中的内容就如下表:

当数据中这11列内容哪个对应内容有缺失或者无效的话就用‘0' 或 `*'代替。

1.QNAME:Read ID,就简单理解成Read的名称,其中会包含一些测序平台信息

2.FALG:可以理解为比对结果的标志,可以根据FLAG值筛选比对结果

FLAG的值代表含义如下表(本图来自http://www.360doc.com/content/18/0329/12/28491187_741231514.shtml,感谢作者,如有侵权请联系本人删除)

FLAG值对应表

FLAG值若记不住可以直接在中进行查询。在SAM文件中出现的flag值是涉及到的value相加得到的值,比如99,177等,可以在上述网站查询

3.RNAME:理解成比对上的染色体号即可

4.POS:比对到参考序列上的位置

5.MAPQ:比对质量值,算法与文中第一部分Phred score算法一致:-10log10(错误率),若值为255表示其比对结果不可用,如果是unmapped read则MAPQ为0

6.CIGAR(Compact Idiosyncratic Gapped Alignment Report)简要描述比对结果(来自http://www.360doc.com/content/18/0329/12/28491187_741231514.shtml,感谢作者,如有侵权请联系本人删除)

简要比对信息表达式(Compact Idiosyncratic Gapped Alignment Report),其以参考序列为基础,使用数字加字母表示比对结果,比如3S6M1P1I4M,前三个碱基被剪切去除了,然后6个比对上了,然后打开了一 个缺口,有一个碱基插入,最后是4个比对上了,是按照顺序的;

M”表示 match或 mismatch;

“I”表示 insertion;

“D”表示 deletion(表示的是READ和ref相匹配时,参考基因组中需要deletion的部分,不是READ);

“N”表示 skipped(跳过这段区域);

“S”表示 soft clipping(被剪切的序列存在于序列中);

“H”表示 hard clipping(被剪切的序列不存在于序列中,已经在除去低质量READ的时候被过滤了);

“P”表示 padding(填充);比如参考基因组序列本来为ATTACGAC,read序列为ATAATACGAC,那么,如果在比对的sam结果文件中,有两种比对方式,如下:

1.如果reference是padded reference,也即参考基因组展示为AT**TACGAC,那么**就是代表了ref考虑了read序列的插入情况,但是这种表示情况下CIGAR中就不会再出现p和I(Insertion)了,因为参考基因组中已经考虑了插入的情况了。而如果有其他read跟这个padded ref比对时,对于ref中pad的区域没有序列的话,就会以D(deletion)来表示这个read。

2.另一种情况就是unpadded reference,也即reference是正常的参考基因组(可参考https://blog.csdn.net/xubo245/article/details/51283022),一般情况下,以这一种为主。

“=”表示 match;

“X”表示 mismatch(错配,位置是一一对应的);

7.RNEXT:mate序列匹配上的染色体号(比如,Read2对应的染色体号),如果第三列和这列都是“*”则说明此Read没有匹配到任何一个染色体,如果第三列有信息,这列是“”或者“=”则代表此Read匹配到第三列对应的染色体上。

8.PNEXT:该列表示与该reads对应的mate pair reads的比对位置,如果这对pair-end reads比对到同一条reference序列上,在sam文件中reads的id出现2次,Read1比对的第4列等于Read2比对的第8列。同样Read1比对的第8列等于Read2比对的第4列。例如:

第1列(Read id)····第4列(Read1比对位置)····第8列(mate-pair reads比对位置)

22699:1759····124057649····124057667

22699:1759····124057667····124057649

相同的reads id一个来自Read1文件,一个来自Read2文件,第4列和第8列是对应的

9.TLEN:signed observed Template LENgth (可以理解为文库插入片段长度)

如果R1端的read和R2端的read能够mapping到同一条Reference序列上(即第三列RNAME相同),则该列的值表示第8列减去第4列加上第6列的值,R1端和R2端相同id的reads其第九列值相同,但该值为一正一负,R1文件的reads和R2文件的reads,相同id的reads要相对来看。在进行该第列值的计算时,如果取第6列的数值,一定要取出现M的值,S或H的值不能取。

(8.9解释来自:作者:genome_denovo,来源:CSDN,原文:https://blog.csdn.net/genome_denovo/article/details/78712972,感谢作者,如有侵权请联系本人删除)

10.序列信息,

11.Read质量信息(ASCII编码)

12.可选区域

举几个例子:NM,MD

NM可以简单的理解为:如果要把READ变为跟reference一致需要几步(一般步骤有三种类型,碱基替换,删除碱基,添加碱基)

MD:简单的理解为匹配结果(跟第10列结合可以看出reference的原序列),下面举个例子:

从上面的图可以看出来这里的NM值为2,MD为3^A42T5,由MD可以知道这个read比对结果中前3个匹配上了,之后需要在read中插入一个A(也就是read相比ref少了一个A),再之后42个碱基是匹配上了,到然后这42个之后的碱基在参考基因组上应该是T,后面5个是匹配上了。这样的话我们根据read的序列就能推断出ref的序列为:CCGATAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCC。

我们检查一下推测的是否正确,来看hg19中这段序列:

hg19中的序列信息

可以发现除了开头四个碱基在hg19中是N以外,其他的是匹配一致的(好像是因为有写序列是在端粒中的,所以呈现的序列信息就是N,也有些是因为还没有测出来,所以用N代替。)、

这样再来看的话,就可以知道read相比ref有两处不同,一个是缺了个A,一个是碱基错配(T-A),所以如果read要和ref一致的话需要编辑两次,一次是添加一个A,一次是把A变成T。

如果觉得有用,就帮我点个赞吧^_^!

bam文件读取_SAM/BAM 格式文件内容解析相关推荐

  1. CDays-3 习题二 (字典及文件读取练习)及相关内容解析。Python 基础教程

    读取某一简单索引文件cdays-3-test.txt,其每行格式为文档序号 关键词,现需根据这些信息转化为倒排索引,即统计关键词在哪些文档中,格式如下:包含该关键词的文档数 关键词 => 文档序 ...

  2. 使用Python读取LabVIEW TDMS 格式文件转成Excel格式+多进程版本

    使用Python读取LabVIEW TDMS 格式文件转成Excel格式+多进程版本 文章目录 使用Python读取LabVIEW TDMS 格式文件转成Excel格式+多进程版本 前言: 背景 tm ...

  3. java将excel文件转换成txt格式文件

    在实际应用中,我们难免会遇到解析excel文件入库事情,有时候为了方便,需要将excel文件转成txt格式文件.下面代码里面提供对xls.xlsx两种格式的excel文件解析,并写入到一个新的txt文 ...

  4. java使用jxl包读写excel表格文件,即xls格式文件

    全栈工程师开发手册 (作者:栾鹏) java教程全解 java使用jxl包读写excel表格文件,即xls格式文件 本实例演示使用jxl包实现对excel文件的操作,下载 测试代码 public st ...

  5. 【ASE+python】实现将poscar格式文件批量转换为xsd格式文件

    将poscar格式文件批量转换为xsd格式文件 ASE介绍 ASE安装 ASE的ase.io.read()与ase.io.write() ase.io.read() ase.io.write() 单份 ...

  6. AAC文件解码成PCM格式文件

    上一篇写到PCM格式文件编码成AAC格式文件,这一步的原因是有利于传输.可以将PCM文件做了很大的压缩力度,使得包变得更小,便于传输.我使用播放器播放了AAC文件听到的是音速明显变快了,声音也变得尖锐 ...

  7. 如何在Photoshop中载入使用pat格式的文件?ps图案pat格式文件载入教程

    PS中有一种为"pat"后缀的ps图案模式,小编今天为大家带来了如何在Photoshop中载入使用pat格式的文件?ps图案pat格式文件载入教程,有需要的小伙伴快来看看吧! 在桌 ...

  8. python打包zip文件_python 解压文件,合并文件 打包成zip格式文件 生成MD5值

    #!/usr/bin/env python #_*_encoding:utf-8 # 2018/05/29 #augustyang #2.0 ''' 解压文件,合并文件 打包成zip格式文件 生成MD ...

  9. 如何将字库生成工具生成的 .DZK/ . bin格式的文件转成.c格式文件

    如何将字库生成工具生成的 .DZK/ . bin格式的文件转成.c格式文件 在我们项目开发的过程中,特别是做界面显示的时候,经常会遇到如 :多国文字点阵字库生成器TS3等软件其生成的.DZK格式文件, ...

最新文章

  1. tensorflow1.14.0  包含了1.x和2.x内容,此后版本要求兼容该版本
  2. 返璞归真 asp.net mvc (2) - 路由(System.Web.Routing)
  3. 套接字编程(VC_Win32)
  4. 面试官问:你讲讲分布式事务问题的几种方案?
  5. C# ckeditor+ckfinder的图片上传配置
  6. mysql 140824,Oracle 12c创建可插拔数据库(PDB)及用户
  7. memcache如何更新mysql_使用MySQL触发器如何实现memcache自动更新
  8. 第一次收到这么用心的感谢信
  9. delphi 生成 超大量xml_用OpenCV4实现图像的超分别率
  10. 洛谷P1217回文质数
  11. html背景图片自适应屏幕
  12. Win11系统开启控制面板会闪退怎么解决?
  13. 三层交换机设置成路由
  14. 20230220学习总结02
  15. python预测子女身高_Python 孩子身高预测
  16. Excel 轻松制作 二级联动 下拉列表清单
  17. CentOS7下安装jmeter5.3
  18. 齐供应TAPPI四碘化5,10,15,20-四(对-N,N,N三甲基苯胺基)卟啉敏化的钛酸盐纳米管(TAPPI-TNTs)高效的可见光催化剂岳
  19. linux系统怎么重启网卡?linux重启网卡的三种教程
  20. Linux tar/rpm/yum命令软件安装

热门文章

  1. CDA LEVEL 1 考试,知识点《机器学习基本概念》
  2. Python之ARP协议探测MAC地址
  3. html body language,关于BODY LANGUAGE
  4. Unhandled exception:java.lang.IllegalAccessException提示报错
  5. 互联网自动化赚钱的方法
  6. 如何添加RichFaces 3.3.x到Maven项目
  7. 〖Python 数据库开发实战 - Python与MySQL交互篇④〗- 数据库连接池技术
  8. Android 页面跳转时发生双击导致app闪退的解决方案
  9. 什么是意志力?如何提高意志力?
  10. 关于大三学生的请教回复