[Bioperl-l] coverage percentage

Rondon Neto rondonbio at yahoo.com.br
Mon Jun 13 17:44:40 UTC 2011


I need the subject's coverage percentage from an alignment (in this case, a BLAT alignment). I saw that exist an coverage map in Bio::Search::Tiling, but it doesn't help..
So, I wrote the script below, but as you can see, It don't consider that can be non-aligned regions between HSPs. CAn anyone help me?
thank you

Rondon

#!/usr/bin/perl
use warnings;
use strict;
use Bio::SearchIO;

my $in = new Bio::SearchIO ( -format => 'psl',
-file   => $ARGV[0] );

my %hash; my $i=0;
while( my $result = $in->next_result ) {
while( my $hit = $result->next_hit ) {
while( my $hsp = $hit->next_hsp ) {
my $subject = $hit->name(); my $start = $hsp->start('hit'); my $end = $hsp->end('hit'); my $length = $hit->length();

if (defined $hash{$subject}){
if ($start < $hash{$subject}{start}){
$hash{$subject}{start} = $start;
} elsif ($end > $hash{$subject}{end}){
$hash{$subject}{end} = $end;
}
} else {
$hash{$subject}{start} = $start;
$hash{$subject}{end} = $end;
}
$hash{$subject}{length} = $length;
}  
}
}

#calc percentage for each subject
foreach my $key (keys %hash){
my $tam = $hash{$key}{end} - $hash{$key}{start};
my $per = $tam * 100 / $hash{$key}{length};
printf "$key\t%4.2f%%\n", $per;
}

exit;



More information about the Bioperl-l mailing list