Method for detecting soft clipped read pileup?

[ndovu9]

I'm looking for a way to detect positions in the genome where there is a pileup of soft clipped reads. I've attached an image with an example of what this situation would look like in IGV. Essentially, I'm looking for a tool similar to samtools mpileup, but with the important difference that I want to see soft-clipped reads. I think the issue is that soft-clipped based are not technically aligned at these positions. I know that I could write a script to parse the CIGAR string for each read and detect locations like these, but I'm wondering, is there a tool that can quickly report the locations where reads start getting soft-clipped?


I'm imagining a version of samtools mpileup that would report something like this:
10 141352 N 105 a$A$A$a$aSASASaSaaaaAAaAAaAaAaaaAaAaAAAaAAaaAAaAAaaAaAaAAaAAAaAAaAaAAAAAAaAAaaAaaaaaaaaaAaaaAaAAaAaaAaaAAaaaaAaA^]a @<@?;?>[email protected]@[email protected]????A>@@[email protected]@[email protected]@@?>@:>[email protected]@[email protected][email protected]@>[email protected]@@[email protected][email protected][email protected][email protected]@?>[email protected][email protected]>>@>[email protected]@@?A>@>?A>[email protected]=?>=??>=?=C>9>

where the "$", as usual, means that reads are ending at this position, whereas the "S" would mean that bases are "aligned" and soft-clipped at this position.

[martinghunt]

This is implemented as part of the REAPR pipeline.
http://www.sanger.ac.uk/resources/software/reapr/

You can run the pipeline with:
reapr pipeline reference.fasta mapped_reads.bam output_directory

Note that it assumes that the orientation of the reads in the input bam is to point towards each other.

One of its early stages is to write a file of stats at each base of the reference:
01.stats.per_base.gz
(If you don't care about other results, then you can kill the pipeline when it starts writing files with names beginning with 02.)

Columns 1-2 are the chromosome name and position.
Columns 17-20 have the number of reads that had their ends soft-clipped at that position. 4 columns because it's broken down into the 4 combinations from reads mapped to the fwd or rev strand and clipped at their start or end (w.r.t to the reference).

Hope that helps,
Martin