Class: ViralSeq::SeqHashPair

Inherits:
Object
  • Object
show all
Defined in:
lib/viral_seq/seq_hash_pair.rb

Overview

Class for paired-end sequences.

Examples:

initialize a new SeqHashPair object from a directory containing paired-end sequences

my_seqhashpair = ViralSeq::SeqHashPair.fa('my_seq_directory')

join the paired-end sequences with an overlap of 100 bp

my_seqhashpair.join1(100)

join the paired-end sequences with unknown overlap, each pair of sequences has its own overlap size

my_seqhashpair.join2(model: :indiv)

Instance Attribute Summary collapse

Class Method Summary collapse

Instance Method Summary collapse

Constructor Details

#initialize(dna_hash = {}, title = "", file = []) ⇒ SeqHashPair

initialize SeqHashPair object with @dna_hash, @title and @file



16
17
18
19
20
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 16

def initialize (dna_hash = {}, title = "", file = [])
  @dna_hash = dna_hash
  @title = title
  @file = file
end

Instance Attribute Details

#dna_hash ⇒ Hash

Returns Hash object for :name => [:r1_sequence_string, :r2_sequence_string].

Returns:

  • (Hash) —

    Hash object for :name => [:r1_sequence_string, :r2_sequence_string]



24
25
26
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 24

def dna_hash
  @dna_hash
end

#file ⇒ String

Returns the r1 and r2 files that are used to initialize SeqHash object, if they exist.

Returns:

  • (String) —

    the r1 and r2 files that are used to initialize SeqHash object, if they exist



33
34
35
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 33

def file
  @file
end

#title ⇒ String

default as the directory basename if SeqHash object is initialized using ::fa

Returns:

  • (String) —

    the title of the SeqHash object.



29
30
31
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 29

def title
  @title
end

Class Method Details

.new_from_fasta(indir) ⇒ ViralSeq::SeqHashPair Also known as: fa

initialize a new ViralSeq::SeqHashPair object from a directory containing paired sequence files in the FASTA format

Examples:

initialize a new SeqHashPair object from a directory containing paired-end sequences

my_seqhashpair = ViralSeq::SeqHashPair.fa('spec/sample_paired_seq')

Parameters:

  • indir (String) —

    directory containing paired sequence files in the FASTA format,

    Paired sequence files need to have “r1” and “r2” in their file names

    Example for the file structure

    The sequence taxa should only differ by last 3 characters to distinguish r1 and r2 sequence.

Returns:



49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 49

def self.new_from_fasta(indir)
  files = Dir[indir + "/*"]
  r1_file = ""
  r2_file = ""
  files.each do |f|
    if File.basename(f) =~ /r1/i
      r1_file = f
    elsif File.basename(f) =~ /r2/i
      r2_file = f
    end
  end

  seq1 = ViralSeq::SeqHash.fa(r1_file).dna_hash
  seq2 = ViralSeq::SeqHash.fa(r2_file).dna_hash

  new_seq1 = seq1.each_with_object({}) {|(k, v), h| h[k[0..-4]] = v}
  new_seq2 = seq2.each_with_object({}) {|(k, v), h| h[k[0..-4]] = v}

  seq_pair_hash = {}

  new_seq1.each do |seq_name,seq|
    seq_pair_hash[seq_name] = [seq, new_seq2[seq_name]]
  end
  seq_hash = ViralSeq::SeqHashPair.new
  seq_hash.dna_hash = seq_pair_hash
  seq_hash.title = File.basename(indir,".*")
  seq_hash.file = [r1_file, r2_file]
  return seq_hash
end

Instance Method Details

#join1(overlap = 0, diff = 0.0) ⇒ ViralSeq::SeqHash

Pair-end join function for KNOWN overlap size.

Examples:

join paired-end sequences with different :diff cut-offs, overlap provided.

paired_seqs = {">pair1"=>["GGGGGGGGGGAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
                          "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAATTTTTTTTTT"],
               ">pair2"=>["GGGGGGGGGGAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
                          "AAAAAAAAAAGAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAATTTTTTTTTT"],
               ">pair3"=>["GGGGGGGGGGAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
                          "AAAAAAAAAAGGAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAATTTTTTTTTT"]}
my_seqhashpair = ViralSeq::SeqHashPair.new(paired_seqs)
my_seqhashpair.join1(100).dna_hash.keys
=> [">pair1"]
my_seqhashpair.join1(100,0.01).dna_hash.keys
=> [">pair1", ">pair2"]
my_seqhashpair.join1(100,0.02).dna_hash.keys
=> [">pair1", ">pair2", ">pair3"]

