bcftools consensus: include the OVERHANG

Question

[Marvin]

I have e.g. a read that overhangs my reference.

It's CIGAR string is e.g. 10S40M.

My goal is to call a consensus that is longer than the original reference, meaning that the 10 soft-clipped bases are taken into account when building the consensus. There might be more than 1 read that overhangs the reference, so I cannot just add the soft-clipped bases to the reference manually.

My approach follows 2 ideas:

1) Modify the CIGAR string, so that 10S40M becomes 50M. I have also tried 10I40M.

2) Add x N's (or gaps "-") at the start of the reference. x is the length of the longest found soft-clip-overhang of all reads.

But using the mpileup-bcftools-call-tabix-bcftools-consensus pipeline results in a reference where I see N's (or gaps) at the beginning instead of the 10 bases that have been soft-clipped.

How do I achieve this?

Answer 1

[Devon Ryan]

You're using the wrong tools for this purpose. What you're trying to do is called "microassembly" and you won't be able to easily trick samtools into ignoring the fact that you don't have a reference sequence for this portion. To just assemble this you might use Mapsembler (this is slow but apparently works OK) or Spades. One of me colleagues has had success with both of these doing something similar.

Answer 2

[JstRoRR]

Is there any way, by using samtools/bcftools or programmatically, to retrieve softclipped bases along with consensus instead of using an assembler??

Any ideas?