Make mpileup file with several BAM files

Question

I'm new to this business, so my question may seem as dumb. I want to make an mpileup file using a reference sequence and several (two to begin with) aligned sorted BAM files. And I want the output file to have certain easy-to-read tab-delimited formatting (which is helpful for subsequent analyses of the output file):

field 1 - chromosome,
field 2 - position,
field 3 - reference base,
field 4 - base(s) in SAMPLE1.BAM for position in field 2,
field 5 - number of reads supporting bases in field 4,
field 6 - base(s) in SAMPLE2.BAM for position in field 2,
field 7 - number of reads supporting bases in field 6.

Something like shown below:

field1   field2   field3   field4   field5   field6   field7
Chr1     15967    T        T,A      15,17    T,_      25,_

where _ means no alternative (or second) base (or zero reads) for this position in the BAM file. '0' or empty space (or nothing) instead of _ is ok. It's ok if information in fields 4-7 is shown not in a single field but in several fields.

Which command is suitable for such a purpose? Thank you.

Answer 1

Have you tried using samtools mpileup in.bam [in2.bam [...]]? Here's information about the output format from the manual page:

The first three columns give the position and reference:

• Chromosome name.
• 1-based position on the chromosome.
• Reference base at this position (this will be “N” on all lines if -f/--fasta-ref has not been used).

The remaining columns show the pileup data, and are repeated for each input BAM file specified:

• Number of reads covering this position.
• Read bases. This encodes information on matches, mismatches, indels, strand, mapping quality, and starts and ends of reads.

This format seems nearly identical to what you've asked for, so hopefully is close enough for you to tweak your downstream analyses.