Как найти индексы идентичных подпоследовательностей в двух строках в Ruby?

Question

Здесь каждый экземпляр класса DNA соответствует строке, такой как 'GCCCAC'. Массивы подстрок, содержащие k-мерс , могут быть построены из этих строк. Для этой струны есть 1-мерный, 2-мерный, 3-мерный, 4-мерный, 5-мерный и один 6-мерный:

• 6 1-мерный: ["G", "C", "C", "C", "A", "C"]
• 5 2-член: ["GC", "CC", "CC", "CA", "AC"]
• 4 3-член: ["GCC", "CCC", "CCA", "CAC"]
• 3 4-член: ["GCCC", "CCCA", "CCAC"]
• 2 5-член:["GCCCA", "CCCAC"]
• 1 6-мер: ["GCCCAC"]

Шаблон должен быть очевидным. Подробнее см. Wiki .

Проблема состоит в том, чтобы написать метод shared_kmers (k, dna2) класса DNA, который возвращает массив всех пар [i, j], где находится эта ДНКобъект (который получает сообщение) делит с dna2 общий k-мер в позиции i в этой ДНК и в позиции j в dna2.

dna1 = DNA.new('GCCCAC')
dna2 = DNA.new('CCACGC')

dna1.shared_kmers(2, dna2)
#=> [[0, 4], [1, 0], [2, 0], [3, 1], [4, 2]]

dna2.shared_kmers(2, dna1)
#=> [[0, 1], [0, 2], [1, 3], [2, 4], [4, 0]]

dna1.shared_kmers(3, dna2)
#=> [[2, 0], [3, 1]]

dna1.shared_kmers(4, dna2)
#=> [[2, 0]]

dna1.shared_kmers(5, dna2)
#=> []

Answer 1

class DNA
  attr_accessor :sequencing

  def initialize(sequencing)
    @sequencing = sequencing
  end

  def kmers(k)
    @sequencing.each_char.each_cons(k).map(&:join)
  end

  def shared_kmers(k, dna)
    kmers(k).each_with_object([]).with_index do |(kmer, result), index|
      dna.kmers(k).each_with_index do |other_kmer, other_kmer_index|
        result << [index, other_kmer_index] if kmer.eql?(other_kmer)
      end
    end
  end
end

dna1 = DNA.new('GCCCAC')
dna2 = DNA.new('CCACGC')

dna1.kmers(2)
#=> ["GC", "CC", "CC", "CA", "AC"]

dna2.kmers(2)
#=> ["CC", "CA", "AC", "CG", "GC"]

dna1.shared_kmers(2, dna2)
#=> [[0, 4], [1, 0], [2, 0], [3, 1], [4, 2]]

dna2.shared_kmers(2, dna1)
#=> [[0, 1], [0, 2], [1, 3], [2, 4], [4, 0]]

dna1.shared_kmers(3, dna2)
#=> [[2, 0], [3, 1]]

dna1.shared_kmers(4, dna2)
#=> [[2, 0]]

dna1.shared_kmers(5, dna2)
#=> []

Answer 2

Я рассмотрю только суть вашей проблемы, без ссылки на класс DNA. Должно быть легко реорганизовать то, что следует, довольно легко.

Код

def match_kmers(s1, s2, k)
  h1 = dna_to_index(s1, k)
  h2 = dna_to_index(s2, k)
  h1.flat_map { |k,_| h1[k].product(h2[k] || []) }
end

def dna_to_index(dna, k)
  dna.each_char.
      with_index.
      each_cons(k).
      with_object({}) {|arr,h| (h[arr.map(&:first).join] ||= []) << arr.first.last}
end

Примеры

dna1 = 'GCCCAC'
dna2 = 'CCACGC'

match_kmers(dna1, dna2, 2)
  #=> [[0, 4], [1, 0], [2, 0], [3, 1], [4, 2]] 
match_kmers(dna2, dna1, 2)
  #=> [[0, 1], [0, 2], [1, 3], [2, 4], [4, 0]] 

match_kmers(dna1, dna2, 3)
  #=> [[2, 0], [3, 1]] 
match_kmers(dna2, dna1, 3)
  #=> [[0, 2], [1, 3]] 

match_kmers(dna1, dna2, 4)
  #=> [[2, 0]] 
match_kmers(dna2, dna1, 4)
  #=> [[0, 2]] 

