]> git.donarmstrong.com Git - biopieces.git/blobdiff - code_ruby/lib/maasha/align.rb
finished refactoring s/has_key?/[]/
[biopieces.git] / code_ruby / lib / maasha / align.rb
index 737f2f8f978c4e28477677939f6961d14df344de..c139fe8d89c30a773df458f23e730d2a73e992c4 100755 (executable)
 require 'pp'
 require 'open3'
 require 'narray'
+require 'maasha/align/pair'
 require 'maasha/fasta'
 
 class AlignError < StandardError; end;
 
-ALPH_DNA    = %w{A T C G}
-ALPH_AMBI   = %w{A T C G M R W S Y K V H D B N}
+ALPH_DNA  = %w{A T C G}
+ALPH_AMBI = %w{A T C G M R W S Y K V H D B N}
 
 BIT_INDEL = 0
-BIT_A = 1 << 0
-BIT_T = 1 << 1
-BIT_C = 1 << 2
-BIT_G = 1 << 3
+BIT_A = 1 << 0
+BIT_T = 1 << 1
+BIT_C = 1 << 2
+BIT_G = 1 << 3
 
 BIT_M = BIT_A | BIT_C
 BIT_R = BIT_A | BIT_G
@@ -91,6 +92,8 @@ class Align
   attr_accessor :options
   attr_reader   :entries
 
+  include PairAlign
+
   # Class method to align sequences in a list of Seq objects and
   # return these as a new list of Seq objects.
   def self.muscle(entries, verbose = false)
@@ -99,9 +102,9 @@ class Align
 
     Open3.popen3("muscle", "-quiet") do |stdin, stdout, stderr|
       entries.each do |entry|
-        raise AlignError, "Duplicate sequence name: #{entry.seq_name}" if index.has_key? entry.seq_name
+        raise AlignError, "Duplicate sequence name: #{entry.seq_name}" if index[entry.seq_name]
 
-        index[entry.seq_name] = entry
+        index[entry.seq_name] = entry.dup
 
         stdin.puts entry.to_fasta
       end
@@ -126,6 +129,17 @@ class Align
 
     self.new(result)
   end
+  
+  # Class method to create a pairwise alignment of two given Seq objects. The
+  # alignment is created by casting a search space the size of the sequences
+  # and save the best scoring match between the sequences and recurse into
+  # the left and right search spaced around this match. When all search spaces
+  # are exhausted the saved matches are used to insert gaps in the sequences.
+  def self.pair(q_entry, s_entry)
+    AlignPair.align(q_entry, s_entry)
+
+    self.new([q_entry, s_entry])
+  end
 
   # Method to initialize an Align object with a list of aligned Seq objects.
   def initialize(entries, options = {})
@@ -148,13 +162,13 @@ class Align
     if self.members != 2
       raise AlignError "Bad number of members for similarity calculation: #{self.members}"
     end
-      
+
     na1 = NArray.to_na(@entries[0].seq.upcase, "byte")
     na2 = NArray.to_na(@entries[1].seq.upcase, "byte")
 
     shared   = (na1 - na2).count_false
     total    = (@entries[0].length < @entries[1].length) ? @entries[0].length : @entries[1].length
-    identity = (shared.to_f / total).round(2)
+    identity = shared.to_f / total
 
     identity
   end
@@ -186,6 +200,15 @@ class Align
     output
   end
 
+  # Method for iterating each of the aligned sequences.
+  def each
+    if block_given?
+      @entries.each { |entry| yield entry }
+    else
+      return @entries
+    end
+  end
+
   private
 
   class Consensus
@@ -239,7 +262,7 @@ class Align
     def na_add_entries
       @entries.each_with_index do |entry, i|
         @na_seq[true, i]  = NArray.to_na(entry.seq.downcase.tr(TR_NUC, TR_HEX), "byte")
-        @na_qual[true, i] = NArray.to_na(entry.qual, "byte") - SCORE_BASE if @has_qual
+        @na_qual[true, i] = NArray.to_na(entry.qual, "byte") - Seq::SCORE_BASE if @has_qual
       end
     end
 
@@ -247,14 +270,14 @@ class Align
     # NArrays.
     def consensus_calc
       if @has_qual
-        if @options.has_key? :quality_min
+        if @options[:quality_min]
           mask = mask_quality_min
 
           @na_seq  *= mask
           @na_qual *= mask
         end
 
-        if @options.has_key? :quality_mean
+        if @options[:quality_mean]
           mask = mask_quality_mean
 
           @na_seq  *= mask
@@ -262,21 +285,21 @@ class Align
         end
       end
 
-      if @options.has_key? :sequence_min
+      if @options[:sequence_min]
         mask = mask_sequence_min
 
         @na_seq  *= mask
         @na_qual *= mask if @has_qual
       end
 
-      if @options.has_key? :gap_max
+      if @options[:gap_max]
         mask = mask_gap_max
 
         @na_seq  *= mask
         @na_qual *= mask if @has_qual
       end
 
-      if @options.has_key? :residue_min
+      if @options[:residue_min]
         mask = mask_residue_min
 
         @na_seq  *= mask
@@ -352,7 +375,7 @@ class Align
 
     # Method to calculate a consensus quality from a Consensus object.
     def consensus_qual
-      (@na_qual.mean(1).round.to_type("byte") + SCORE_BASE).to_s
+      (@na_qual.mean(1).round.to_type("byte") + Seq::SCORE_BASE).to_s
     end
   end
 end