I'm trying to call SNPs from a pooled DNA sample using the PoPoolation2 pipeline. I'm stuck on a step that requires me to make a mpileup file from my BAM files. According to multiple tutorials this should work;
samtools mpileup -B input1.bam input2.bam > out.mpileup
However, when I try it makes an empty file.
I've tried adding a reference, re-sorting my BAM files and I checked for errors using Picardtool's ValidateSamFile. Even weirder, I can generate a file using GATK's Pileup command. Any ideas why samtools mpileup
doesn't work for me?
version of samtools ? any error message ? number of reads per bam ? average mapping quality of the reads ? etc...
Samtools v1.9. No error messages. Average 470 million reads per bam. MAPQ = 48.56.
i do have a follow up question. Is using the -A option a good idea? I have same problem as OP (i.e. empty mpileup), but can't find any useful info. of course with -A it seems to work but...how reliable will all downstream analysis be? I asked the same question here but wondering whether you have other ideas...