Hello everyone,
I need to get the precise 5' position of my mapped directional RNA-seq reads. Problem is, for reads that map to the reverse strand, the BAM/SAM format output by my mapping software (TMAP) seems to be the left-most position of the last *aligning* base relative to the reference. In other words, it's the 3'-most base (using the orientation of the read) after clipping non-aligning bases.
I've tried just adding the length of the read onto the position, but this fails for reads with bad 3' ends or several indels in the alignment. I suppose I could write a script to count up the number of "S"s etc. in the cigar, but I am a wet-lab person and the thought of that just makes me sad.
Is there a better way to recover the original 5' position of the seed alignment for reverse reads? Thanks in advance.
I need to get the precise 5' position of my mapped directional RNA-seq reads. Problem is, for reads that map to the reverse strand, the BAM/SAM format output by my mapping software (TMAP) seems to be the left-most position of the last *aligning* base relative to the reference. In other words, it's the 3'-most base (using the orientation of the read) after clipping non-aligning bases.
I've tried just adding the length of the read onto the position, but this fails for reads with bad 3' ends or several indels in the alignment. I suppose I could write a script to count up the number of "S"s etc. in the cigar, but I am a wet-lab person and the thought of that just makes me sad.
Is there a better way to recover the original 5' position of the seed alignment for reverse reads? Thanks in advance.
Comment