[Bioperl-l] Clustalw problems!!!
duxroq
duxroq at hotmail.com
Fri Apr 29 22:45:01 UTC 2011
Hi, I'm not sure whether my program is not finding the clustal w exceutable
or if it is having trouble with the module itself. Here is my error, in the
image below:
http://old.nabble.com/file/p31509401/Untitled.png Untitled.png
here is my code:
#!/usr/bin/local/perl -w
use Bio::Perl;
use Bio::AlignIO;
use Bio::Root::IO;
use Bio::Seq;
use Bio::SeqIO;
use Bio::SimpleAlign;
use Bio::TreeIO;
#--------------------------------------------------------------------------------#
# Main
unless(($#ARGV + 1) > 1) {
print_start_message();
exit;
}
BEGIN { $ENV{CLUSTALDIR} = 'C:\Program Files\clustalw.exe'} #-Set
CLUSTALDIR to correct directory
use Bio::Tools::Run::Alignment::Clustalw;
my $file_name1 = $ARGV[0]; #-name of file containing all sequences
my @sequences = read_all_sequences($file_name1,'fasta');
# alt_revcom(\@sequences); #-reverse complement of every other seq
# print_sequences(@sequences);
my @params = ('outfile' => 'mult_aln.aln'); #-sets parameters for
alignment factory
###################################################
#
#
#
#The error is probably in the next few lines! agh!#
#
#
#
###################################################
my $factory = Bio::Tools::Run::Alignment::Clustalw->new(@params);
$clustalfound = Bio::Tools::Run::Alignment::Clustalw->exists_clustal();
if ($clustalfound) {
print "\n we found it!!!\n" }
my $aln = $factory->align(\@sequences);
print "\nDevins name is bob. \n";
#-creates alignment
print "\nPercentage Identity:\n",$aln->percentage_identity,"\n\n";
my $cons_str = $aln->consensus_iupac(); #-configures consensus using
IUPAC codes
my $cons_name = "Consensus_".save_id(@sequences);
my $cons_seq = new_sequence($cons_str, $cons_name);
write_seq_to_file(">cons.fa",$cons_seq); #-writes consensus to file
my $file_name2 = $ARGV[1];
my $lead_seq = read_sequence($file_name2,'fasta'); #-compare consensus
sequence to leader
# print_sequences($lead_seq);
# print_sequences($cons_seq);
my @lead_cons_seqs = ($lead_seq, $cons_seq); #-forms array for alignment
of leader and consensus
@params = ('pairgap' => 50);
print "\n bobisnotmyname\n";
$factory = Bio::Tools::Run::Alignment::Clustalw->new(@params);
my $aln_lead_cons = $factory->align(\@lead_cons_seqs);
print "\nPercentage
Identity:\n",$aln_lead_cons->percentage_identity,"\n\n";
my $seq_len = $aln_lead_cons->length();
my @aln_seqs = (); #-array of aligned sequences, including gaps
my $i = 0;
foreach $seq ($aln_lead_cons->each_seq() ) {
$aln_seqs[$i] = $seq;
$i++;
}
# print_sequences(@aln_seqs);
my $l_aln_str = ''; #-str of leader sequence from alignment
my $c_aln_str = ''; #-str of consensus sequence from alignment
$found = 0;
$i = 1;
while ($found == 0 && $i < $seq_len) {
$l_aln_str = substr($aln_seqs[0]->seq(),$i,1); #-gets a substring from
l_aln_str
if ($l_aln_str !~ m/\./i) {
#-checks if substring has gap characters
$cons_slice_str = substr($aln_seqs[1]->seq(),$i,490);
$found = 1; #-retrieves 490 characters of consensus where
}
$i++; #-leader begins in alignment
}
$cons_slice_seq = new_sequence($cons_slice_str,"Sliced_".$cons_name);
# print_sequences($cons_slice_seq);
write_seq_to_file(">cons_slice.fa",$cons_slice_seq);
# End Main
#--------------------------------------------------------------------------------#
# Subroutines
sub alt_revcom {
my @sequences = @_;
for ( $i=0; $i <= $#sequences; $i++) {
if ($i%2==1) {
$sequences[$i] = reverse_complement($sequences[$i]);
}
}
}
sub save_id {
my @sequences = @_;
$id = $sequences[0]->display_id;
return $id;
}
sub print_sequences {
my @sequences = @_;
for ($i = 0; $i <= $#sequences; $i++) {
print "Sequence name:",$sequences[$i]->display_id,"\n";
print "Sequence acc:",$sequences[$i]->accession_number,"\n";
print $sequences[$i]->seq(),"\n";
}
}
sub write_seq_to_file {
my ($file_name,$seq) = @_;
write_sequence($file_name,'fasta',$seq);
print "\n",$seq->display_id," written to file.\n\n";
}
--
View this message in context: http://old.nabble.com/Clustalw-problems%21%21%21-tp31509401p31509401.html
Sent from the Perl - Bioperl-L mailing list archive at Nabble.com.
More information about the Bioperl-l
mailing list