SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
bitwise flag in sam format and others aurelielaugraud Bioinformatics 48 11-24-2012 11:01 PM
BWA output bitwise flag for mapped/unmapped reads wenhuang Bioinformatics 1 08-29-2011 03:54 PM
Flag=4 in SAM Rachelly Bioinformatics 2 12-22-2010 02:54 AM
Aligner that outputs H2 flag in SAM file phalaenopsis Bioinformatics 0 08-04-2010 05:14 AM
SAM file flag problem ptong7 Bioinformatics 4 07-30-2009 03:32 AM

Reply
 
Thread Tools
Old 06-02-2010, 06:04 PM   #1
win804
Member
 
Location: Singapore

Join Date: Apr 2010
Posts: 18
Default bitwise flag in SAM file

Hi Everyone,
Recently I have aligned my sequencing results using BWA and get the aligned result in SAM format.
My sequence is a single end read and therefore I would expect the bitwise flag (second field) in the SAM format contains only 0, 4, 16 (correct me if i m wrong).

However, I notice one aligned result with strange value:
id1 20 chrX 154913749 25 36M * 0 0 TGCGGTCCAACCCTAACCCTAACCCTAACCCTAACC SVVVVV[[[[[[[[[[Z[Z[[[[ZZ[[[[[Z[[[[Z XT:A:U NM:i:2X0:i:1 X1:i:0 XM:i:2 XO:i:0 XG:i:0 MD:Z:2T4T28

I cannot explain why there is a bitwise flag 20 in the second field. I expect it to come from 4+16, however, since it is mapped in reverse strand (give 16), then why there is 4 (which means unmapped) contribute to the bitwise flag?

Thank you very much for the help.
win804 is offline   Reply With Quote
Old 06-02-2010, 07:16 PM   #2
nilshomer
Nils Homer
 
nilshomer's Avatar
 
Location: Boston, MA, USA

Join Date: Nov 2008
Posts: 1,285
Default

Quote:
Originally Posted by win804 View Post
Hi Everyone,
Recently I have aligned my sequencing results using BWA and get the aligned result in SAM format.
My sequence is a single end read and therefore I would expect the bitwise flag (second field) in the SAM format contains only 0, 4, 16 (correct me if i m wrong).

However, I notice one aligned result with strange value:
id1 20 chrX 154913749 25 36M * 0 0 TGCGGTCCAACCCTAACCCTAACCCTAACCCTAACC SVVVVV[[[[[[[[[[Z[Z[[[[ZZ[[[[[Z[[[[Z XT:A:U NM:i:2X0:i:1 X1:i:0 XM:i:2 XO:i:0 XG:i:0 MD:Z:2T4T28

I cannot explain why there is a bitwise flag 20 in the second field. I expect it to come from 4+16, however, since it is mapped in reverse strand (give 16), then why there is 4 (which means unmapped) contribute to the bitwise flag?

Thank you very much for the help.
Probably a bug in BWA. Why is there a chr/position and non-zero mapping quality for this read? Also, I "blatted" the base sequence and it did not map to the X chromosome. Anyone else?
nilshomer is offline   Reply With Quote
Old 06-02-2010, 07:55 PM   #3
win804
Member
 
Location: Singapore

Join Date: Apr 2010
Posts: 18
Default

I guess it's a bug also. I modified the sequence according to MD:Z:2T4T28 and blatted against hg18 and I got this hit:
QUERY SCORE START END QSIZE IDENTITY CHRO STRAND START END SPAN
id1 35 1 36 36 100.0% 2 + 181848646 181848793 148

Thanks Nilshomer.
I wonder whether anybody has encountered this kind of problem before.
win804 is offline   Reply With Quote
Old 06-03-2010, 11:30 AM   #4
drio
Senior Member
 
Location: 41°17'49"N / 2°4'42"E

Join Date: Oct 2008
Posts: 323
Default

Let see what Heng saids. If you can, share the reference and the read that is causing the issue so the bug can be reproduced. Actually, do you get the same results when aligning that single read against your genome?
__________________
-drd
drio is offline   Reply With Quote
Old 06-03-2010, 12:34 PM   #5
lh3
Senior Member
 
Location: Boston

Join Date: Feb 2008
Posts: 693
Default

Nothing is wrong. Flag 16 only means the read sequence is on the reverse strand. It may or may not be mapped at all. You see coordinate because it stands out of chrX. The BLAT result is not informative. Note the "148"bp span. Maybe you are mapping RNA-seq reads, but bwa is unable to map reads across exon juctions.
lh3 is offline   Reply With Quote
Old 06-03-2010, 01:15 PM   #6
yh253
Member
 
Location: Ireland

Join Date: Jul 2009
Posts: 16
Default

Hi Heng,

Would you please help to clarify one point that if a read mapped to the junction of two chromosomes, I understand that its Flag will be 4 and I can still see the position, CIGAR and all other tags, but would the MAPQ be set at zero?

I saw some chip-seq single end alignments from BWA in my data with Flag=20. I took it as unmapped, but the MAPQ is not equal zero. Is this legal in terms of sam format?

Thanks in advance!

Yuan
yh253 is offline   Reply With Quote
Old 06-03-2010, 05:04 PM   #7
win804
Member
 
Location: Singapore

Join Date: Apr 2010
Posts: 18
Default

My data is a ChIP-seq data.
To test again, I have extracted the single read and aligned towards hg18.
The result becomes:
id1 4 chrX 154913749 25 36M * 0 0 TGCGGTCCAACCCTAACCCTAACCCTAACCCTAACC SVVVVV[[[[[[[[[[Z[Z[[[[ZZ[[[[[Z[[[[Z XT:A:U NM:i:2 X0:i:1 X1:i:0 XM:i:2 XO:i:0
XG:i:0 MD:Z:2T4T28

So, this time it tells me that it is unmapped (with bitwise flag = 4). But strangely it still contains the coordinate. I noticed that the length of chrX is @SQ SN:chrX LN:154913754, which means the mapping coordinate may be wrong? Since the POS=154913749 and if I add the length of my sequence, it is more than chrX length.
win804 is offline   Reply With Quote
Old 06-04-2010, 01:32 AM   #8
yh253
Member
 
Location: Ireland

Join Date: Jul 2009
Posts: 16
Default

If a read mapped to the junction of two chromosomes, it is shown as unmapped, but you can still see all the tags. Your read may mapped to the end of chrX and the beginning the following chromosome.

Yuan
yh253 is offline   Reply With Quote
Old 06-04-2010, 03:09 AM   #9
win804
Member
 
Location: Singapore

Join Date: Apr 2010
Posts: 18
Default

Hi Yuan,
Thanks for your explanation. Now I understand already.
win804 is offline   Reply With Quote
Old 07-05-2010, 09:56 AM   #10
Xi Wang
Senior Member
 
Location: MDC, Berlin, Germany

Join Date: Oct 2009
Posts: 317
Default

Code:
GA-C_0001:1:1:4471:3450#0 117 hg19_refGene_NR_033322 1558 0 * = 1558 0 NCCAATGTGCCCTGTGTGGTCCCGTTACGTGTCCCATCATGTTTTGTGC B````ba`bb`bbbbbbbbbbbb`bbbbbbbabbbbbbbbbbbbbbbbb  
GA-C_0001:1:1:4471:3450#0 157 hg19_refGene_NR_033322 1558 37 49M = 1558 0 CGCTGAGCATGATGGGGCATGTGCGGGAGCGCCAGGCGGGGCATGTAAC bbaaa^b_abcbbbbObbbbcbbbbbbbbbbbbaabbbbbbbbbbbbbb XT:A:U NM:i:1 SM:i:37 AM:i:0 X0:i:1 X1:i:0 XM:i:1 XO:i:0 XG:i:0 MD:Z:1T47
Two lines for two paired-end reads' mapping results are shown here. It was done by BWA. Please note the flags of the two reads: 117 and 157. The third and the fourth bits (from right) of 117 and 157 are (0,1) and (1,1), respectively. The third bit indicates whether the current read is unmapped and the fourth bit indicates whether the mate is unmapped. So it is clear that the third and fourth bits (0,1) and (1,1) for the mate reads are conflicted with each other. Is this a bug of BWA? Any suggestions? Thanks in advance.
__________________
Xi Wang
Xi Wang is offline   Reply With Quote
Old 07-06-2010, 02:54 PM   #11
misko
Junior Member
 
Location: Toronto

Join Date: Jun 2010
Posts: 7
Default

I am not sure if this is a bug, but it does not look like proper SAM format?
The 157 flag -> mate on positive, query on negative
The 117 flag -> mate on negative, query on negative
Not sure if this helps but i was having trouble parsing SAM output and checking it a while back so i made a small tool to help out with making the SAM output a bit more verbose. If you want i posted the source and binaries to,
http://misko.ca/sam2pretty.html
Hope this helps?

Misko
misko is offline   Reply With Quote
Old 07-06-2010, 05:53 PM   #12
Xi Wang
Senior Member
 
Location: MDC, Berlin, Germany

Join Date: Oct 2009
Posts: 317
Default

A nice tool. But how did you think it's not in a proper SAM format? It is the output of BWA.
__________________
Xi Wang
Xi Wang is offline   Reply With Quote
Old 07-07-2010, 05:21 AM   #13
lh3
Senior Member
 
Location: Boston

Join Date: Feb 2008
Posts: 693
Default

Sometimes bwa flags read1 as unmapped, but for read2, the mate-unmapped flag is not set. This is a bug.
lh3 is offline   Reply With Quote
Old 07-07-2010, 05:26 AM   #14
Xi Wang
Senior Member
 
Location: MDC, Berlin, Germany

Join Date: Oct 2009
Posts: 317
Default

Quote:
Originally Posted by lh3 View Post
Sometimes bwa flags read1 as unmapped, but for read2, the mate-unmapped flag is not set. This is a bug.
Thanks, Heng. Hope next release will solve this:-)
__________________
Xi Wang
Xi Wang is offline   Reply With Quote
Old 07-07-2010, 08:36 AM   #15
Ting Zhang
Junior Member
 
Location: China

Join Date: Jul 2010
Posts: 2
Smile Questions about BWA

Quote:
Originally Posted by lh3 View Post
Sometimes bwa flags read1 as unmapped, but for read2, the mate-unmapped flag is not set. This is a bug.
BWA is quite a good software, but I met some problems when using.
My first question is: as Doc. Heng Li said, “Sometimes BWA flags read1 as unmapped, but for read2, the mate-unmapped flag is not set”, flag may make mistakes. So which field should be checked to see if the read has been mapped? Is that only “NM” or something else?
The second question is: when there are several best hits (X0: i: n, n>1), why are annotation fields given for only one best hit, but not for the rest hits in the “XA” record? How do you choose the fully annotated hit? And why are hits in the “XA” sometimes the same as that hit?
Waiting for the answers eagerly
Ting Zhang is offline   Reply With Quote
Old 07-07-2010, 10:09 AM   #16
lh3
Senior Member
 
Location: Boston

Join Date: Feb 2008
Posts: 693
Default

The flag. It is the mate-unmap bit that can be mislabeled. The unmap bit is always set correctly.
lh3 is offline   Reply With Quote
Old 07-07-2010, 11:13 PM   #17
Ting Zhang
Junior Member
 
Location: China

Join Date: Jul 2010
Posts: 2
Smile

Quote:
Originally Posted by lh3 View Post
The flag. It is the mate-unmap bit that can be mislabeled. The unmap bit is always set correctly.
Ok, but for the other best hits listed in "XA", no flags are given, so how can we judge if they are also mapped?
And sometimes, the first best hit appear in in "XA" again, but sometimes not, why? What's the difference between the first best hit and best hits in "XA"?
An example is as follows, where the "hg19_refGene_NM_001190241" appear in "XA" again:
GA-C_0001:1:1:1519:3751#0 83 hg19_refGene_NM_001190241 2021 0 49M = 1932 -138 NCTTTCTCATAAGAATGAAATCTTGGAAATTGCTCTGGATCAAAAAGGA BB
BBBbb_bb_abbbabbb`bbbbbaabbbbb`bbbbbbbabbabbbbb XT:A:R NM:i:1 SM:i:0 AM:i:0 X0:i:3 X1:i:0 XM:i:1 XO:i:0 XG:i:0 MD:Z:0T48 XA:Z:hg19_refGene_NM_020800,-1801,49M,1;hg
19_refGene_NM_001190241,-2021,49M,1;
Ting Zhang is offline   Reply With Quote
Old 07-09-2010, 03:29 PM   #18
lh3
Senior Member
 
Location: Boston

Join Date: Feb 2008
Posts: 693
Default

If you map to human, 99.999% are mapped. Please read the manual page for the other question.
lh3 is offline   Reply With Quote
Old 01-27-2013, 11:49 PM   #19
edge
Senior Member
 
Location: China

Join Date: Sep 2009
Posts: 199
Default

Hi Yuan,

Do you mind to explain what is the exact meaning of 0, 4, 16 if I map my single-end read against reference transcript by using bowtie?
As I know, 4 is refer to unmapped.
How about 0 and 16 represent?
Thanks.
edge is offline   Reply With Quote
Old 01-28-2013, 12:30 AM   #20
win804
Member
 
Location: Singapore

Join Date: Apr 2010
Posts: 18
Default

0 means map in forward strand.
16 means map in reverse strand.
win804 is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 09:30 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO