This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Extract base information from aligned reads in a sam file

I generated a sam (bam) file by aligning fastq to reference using BWA. Now I would like to retrieve the the basepair at a given chromosome position, for example, chr12:11010100:11010100, for all reads aligned to this region, how can I do that? Thanks! (My expected workflow would be: (1) extract reads from sam that cover the querying position; (2) get the basepair from these reads at this position)

sam basepair reads
• 2,568 views • link
written by yliueagle yliueagle

Look into using igvtools which is part of IGV. I assume you are looking to get counts of basepairs (alleles) at a particular location.

• 0 views • link
written by GenoMax GenoMax

Not only counts, but the basepair at that position for each read, with read name included in the output. (Of course the counts can be generated from this table as a downstream step)

• 0 views • link
written by yliueagle yliueagle

Is this good enough? samtools view <your bam> chr:start-end > reads.sam

• 0 views • link
written by JC JC

That will get the reads but not specific base at that position. Not sure what is the best way to do that.

• 0 views • link
written by GenoMax GenoMax

Thank you. This do obtain the reads but I also need to retrieve the basepair information from reads.sam at chr:start-end

• 0 views • link
written by yliueagle yliueagle

1 answer

I found http://lindenb.github.io/jvarkit/Sam2Tsv.html can do the job. Thanks all!

• 0 views • link
written by yliueagle yliueagle

Log in to answer this question.

More posts like this