Calculate some statistics for a DNA alignment with Perl
By Abhimanyu Singh 2633 days ago
use Bio::AlignIO;
use Bio::Align::DNAStatistics;
my $stats = Bio::Align::DNAStatistics->new();
my $alignin = Bio::AlignIO->new(-format => 'emboss',
-file => 't/data/insulin.water');
my $aln = $alignin->next_aln;
my $jcmatrix = $stats->distance(-align => $aln,
-method => 'Jukes-Cantor');
print $jcmatrix->print_matrix;
## and for measurements of synonymous /nonsynonymous substitutions ##
my $in = Bio::AlignIO->new(-format => 'fasta',
-file => 't/data/nei_gojobori_test.aln');
my $alnobj = $in->next_aln;
my ($seq1id,$seq2id) = map { $_->display_id } $alnobj->each_seq;
my $results = $stats->calc_KaKs_pair($alnobj, $seq1id, $seq2id);
print "comparing ".$results->[0]{'Seq1'}." and ".$results->[0]{'Seq2'}."\n";
for (sort keys %{$results->[0]} ){
next if /Seq/;
printf("%-9s %.4f \n",$_ , $results->[0]{$_});
}
my $results2 = $stats->calc_all_KaKs_pairs($alnobj);
for my $an (@$results2){
print "comparing ". $an->{'Seq1'}." and ". $an->{'Seq2'}. " \n";
for (sort keys %$an ){
next if /Seq/;
printf("%-9s %.4f \n",$_ , $an->{$_});
}
print "\n\n";
}
my $result3 = $stats->calc_average_KaKs($alnobj, 1000);
for (sort keys %$result3 ){
next if /Seq/;
printf("%-9s %.4f \n",$_ , $result3->{$_});
}
