X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=code_ruby%2Flib%2Fmaasha%2Fseq.rb;h=c22ecec8948a0230079aacd05af1e552a168819e;hb=eb15b4040e0ceba1c3df40b7f0a597c469771cee;hp=2788e5dcc02d78cdb9e3fd38ef83749bad158b79;hpb=740804b81930af277b29c64f8b8276daff633f44;p=biopieces.git diff --git a/code_ruby/lib/maasha/seq.rb b/code_ruby/lib/maasha/seq.rb index 2788e5d..c22ecec 100644 --- a/code_ruby/lib/maasha/seq.rb +++ b/code_ruby/lib/maasha/seq.rb @@ -27,9 +27,12 @@ require 'maasha/seq/digest' require 'maasha/seq/trim' require 'narray' -autoload :BackTrack, 'maasha/seq/backtrack.rb' -autoload :Dynamic, 'maasha/seq/dynamic.rb' -autoload :Homopolymer, 'maasha/seq/homopolymer.rb' +autoload :BackTrack, 'maasha/seq/backtrack' +autoload :Dynamic, 'maasha/seq/dynamic' +autoload :Homopolymer, 'maasha/seq/homopolymer' +autoload :Hamming, 'maasha/seq/hamming' +autoload :Levenshtein, 'maasha/seq/levenshtein' +autoload :Ambiguity, 'maasha/seq/ambiguity' # Residue alphabets DNA = %w[a t c g] @@ -68,7 +71,6 @@ TRANS_TAB11 = { "GTG" => "V", "GCG" => "A", "GAG" => "E", "GGG" => "G" } - # Error class for all exceptions to do with Seq. class SeqError < StandardError; end @@ -91,7 +93,7 @@ class Seq type = record[:SEQ_TYPE].to_sym if record[:SEQ_TYPE] qual = record[:SCORES] - self.new(seq_name, seq, type, qual) + self.new(seq_name: seq_name, seq: seq, type: type, qual: qual) end # Class method that generates all possible oligos of a specifed length and type. @@ -123,16 +125,20 @@ class Seq oligos end - # Initialize a sequence object with the following arguments: - # - seq_name: Name of the sequence. - # - seq: The sequence. - # - type: The sequence type - DNA, RNA, or protein - # - qual: An Illumina type quality scores string. - def initialize(seq_name = nil, seq = nil, type = nil, qual = nil) - @seq_name = seq_name - @seq = seq - @type = type - @qual = qual + # Initialize a sequence object with the following options: + # - :seq_name Name of the sequence. + # - :seq The sequence. + # - :type The sequence type - DNA, RNA, or protein + # - :qual An Illumina type quality scores string. + def initialize(options = {}) + @seq_name = options[:seq_name] + @seq = options[:seq] + @type = options[:type] + @qual = options[:qual] + + if @seq and @qual and @seq.length != @qual.length + raise SeqError, "Sequence length and score length mismatch: #{@seq.length} != #{@qual.length}" + end end # Method that guesses and returns the sequence type @@ -288,8 +294,8 @@ class Seq raise SeqError, "Missing seq_name" if self.seq_name.nil? or self.seq_name == '' raise SeqError, "Missing seq" if self.seq.nil? or self.seq.empty? - seq_name = self.seq_name.to_s - seq = self.seq.to_s + seq_name = self.seq_name + seq = self.seq.dup unless wrap.nil? seq.gsub!(/(.{#{wrap}})/) do |match| @@ -337,7 +343,14 @@ class Seq # Method to reverse the sequence. def reverse - Seq.new(self.seq_name, self.seq.reverse, self.type, self.qual ? self.qual.reverse : self.qual) + entry = Seq.new( + seq_name: self.seq_name, + seq: self.seq.reverse, + type: self.type, + qual: (self.qual ? self.qual.reverse : self.qual) + ) + + entry end # Method to reverse the sequence. @@ -351,10 +364,11 @@ class Seq def complement raise SeqError, "Cannot complement 0 length sequence" if self.length == 0 - entry = Seq.new - entry.seq_name = self.seq_name - entry.type = self.type - entry.qual = self.qual + entry = Seq.new( + seq_name: self.seq_name, + type: self.type, + qual: self.qual + ) if self.is_dna? entry.seq = self.seq.tr('AGCUTRYWSMKHDVBNagcutrywsmkhdvbn', 'TCGAAYRWSKMDHBVNtcgaayrwskmdhbvn') @@ -384,8 +398,18 @@ class Seq # Method to determine the Hamming Distance between # two Sequence objects (case insensitive). - def hamming_distance(seq) - self.seq.upcase.hamming_distance(seq.seq.upcase) + def hamming_distance(entry, options = nil) + if options and options[:ambiguity] + Hamming.distance(self.seq, entry.seq) + else + self.seq.upcase.hamming_distance(entry.seq.upcase) + end + end + + # Method to determine the Edit Distance between + # two Sequence objects (case insensitive). + def edit_distance(entry) + Levenshtein.distance(self.seq, entry.seq) end # Method that generates a random sequence of a given length and type. @@ -408,7 +432,12 @@ class Seq # Method to return a new Seq object with shuffled sequence. def shuffle - Seq.new(self.seq_name, self.seq.split('').shuffle!.join, self.type, self.qual) + Seq.new( + seq_name: self.seq_name, + seq: self.seq.split('').shuffle!.join, + type: self.type, + qual: self.qual + ) end # Method to shuffle a sequence randomly inline. @@ -475,7 +504,7 @@ class Seq seq_name = self.seq_name.nil? ? nil : self.seq_name.dup - Seq.new(seq_name, seq, self.type, qual) + Seq.new(seq_name: seq_name, seq: seq, type: self.type, qual: qual) end # Method that replaces a sequence with a subsequence from a given start position @@ -640,22 +669,26 @@ class Seq regex_start = Regexp.new(start_codons.join('|'), true) regex_stop = Regexp.new(stop_codons.join('|'), true) - while pos_beg and pos_beg < self.length - size_min - if pos_beg = self.seq.index(regex_start, pos_beg) - if pos_end = self.seq.index(regex_stop, pos_beg) - length = (pos_end - pos_beg) + 3 + while pos_beg = self.seq.index(regex_start, pos_beg) + pos_end = pos_beg + 3 + + while pos_end = self.seq.index(regex_stop, pos_end) + length = (pos_end - pos_beg) + 3 - if (length % 3) == 0 - if size_min <= length and length <= size_max - subseq = self.subseq(pos_beg, length) + if (length % 3) == 0 + if size_min <= length and length <= size_max + subseq = self.subseq(pos_beg, length) - orfs << [subseq, pos_beg, pos_end + 3] - end + orfs << [subseq, pos_beg, pos_end + 3] end + + break end - pos_beg += 1 + pos_end += 1 end + + pos_beg += 1 end if pick_longest