This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Get single cigar character for every position from reference start position to end position for a read ?

HI, ref_start_position is 100 and ref_end_position is 175, so is there a function to return every cigar character for each ref position from 100 to 175 . eg . the cigar character of position 110 is 'D', the cigar character of position 172 is 'I' .not just cigar character for ref position matching to referece . the return result can be a tuple ,like [(single cigar character,ref position start ) ...... (single cigar character,ref position end ) ]

Best Wishes!

snp sequencing alignment
• 2,541 views • link
updated by d-cameron d-cameron • written by winter_li • winter_li

2 answers

Given your phrasing I assume you're using pysam, so just write a function that takes the cigartuples tuple and start position expands it accordingly.

• 0 views • link
written by Devon Ryan Devon Ryan

Unfortunately, such a mapping losing information as there can be read sequence between aligned bases. Consider the following alignment against a poly A reference sequence:

1S1M2I1D2I1M
-A--A--A <- reference sequence
SMIIDIIM

and the one-to-one mapping function you are requesting would return [(M,100), (D, 101), (M, 102)]. This loses information about the soft clipped base as well as the two insertions since these read bases have no corresponding reference base thus no corresponding position.

If, for some reason, you want to consider unmapped bases to align to a particular base, you will need a convention as to whether they are left assigned or right assigned.

• 0 views • link
written by d-cameron d-cameron

Log in to answer this question.

More posts like this