match_kmers(dna1, dna2, 5)
  #=> [] 
match_kmers(dna2, dna1, 5)
  #=> [] 

match_kmers(dna1, dna2, 6)
  #=> [] 
match_kmers(dna2, dna1, 6)
  #=> [] 

Объяснение

Рассмотрим dna1 = 'GCCCAC'. Он содержит 5 2-мер (k = 2):

dna1.each_char.each_cons(2).to_a.map(&:join)
  #=> ["GC", "CC", "CC", "CA", "AC"] 

Аналогично для dna2 = 'CCACGC':

dna2.each_char.each_cons(2).to_a.map(&:join)
  #=> ["CC", "CA", "AC", "CG", "GC"]

Это ключи хэшей, которые dna_to_index создает дляdna1 и dna2 соответственно. Значения хэша - это массивы индексов, где соответствующий ключ начинается в строке ДНК. Давайте вычислим эти хэши для k = 2:

h1 = dna_to_index(dna1, 2)
  #=> {"GC"=>[0], "CC"=>[1, 2], "CA"=>[3], "AC"=>[4]} 
h2 = dna_to_index(dna2, 2)
  #=> {"CC"=>[0], "CA"=>[1], "AC"=>[2], "CG"=>[3], "GC"=>[4]} 

h1 показывает, что:

• "GC" начинается с индекса 0 dna1
• "CC" начинается с индексов 1 и 2 из dna1
• "CA" начинается с индекса 3 из dna1
• "CC" начинается с индекса 4 из dna1

h2 имеет аналогичную интерпретацию. См. Enumerable # flat_map и Array # product .

Затем метод match_kmers используется для построения желаемого массива пар индексов [i, j], таких, что h1[i] = h2[j].

Теперь давайте посмотрим на хэши, созданные для 3-мер (k = 3):

h1 = dna_to_index(dna1, 3)
  #=> {"GCC"=>[0], "CCC"=>[1], "CCA"=>[2], "CAC"=>[3]} 
h2 = dna_to_index(dna2, 3)
  #=> {"CCA"=>[0], "CAC"=>[1], "ACG"=>[2], "CGC"=>[3]} 

Мы видим, что первый 3-элемент в dna1 равен "GCC", начиная с индекса 0. Однако этот 3-мер отсутствует в dna2, поэтому в возвращаемом массиве нет элементов [0, X] (X является просто заполнителем). "CCC" также не является ключом ко второму хешу. "CCA" и "CAC" присутствуют во втором хеше, поэтому возвращаемый массив:

h1["CCA"].product(h2["CCA"]) + h1["CAC"].product(h2["CAC"]) 
  #=> [[2, 0]] + [[3, 1]]
  #=> [[2, 0], [3, 1]]

Answer 3

Я бы начал с написания метода для перечисления подпоследовательностей заданной длины (т.е. k-мер):

class DNA
  def initialize(sequence)
    @sequence = sequence
  end

  def each_kmer(length)
    return enum_for(:each_kmer, length) unless block_given?

    0.upto(@sequence.length - length) { |i| yield @sequence[i, length] }
  end
end

DNA.new('GCCCAC').each_kmer(2).to_a
#=> ["GC", "CC", "CC", "CA", "AC"]

Кроме того, вы можете легко собрать индексы идентичных k-мер, используявложенный цикл:

class DNA
  # ...

  def shared_kmers(length, other)
    indices = []
    each_kmer(length).with_index do |k, i|
      other.each_kmer(length).with_index do |l, j|
        indices << [i, j] if k == l
      end
    end
    indices
  end
end

dna1 = DNA.new('GCCCAC')
dna2 = DNA.new('CCACGC')

dna1.shared_kmers(2, dna2)
#=> [[0, 4], [1, 0], [2, 0], [3, 1], [4, 2]]

К сожалению, приведенный выше код пересекает other.each_kmer для каждого k-mer в приемнике. Мы можем оптимизировать это, создав хеш, содержащий все индексы для каждого k-mer в other заранее:

class DNA
  # ...

  def shared_kmers(length, other)
    hash = Hash.new { |h, k| h[k] = [] }
    other.each_kmer(length).with_index { |k, i| hash[k] << i }

    indices = []
    each_kmer(length).with_index do |k, i|
      hash[k].each { |j| indices << [i, j] }
    end
    indices
  end
end