Difference between Hard Clip(H) and Soft Clip(S) in Samtools CIGAR string

一般人都知道 H 和 S 的表面上的区别,即 S 就是 soft, H 就是 hard,S 后,序列里还是会保留序列的信息,而 H 则不会。

但这只是表面上的,在深层次的意义上, H 和 S 又有什么本质的不同呢?

首先要了解嵌合体的概念:

嵌合体就是两个不同的序列错误的拼接到了一起,也就是一条序列分别比对到了 ref 的两个地方(这和多重比对、次级比对之间又有区别)

Example of extended CIGAR and the pileup output.

(a) Alignments of one pair of reads and three single-end reads.

(b) The corresponding SAM file. The ‘@SQ’ line in the header section gives the order of reference sequences. Notably, r001 is the name of a read pair. According to FLAG 163 (=1 + 2 + 32 + 128), the read mapped to position 7 is the second read in the pair (128) and regarded as properly paired (1 + 2); its mate is mapped to 37 on the reverse strand (32). Read r002 has three soft-clipped (unaligned) bases. The coordinate shown in SAM is the position of the first aligned base. The CIGAR string for this alignment contains a P (padding) operation which correctly aligns the inserted sequences. Padding operations can be absent when an aligner does not support multiple sequence alignment. The last six bases of read r003 map to position 9, and the first five to position 29 on the reverse strand. The hard clipping operation H indicates that the clipped sequence is not present in the sequence field. The NM tag gives the number of mismatches. Read r004 is aligned across an intron, indicated by the N operation.

(c) Simplified pileup output by SAMtools. Each line consists of reference name, sorted coordinate, reference base, the number of reads covering the position and read bases. In the fifth field, a dot or a comma denotes a base identical to the reference; a dot or a capital letter denotes a base from a read mapped on the forward strand, while a comma or a lowercase letter on the reverse strand.



clipped alignment因为着在比对过程中,并没有用到全部的read的序列,read两段的序列被截取了(clip or trim)。如下表示,即为clip alignment。

Alignment:
Read:          ACGGTTGCGTTAA-TCCGCCACG
|                           ||||||||| ||||||
Reference: TAACTTGCGTTAAATCCGCCTGG

与clipped alignment对应的是spliced alignment,即read的中间没有比对到而两段比对上了。对应的表示如下:

Alignment:
Read:          ACGGTTGCGTTAAGCTCATCCGCCACG
|                 |||||||||||||         |||||||||
Reference: ACGGTTGCGTTAA…..TCCGCCACG

clip alignment对应的CIGAR表示有两种S (soft clip) 和H (hard clip)。
BWA提到If the read has a chimeric alignment, the paired or the top hit uses soft clipping and is marked with neither 0x800 nor 0x100 bits. All the other hits part of the chimeric alignment will use hard clipping and be marked with 0x800 if option “-M” is not in use, or marked with 0x100 otherwise.

即如果发现嵌合比对,最好的比对top hit标记为soft clipping,其余的则标记为hard clipping。

如果是hard clip,则截取的部分不会在SAM文件对应的read中出现 (clipped sequences not present in SEQ),如果是soft clip (clipped sequences present in SEQ),则会出现。

Understand?

Ref:https://github.com/lh3/bwa/blob/master/NEWS.md

转自:http://wp.zxzyl.com/?p=131


理解1:

Hard masked bases do not appear in the SEQ string, soft masked bases do.

So, if your cigar is: `10H10M10H` then the SEQ will only be 10 bases long.

if your cigar is 10S10M10S then the SEQ and base-quals will be 30 bases long.

首先,结果展示方式有区别:?比如说10H10M10H,第10列的碱基序列只显示10bp;而如果是10S10M10S的话,就会显示30bp的序列,尽管开头和结尾的20bp也没比上。

In the case of soft-masking, even though the SEQ is present, it is not used by variant callers and not displayed when you view your data in a viewer. In either case, masked bases should not be used in calculating coverage.

在soft中,即使显示的序列比hard的要长,但是计算变异或可视化比对结果时,这些序列也不会被考虑。而且,2种情况计算覆盖度?时,mask的碱基都不会考虑。

例子:

20692128    97    viral_genome    21417    60    69M32S    chr7    101141242    0    TACATCTTCTCCCTCTCTCACGACACAAGAATTAGTCACATAGGGATGTTCTCGTAAATCTACATTATCTTACAAAAACATTTTTTAAAAATTTGCTAGGT (101bp)    GGGGGGGGGGGGGGEGGEGGGGGGGGGFGGGGGGGGG[email protected]    NM:i:4    MD:Z:6G34G6C5C14    AS:i:49    XS:i:0    SA:Z:chr7,101141091,+,66S35M,60,0;
20692128    353    chr7    101141091    60    66H35M    =    101141242    252    ATCTTACAAAAACATTTTTTAAAAATTTGCTAGGT  (35bp)[email protected]    NM:i:0    MD:Z:35    AS:i:35    XS:i:23    SA:Z:gi|224020395|ref|NC_001664.2|,21417,+,69M32S,60,4;
20692128    145    chr7    101141242    60    101M    gi|224020395|ref|NC_001664.2|    21417    0    GCAACAGAGCGAGACCCTATATTCATGAGTGTTGCAATGAGCCAAGTAGTGGAGGTTGGCTTTTGAAGGCAGAAAAGGACTGAGAAAAGCTAACACAGAGA    FEGCGGGGGCGEFCDEEEEGGGGGGGGGGGGGGGEGGGGGGFGGGEGGG

?理解2:

当同一条reads比对到不同chr时(嵌合reads),会以hard clip的显示显示。比如上面的例子,R1分别比到了viral基因组和chr7上(前面69bp比到viral,后面35比到chr7),R2比到了chr7。

