[Bioperl-l] Error calling alignment method
Lorenzo Carretero Paulet
locarpau at upvnet.upv.es
Wed May 25 17:58:15 UTC 2011
Hi all,
I'm trying to create the following subroutine but I get the error message:
"Can't call method "alignment" on an undefined value at
/Users/marioafares/Documents/workspace/PlantEvolGen/test.pl line 80,
<GEN1> line 2."
which indicates $dna_aln is undefined. However, I manage to print it
using Dumper, so I guess it is actually defined.
Can anyone see where the problem is?
I attach the two files I'm using in the test.
Thanks in advance,
Lorenzo
#!/usr/local/bin/perl -w
use 5.010;
use strict;
use Data::Dumper;
# definition of the environmental variable CLUSTALDIR
BEGIN {$ENV{CLUSTALDIR} =
'/Applications/Bioinformatics/clustalw-2.1-macosx/'}
use Bio::Tools::Run::Alignment::Clustalw;
use Bio::Align::Utilities qw(aa_to_dna_aln);
BEGIN {$ENV{PAMLDIR} = '/Applications/Bioinformatics/paml44/bin/'}
use Bio::Tools::Run::Phylo::PAML::Codeml;
my $sequencesfilenameAA =
"/Users/marioafares/Documents/SequenceDatabase/plaza_public_02_Apr27/plaza_public_02/BLAST_Parsed_results/PerSpecies/test_vs_test.par.aa.1.fas";
my $sequencesfilenameNT =
"/Users/marioafares/Documents/SequenceDatabase/plaza_public_02_Apr27/plaza_public_02/BLAST_Parsed_results/PerSpecies/test_vs_test.par.nt.1.fas";
my $format = 'fasta';
GettingBioperlAlignmentAAtoDNAplusPAMLcalculation
($sequencesfilenameAA, $sequencesfilenameNT, $format);
sub GettingBioperlAlignmentAAtoDNAplusPAMLcalculation
{
my ( $sequencesfilenameAA, $sequencesfilenameNT, $format, $ktuple,
$matrix ) = @_;
my %hashNTseqs = ();
my $likelihood;
my $Ks;
my $Ka;
my $omega;
my $factory = Bio::Tools::Run::Alignment::Clustalw->new (); #use
default parameters in here
my $pep_aln = $factory->align($sequencesfilenameAA);
my $inseq = Bio::SeqIO->new(-file => "<$sequencesfilenameNT",
-format => $format );
while (my $seq = $inseq->next_seq)
{
my $seq_id = $seq->display_id();
$hashNTseqs{$seq_id} = $seq;
}
my $dna_aln = aa_to_dna_aln($pep_aln, \%hashNTseqs);
say Dumper \%hashNTseqs;
say Dumper $dna_aln;
###########################################################
#my $codeml_runs = 5;
my $ks_factory = new Bio::Tools::Run::Phylo::PAML::Codeml
( -params => {
#'verbose' => 0,
#'noisy' => 0,
'runmode' => -2,
'seqtype' => 1,
#'model' => 0,
#'NSsites' => 0,
#'icode' => 0,
#'fix_alpha' => 0,
#'fix_kappa' => 0,
#'RateAncestor' => 0,
'CodonFreq' => 2,
'cleandata' => 1,
'ndata' => 1
},
);
my $kaks_factory->alignment($dna_aln);
say "\nCalculating Ks-values of current cluster ...";
# The estimation of the Ks-values is repeated $codeml_runs
times ...
# for(my $k=1;$k <= $codeml_runs;$k++)
# {
# print "\nCodeml-Run $k:\n\n";
# Ka and Ks-vlaues are calculated using codeml
my ($rc,$parser) = $kaks_factory->run();
#$kaks_factory->cleanup();
# If the calculation was succsessful ...
# if($rc != 0)
# {
my $result = $parser->next_result;#not sure what it
does here
#my $NGmatrix = $result->get_MLmatrix();
my $MLmatrix = $result->get_MLmatrix();
$likelihood = $MLmatrix->[0]->[1]->{'lnL'};
$Ks = $MLmatrix->[0]->[1]->{'dS'};
$Ka = $MLmatrix->[0]->[1]->{'dA'};
$omega = $MLmatrix->[0]->[1]->{'omega'};
print " likelihood = $likelihood, Ka = $Ka, Ks =
$Ks, Ka/Ks = $omega\n";
# }
# else # If an error occured during the Ks-calculation ...
# {
# return (-1);
# }
# }
return ( $likelihood, $Ks, $Ka, $omega );
}
--
*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
Lorenzo Carretero Paulet
Institute for Plant Molecular and Cell Biology - IBMCP (CSIC-UPV)
Integrative Systems Biology Group
C/ Ingeniero Fausto Elio s/n.
46022 Valencia, Spain
Phone: +34 963879934
Fax: +34 963877859
e-mail:locarpau at upvnet.upv.es
*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: test_vs_test.par.nt.1.fas
URL: <http://lists.open-bio.org/pipermail/bioperl-l/attachments/20110525/469cb211/attachment-0002.ksh>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: test_vs_test.par.aa.1.fas
URL: <http://lists.open-bio.org/pipermail/bioperl-l/attachments/20110525/469cb211/attachment-0003.ksh>
More information about the Bioperl-l
mailing list