Как найти перекрытие между двумя последовательностями и вернуть его

Я новичок в Python, и уже потратил много часов с этой проблемой, надеюсь, кто-нибудь может мне помочь. Мне нужно найти перекрытие между двумя последовательностями. Перекрытие находится в левом конце первых последовательностей и в правом конце второго. Я хочу, чтобы функция нашла перекрытие и вернула его.

Мои последовательности:

s1 = "CGATTCCAGGCTCCCCACGGGGTACCCATAACTTGACAGTAGATCTC"
s2 = "GGCTCCCCACGGGGTACCCATAACTTGACAGTAGATCTCGTCCAGACCCCTAGC"

Моя функция должна быть названа

def getOverlap(left, right)

С s1 будучи левой последовательностью, а s2 быть правильным.

Результат должен быть

‘GGCTCCCCACGGGGTACCCATAACTTGACAGTAGATCTC’

Любая помощь приветствуется

4 ответа

Посмотрите на difflib библиотека, а точнее в find_longest_match():

import difflib

def get_overlap(s1, s2):
    s = difflib.SequenceMatcher(None, s1, s2)
    pos_a, pos_b, size = s.find_longest_match(0, len(s1), 0, len(s2)) 
    return s1[pos_a:pos_a+size]

s1 = "CGATTCCAGGCTCCCCACGGGGTACCCATAACTTGACAGTAGATCTC"
s2 = "GGCTCCCCACGGGGTACCCATAACTTGACAGTAGATCTCGTCCAGACCCCTAGC"

print(get_overlap(s1, s2)) # GGCTCCCCACGGGGTACCCATAACTTGACAGTAGATCTC

Вы могли бы использовать difflib.SequenceMatcher:

d = difflib.SequenceMatcher(None,s1,s2)
>>> match = max(d.get_matching_blocks(),key=lambda x:x[2])
>>> match
Match(a=8, b=0, size=39)
>>> i,j,k = match
>>> d.a[i:i+k]
'GGCTCCCCACGGGGTACCCATAACTTGACAGTAGATCTC'
>>> d.a[i:i+k] == d.b[j:j+k]
True
>>> d.a == s1
True
>>> d.b == s2
True

Алгоритм Кнута-Морриса-Пратта - хороший метод для нахождения одной строки внутри другой (так как я видел ДНК, я предполагаю, что вы хотели бы, чтобы это масштабировалось до... миллиардов?).

# Knuth-Morris-Pratt string matching
# David Eppstein, UC Irvine, 1 Mar 2002

from __future__ import generators

def KnuthMorrisPratt(text, pattern):

    '''Yields all starting positions of copies of the pattern in the text.
Calling conventions are similar to string.find, but its arguments can be
lists or iterators, not just strings, it returns all matches, not just
the first one, and it does not need the whole text in memory at once.
Whenever it yields, it will have read the text exactly up to and including
the match that caused the yield.'''

    # allow indexing into pattern and protect against change during yield
    pattern = list(pattern)

    # build table of shift amounts
    shifts = [1] * (len(pattern) + 1)
    shift = 1
    for pos in range(len(pattern)):
        while shift <= pos and pattern[pos] != pattern[pos-shift]:
            shift += shifts[pos-shift]
        shifts[pos+1] = shift

    # do the actual search
    startPos = 0
    matchLen = 0
    for c in text:
        while matchLen == len(pattern) or \
              matchLen >= 0 and pattern[matchLen] != c:
            startPos += shifts[matchLen]
            matchLen -= shifts[matchLen]
        matchLen += 1
        if matchLen == len(pattern):
            yield startPos

Ссылка, где я получил код Python KMP (и встроенный, который будет быстрее для небольших проблем из-за постоянной времени выполнения).

Для передовой производительности используйте таблицу префиксов и окна хеш-функции вашей строки в качестве целых чисел с основанием 4 (в биологии вы бы назвали их k-mers или oligos).;)

Удачи!

РЕДАКТИРОВАТЬ: Есть также хороший трюк, где вы сортируете список, содержащий каждый префикс (всего n) в первой строке и каждый префикс (всего n) во второй строке. Если они имеют наибольшую общую подпоследовательность, то они должны быть смежными в отсортированном списке, поэтому найдите элемент из другой строки, которая находится ближе всего в отсортированном списке, а затем возьмите самый длинный префикс, который полностью соответствует.:)

Самая длинная общая подстрока

def LongestCommonSubstring(S1, S2):
  M = [[0]*(1+len(S2)) for i in xrange(1+len(S1))]
  longest, x_longest = 0, 0
  for x in xrange(1,1+len(S1)):
    for y in xrange(1,1+len(S2)):
        if S1[x-1] == S2[y-1]:
            M[x][y] = M[x-1][y-1] + 1
            if M[x][y]>longest:
                longest = M[x][y]
                x_longest  = x
        else:
            M[x][y] = 0
  return S1[x_longest-longest: x_longest]
Другие вопросы по тегам