sam格式详细说明见:http://davetang.org/wiki/tiki-index.php?page=SAM

时间: 2024-10-11 01:08:50

Difference between Hard Clip(H) and Soft Clip(S) in Samtools CIGAR string的相关文章

(转)开源项目miaosha(下)

石墨文档:https://shimo.im/docs/2XlwliBQAYsKCHbq/ (二期)20.开源秒杀项目miaosha解读(下) [课程20]jmeter.xmind81.5KB [课程20]秒杀项目杂谈.xmind0.2MB 代码片段: RedisTool.java3.2KB MiaoshaInterface.java4.7KB 测试过程中,发现了一个bug. wang.moshu.cache.MiaoshaSuccessTokenCache#validateTokenByKey(

基于RTP的h.264视频传输系统设计(一)

一.H.264 的层次介绍 H.264 定义三个层次,每个层次支持一组特定的编码功能,并且依照各个层次指定所指定的功能.基础层次(baselineprofile)支持I 帧和 P 帧[1]的帧内和帧间编码,支持自适应的可变长度的熵编码(CAVLC).主要层次(main profile)支持隔行扫描视频,B帧[2]的帧内编码,使用加权预测的帧内编码和使用上下文的算术编码(CABAV).扩展层次(extendedprofile)不支持隔行扫描视频和CABAC,但增加了码流之间高效的转化模式(SP 和

【C/S通信交互之Http篇】Cocos2dx(Client)使用Curl与Jetty(Server)实现手机网游Http通信框架(内含解决curl.h头文件找不到问题)

之前已经分享过一篇基于Cocos2dx与服务器使用Socket进行通信的框架,还不太熟悉的请移步到如下博文中: [C/S通信交互之Socket篇]Cocos2dx(Client)使用BSD Socket与Mina(Server)手机网游通信框架! 那么今天Himi来分享如何在cocos2dx中使用Http来访问Server端并且获取数据: 这里对于Server端,Himi选用,Jetty,对于Jetty不太熟悉的可以先自行baidu-google-是个servlet的容器.类似JSP. 什么是s

网站被挂暗链、点开同一链接进入不同页面(博彩页面)、恶意脚本(INCLUDE(pack('H*')……)之类

原文链接 论坛被挂暗链问题分析与解决http://blog.kankanan.com/posts/2014/04/01_8bba575b88ab6302669794fe95ee9898520667904e0e89e351b3.html 发现问题 有网友反映我们的论坛被挂了暗链,具体表现为从 google 搜索论坛名称结果如下图所示: 直接搜索论坛网址出现的一些热门帖子也被挂了暗链,通过 google 搜索结果访问会跳到恶意网站, 解决问题 直接通过网址访问论坛则没有任务问题,应该是论坛被注入了恶

黄聪:FFmpeg视频转码技巧之-crf参数(H.264篇)

昨天,有个朋友给我出了个难题:他手上有一个视频,1080P的,49秒,200多兆:要求在确保质量的情况下把文件压缩到10M以内. 这是什么概念呢?按照文件大小10M来计算,码率是:10 x 8 / 49 = 1.6 Mbps.也就比VCD的质量略好一点(注:VCD的标准码率是1150 Kbps).谈何“确保质量”?mission impossible啊! 咱还是现实一点吧.在不明显损失画质的前提下,看看使用FFmpeg能够帮到多少忙.用iPhone拍了一个1920 x 1080的视频,33秒,4

C语言中.h和.c文件解析(转载)

转载:http://www.cnblogs.com/laojie4321/archive/2012/03/30/2425015.html   简单的说其实要理解C文件与头文件(即.h)有什么不同之处,首先需要弄明白编译器的工作过程,一般说来编译器会做以下几个过程:       1.预处理阶段 2.词法与语法分析阶段 3.编译阶段,首先编译成纯汇编语句,再将之汇编成跟CPU相关的二进制码,生成各个目标文件 (.obj文件) 4.连接阶段,将各个目标文件中的各段代码进行绝对地址定位,生成跟特定平台相

关于PROFIBUS Master(H)不能正确识别并处理 DP-Slave 回复的RS帧的一些思考

图1.是在测试过程中,发现PROFIBUS Master(H)不能正确识别并处理 DP-Slave 回复的RS帧,引起Slave回复 RS 帧的操作是"断开Slave与Master之间的PROFIBUS电缆,然后再恢复物理链路",但是Master却一直在 请求Data,而Slave却回复RS通知Master"服务不匹配".图1.中分析了这种RS交互过程,并提出了可能的解决办法. 图1 图2.是在测试过程中,PROFIBUS Master(A)正确识别并处理 DP-S

C++中头文件(.h)和源文件(.cpp)都应该写些什么

头文件(.h): 写类的声明(包括类里面的成员和方法的声明).函数原型.#define常数等,但一般来说不写出具体的实现. 在写头文件时需要注意,在开头和结尾处必须按照如下样式加上预编译语句(如下): #ifndef CIRCLE_H #define CIRCLE_H //你的代码写在这里 #endif 这样做是为了防止重复编译,不这样做就有可能出错. 至于CIRCLE_H这个名字实际上是无所谓的,你叫什么都行,只要符合规范都行.原则上来说,非常建议把它写成这种形式,因为比较容易和头文件的名字对

MVC5+(H+)+EF6+Linq+Cache 通用管理系统,项目包括(自定义表单,工作流,代码生成器,数据字典,权限管理)【1】

经过1年多的沉淀.新的项目即将发布上线. MVC5+(H+)+EF6+Linq+Cache 通用管理系统,项目包括(自定义表单,工作流,代码生成器,数据字典,权限管理)