view release on metacpan or  search on metacpan

lib/Bio/Tools/Run/Alignment/Clustalw.pm  view on Meta::CPAN

  # For our last example, we'll test a way to perform a local alignment by
  # "anchoring" the alignment to a regular expression.  This is similar
  # to the approach taken in the recent dbclustal program.
  # In principle, we could write a script to search for a good regular expression
  # to use. Instead, here we'll simply choose one manually after looking at the
  # previous alignments.

  sub anchored_align {

      my @local_array = ();
      my @seqs_not_matched = ();

      print "\n Beginning anchored-alignment example... \n";

      for ($seq_num = 0; $seq_num < scalar(@seq_array); $seq_num++) {
          my $seqobj = $seq_array[$seq_num];
          my $seq =  $seqobj->seq();
          my $id =  $seqobj->id();
          # if $regex is not found in the sequence, save sequence id name and set
          # array value =0 for later
          unless ($seq =~/$regex/) {
              $local_array[$seq_num] = 0;
              push (@seqs_not_matched, $id) ;
              next;
          }
          # find positions of start and of subsequence to be aligned
          my $match_start_pos = length($`);
          my $match_stop_pos = length($`) + length($&);
          my $start =  ($match_start_pos - $extension) > 1 ?
              ($match_start_pos - $extension) +1 : 1;
          my $stop =  ($match_stop_pos + $extension) < length($seq) ?
              ($match_stop_pos + $extension) : length($seq);
          my $string = $seqobj->subseq($start, $stop);

          $local_array[$seq_num] = Bio::Seq->new(-id=>$id, -seq=>$string);
      }
      @local_array = grep $_ , @local_array; # remove array entries with no match

      # Perform alignment on the local segments of the sequences which match "anchor"
      $aln = $factory->align(\@local_array);
      my $consensus  = $aln->consensus_string(); # Get consensus of local alignment

      if (scalar(@seqs_not_matched) ) {
          print " Sequences not matching $regex : @seqs_not_matched \n"
      } else {
          print " All sequences match $regex : @seqs_not_matched \n"
  }
      print "Consensus sequence of local alignment: $consensus \n";

      return 1;
  }

  #----------------
  sub clustalw_usage {
  #----------------