Samtools mpileup. Which reads are excluded during SNP calling?

[wwmm933]

Hi All,

I have two questions about Samtools mpileup. Hope somebody can help me!

First, I found 'mpileup' filters some reads/nts out before SNP calling (compared to command 'pileup'). How can I know which reads are actually filtered?

Second, in 'pileup' output, we can directly see the nucleotides and corresponding base quality scores mapped to one position. Does mpileup provide similar output? Or mpileup just generates the VCF file?

Many thanks!

[dagarfield]

As I understand it, we shouldn't be using pileup anymore. Yes, the output is nice, but all of the same information is (somewhere) in the VCF file.

As for what's being filtered out, however, I don't know. One of the default options is

-Q INT skip bases with baseQ/BAQ smaller than INT [13]

This suggests to me that you might be losing some bases with low baseQs or those with low BAQ scores (explained here: http://samtools.sourceforge.net/mpileup.shtml)

I wish I could be of more help.

[wwmm933]

As I understand it, we shouldn't be using pileup anymore. Yes, the output is nice, but all of the same information is (somewhere) in the VCF file.

As for what's being filtered out, however, I don't know. One of the default options is

-Q INT skip bases with baseQ/BAQ smaller than INT [13]

This suggests to me that you might be losing some bases with low baseQs or those with low BAQ scores (explained here: http://samtools.sourceforge.net/mpileup.shtml)

I wish I could be of more help.

Thank you! So we cannot get the specific nucleotides and corresponding base qualities from samtools anymore, right?

[dagarfield]

It is not obvious to me where that information is in the VCF file, if it is in there at all. However, you might be able to get something in the file generated by mpileup.

Rather than generating a BCF formatted file with mpileup (as is outlined on the man page for mpileup), have you tried running mpileup without the -u and -g options? The output looks a whole lot like the output from old pileup.

--DG

[oiiio]

I've tried to run simply

samtools mpileup -f ref.fasta -b bam > out

but I get "Segmentation Fault" almost immediately.
@dagarfield
Any ideas on what is going on? I though I'd drop you a question here because you said you have run it without the -u and -g options.

Thanks,
Gareth

[dagarfield]

I just ran mine with the following syntax

samtools mpileup -f myGenome.fasta myBam.bam > myoutput.txt

Where myGenome.fasta is a fasta file on which I have run the command

samtools faidx myGenome.fasta

In the same directory in which myGenome.fasta lives.

This looks pretty much like what you did except for the -b option. I think you can (and maybe should) leave that out when you are running just a single BAM file. For the -b option, you'd specify not a BAM file but rather a file that contains a list of the BAM files you want to analyze.

How'd that work?

[oiiio]

You are right about the -b parameter being a list rather than a file. I think its time for some coffee.

Thanks!