[Bioperl-l] Getting pileup consensus from BAM files using Bio::DB::Sam

[Gowthaman Ramasamy]

Hi Lincoln,
Thats a good lead. I will try to use MAQ in perl rather than using my simple majority rule.

-gowtham
________________________________________
From: Lincoln Stein [lincoln.stein at gmail.com]
Sent: Tuesday, August 03, 2010 9:57 AM
To: Gowthaman Ramasamy
Cc: bioperl-l
Subject: Re: [Bioperl-l] Getting pileup consensus from BAM files using  Bio::DB::Sam

Samtools is running MAQ on the pileup. You could either implement MAQ in perl, or come up with your own consensus caller.

Lincoln

On Tue, Aug 3, 2010 at 1:29 AM, Gowthaman Ramasamy <gowthaman.ramasamy at seattlebiomed.org<mailto:gowthaman.ramasamy at seattlebiomed.org>> wrote:
Hi List,
I am trying to find out the consensus using pileup via Bio::DB::Sam. Using the following script I could parse out the ref_base and different bases from reads at that position. Though, I am not able to find a method to derive consensus. Similar to the values produced by "samtools pileup -c -f xxxxxx.fasta yyyyyyy.bam".

The script I use now retrives ref base, query bases for each position. How do I improve it to get the consensus?

Thanks very much in advance,
Gowthaman


use Bio::DB::Sam;

my $bam = Bio::DB::Sam->new(-bam => 'something.bam',
                           -fasta => 'something.fasta'
                          );

my $cb = sub {
                       my ($seqid, $pos, $pileups) = @_;
                       my $refBase = $bam->segment($seqid, $pos, $pos)->dna;
                       print "\n$pos\t$refBase=>";
                       for my $pileup (@$pileups){
                               my $al = $pileup->alignment;
                               my $qBase = substr($al->qseq, $pileup->qpos, 1);
                               print "$qBase,";
                               }
                       };

$bam->pileup('Lin.chr10i', $cb);

_______________________________________________
Bioperl-l mailing list
Bioperl-l at lists.open-bio.org<mailto:Bioperl-l at lists.open-bio.org>
http://lists.open-bio.org/mailman/listinfo/bioperl-l



--
Lincoln D. Stein
Director, Informatics and Biocomputing Platform
Ontario Institute for Cancer Research
101 College St., Suite 800
Toronto, ON, Canada M5G0A3
416 673-8514
Assistant: Renata Musa <Renata.Musa at oicr.on.ca<mailto:Renata.Musa at oicr.on.ca>>