Background: I've been using bowtie and samtools for a couple years. In the past 6 months have I converted to outputing SAM directly from bowtie (using -S) instead of piping it through bowtie2sam.pl. I just realized recently that bowtie doesn't output number of possible mappings (e.g. using tag X0 or X1) in the output.
Goal: I would like bowtie's ouput to SAM option -S to have the total number of loci to which a read can map. Specifically, tags X0 and/or X1 would be nice to see.
Question: Is the below proposal a reasonable path to accomplishing the goal? Is there some reason that Bowtie didn't do this by default? Is h.oms appropriate?
Proposed Fix: When I examined bowtie 0.12.8 source code, I found that an object of class Hit (variable is typically named 'h') has a field 'oms' that is used to track the number of mappable loci. In the bowtie output, there is output of h.oms, but not in SAM output.
Bowtie's default output in hit.cpp:
Code:
if(!suppress.test(field++)) { if(firstfield) firstfield = false; else ss << '\t'; ss << h.oms; }
Bowtie's -S SAM output in sam.cpp:
Code:
// SEQ ss << '\t' << h.patSeq; // QUAL ss << '\t' << h.quals; // // Optional fields // // Always output stratum ss << "\tXA:i:" << (int)h.stratum; // Output MD field size_t len = length(h.patSeq); int nm = 0; int run = 0; ss << "\tMD:Z:";
The modified bowtie -S SAM output in sam.cpp to output the X0 tag for total number of loci that read can map to:
Code:
// Always output stratum ss << "\tXA:i:" << (int)h.stratum; ss << "\tX0:i:" << (int)h.oms;
Thank you very much for any suggestions or comments you may have on this fix/hack.
Best, Aaron