Parameters:

  • overlap (Integer) (defaults to: 0) —

    how many bases are overlapped. ‘0` means no overlap, R1 and R2 will be simply put together.

  • diff (Integer, Float) (defaults to: 0.0) —

    the maximum mismatch rate allowed for the overlapping region. default at 0.0, i.e. no mis-match allowed.

Returns:

Raises:

  • (ArgumentError) 


108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 108

def join1(overlap = 0, diff = 0.0)
  seq_pair_hash = self.dna_hash
  raise ArgumentError.new(":overlap has to be Integer, input #{overlap} invalid.") unless overlap.is_a? Integer
  raise ArgumentError.new(":diff has to be float or integer, input #{diff} invalid.") unless (diff.is_a? Integer or diff.is_a? Float)
  joined_seq = {}
  seq_pair_hash.each do |seq_name,seq_pair|
    r1_seq = seq_pair[0]
    r2_seq = seq_pair[1]
    if overlap.zero?
      joined_sequence = r1_seq + r2_seq
    elsif diff.zero?
      if r1_seq[-overlap..-1] == r2_seq[0,overlap]
        joined_sequence= r1_seq + r2_seq[overlap..-1]
      end
    elsif r1_seq[-overlap..-1].compare_with(r2_seq[0,overlap]) <= (overlap * diff)
      joined_sequence= r1_seq + r2_seq[overlap..-1]
    else
      next
    end
    joined_seq[seq_name] = joined_sequence if joined_sequence
  end

  joined_seq_hash = ViralSeq::SeqHash.new
  joined_seq_hash.dna_hash = joined_seq
  joined_seq_hash.title = self.title + "_joined"
  joined_seq_hash.file = File.dirname(self.file[0]) if self.file.size > 0
  return joined_seq_hash
end

#join2(model: :con, diff: 0.0) ⇒ ViralSeq::SeqHash

Pair-end join function for UNKNOWN overlap.

Examples:

join paired-end sequences, overlap NOT provided

paired_seq2 = {">pair4" => ["AAAGGGGGGG", "GGGGGGGTT"],
               ">pair5" => ["AAAAAAGGGG", "GGGGTTTTT"],
               ">pair6" => ["AAACAAGGGG", "GGGGTTTTT"] }
my_seqhashpair = ViralSeq::SeqHashPair.new(paired_seq2)
my_seqhashpair.join2.dna_hash
=> {">pair4"=>"AAAGGGGGGGGGGTT", ">pair5"=>"AAAAAAGGGGTTTTT", ">pair6"=>"AAACAAGGGGTTTTT"}
my_seqhashpair.join2(model: :indiv).dna_hash
=> {">pair4"=>"AAAGGGGGGGTT", ">pair5"=>"AAAAAAGGGGTTTTT", ">pair6"=>"AAACAAGGGGTTTTT"}

Parameters:

  • model (Symbol) (defaults to: :con) —

    models used to determine the overlap, ‘:con`, `:indiv`

    model ‘:con`: overlap is determined based on consensus, all sequence pairs are supposed to have the same overlap size

    note: minimal overlap as 4 bases.

    model ‘:indiv`: overlap is determined for each sequence pair, sequence pairs can have different size of overlap

  • diff (Integer, Float) (defaults to: 0.0) —

    the maximum mismatch rate allowed for the overlapping region. default at 0.0, i.e. no mis-match allowed.

Returns:



157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 157

def join2(model: :con, diff: 0.0)
  seq_pair_hash = self.dna_hash
  begin
    raise ArgumentError.new(":diff has to be float or integer, input #{diff} invalid.") unless (diff.is_a? Integer or diff.is_a? Float)
    if model == :con
      overlap = determine_overlap_pid_pair(seq_pair_hash, diff)
      return self.join1(overlap, diff)
    elsif model == :indiv
      joined_seq = {}
      seq_pair_hash.each do |seq_name, seq_pair|
        overlap_list = []
        overlap_matrix(seq_pair[0], seq_pair[1]).each do |overlap1, diff_nt|
          cut_off_base = overlap1 * diff
          overlap_list << overlap1 if diff_nt <= cut_off_base
        end
        if overlap_list.empty?
          joined_seq[seq_name] = seq_pair[0] + seq_pair[1]
        else
          overlap = overlap_list.max
          joined_seq[seq_name] = seq_pair[0] + seq_pair[1][overlap..-1]
        end
      end
      joined_seq_hash = ViralSeq::SeqHash.new
      joined_seq_hash.dna_hash = joined_seq
      joined_seq_hash.title = self.title + "_joined"
      joined_seq_hash.file = File.dirname(self.file[0]) if self.file.size > 0
      return joined_seq_hash
    else
      raise ArgumentError.new("Error::Wrong Overlap Model Argument. Given \`#{model}\`, expected `:con` or `:indiv`.")
    end
  rescue ArgumentError => e
    puts e
    return nil
  end
end

#size ⇒ Integer

the size of nt sequence hash of the SeqHashPair object

Returns:

  • (Integer) —

    size of nt sequence hash of the SeqHash object



85
86
87
 
# File 'lib/viral_seq/seq_hash_pair.rb', line 85

def size
  self.dna_hash.size
end