I would like to perform a
samtools mpileup from a single file that contains thousands of read groups with different SM tags. I could split the bam by read group using
samtools split, and then perform the mpileup, but splitting into thousands of files and then performing an mpileup from thousands of files is crazy slow on my machine.
I could imagine an implementation using pysam, but I was hoping someone has a solution already.
Bcftools mpileup already has this feature, but the feature was never implemented for the original mpileup format, which is what I need.
@HD VN:1.6 SO:unknown @SQ SN:chr1 LN:1000000 @RG ID:0 SM:sample_0 @RG ID:1 SM:sample_1 r0 0 chr1 24 0 1M * 0 0 G I RG:Z:0 r1 0 chr1 24 0 1M * 0 0 G I RG:Z:1
When I run samtools mpileup I get this behavior:
> samtools mpileup input.sam [mpileup] 2 samples in 1 input files chr1 24 N 2 ^!G$^!G$ II
Whereas, what I actually want is this output:
[mpileup] 2 samples in 1 input files chr1 24 N 1 ^!G$ I 1 ^!G$ I
For a long time I thought this was a bug in samtools mpileup, but it appears that it is a feature that was described/hinted at in the samtools documentation, but was never implemented: https://github.com/samtools/samtools/issues/599#issuecomment-604941884