Консолидация аналогичных шаблонов в единый шаблон консенсуса

В предыдущем посте я не разъяснил вопросы должным образом, поэтому я хотел бы начать новую тему здесь.

У меня есть следующие предметы:

  1. отсортированный список из 59 000 шаблонов белка (от 3 символов "FFK" до 152 символов);

  2. некоторые длинные белковые последовательности, иначе моя ссылка.

Я собираюсь сопоставить эти шаблоны с моей ссылкой и найти место, где найдено совпадение. (Мой друг помог написать сценарий для этого.)

import sys
import re
from itertools import chain, izip

# Read input
with open(sys.argv[1], 'r') as f:
    sequences = f.read().splitlines()
with open(sys.argv[2], 'r') as g:
    patterns = g.read().splitlines()

# Write output
with open(sys.argv[3], 'w') as outputFile:
    data_iter = iter(sequences)
    order = ['antibody name', 'epitope sequence', 'start', 'end', 'length']
    header = '\t'.join([k for k in order])
    outputFile.write(header + '\n')
    for seq_name, seq in izip(data_iter, data_iter):
        locations = [[{'antibody name': seq_name, 'epitope sequence': pattern, 'start': match.start()+1, 'end': match.end(), 'length': len(pattern)} for match in re.finditer(pattern, seq)] for pattern in patterns]
        for loc in chain.from_iterable(locations):
            output = '\t'.join([str(loc[k]) for k in order])
            outputFile.write(output + '\n')

f.close()
g.close()
outputFile.close()

Проблема в том, что в этих 59 000 шаблонах после сортировки я обнаружил, что некоторая часть одного шаблона совпадает с частью других шаблонов, и я хотел бы объединить их в один большой шаблон "консенсуса" и просто сохранить консенсус (см. Примеры ниже).):

TLYLQMNSLRAED
TLYLQMNSLRAEDT
  YLQMNSLRAED
  YLQMNSLRAEDT
  YLQMNSLRAEDTA
  YLQMNSLRAEDTAV

даст

TLYLQMNSLRAEDTAV

другой пример:

    APRLLIYGASS
    APRLLIYGASSR
    APRLLIYGASSRA
    APRLLIYGASSRAT
    APRLLIYGASSRATG
    APRLLIYGASSRATGIP
    APRLLIYGASSRATGIPD
  GQAPRLLIY
KPGQAPRLLIYGASSR
KPGQAPRLLIYGASSRAT
KPGQAPRLLIYGASSRATG
KPGQAPRLLIYGASSRATGIPD
       LLIYGASSRATG
       LLIYGASSRATGIPD
   QAPRLLIYGASSR

даст

KPGQAPRLLIYGASSRATGIPD

PS: я выравниваю их здесь, чтобы их было проще визуализировать. Изначально 59 000 шаблонов не отсортированы, поэтому трудно увидеть консенсус в реальном файле.

В моей конкретной проблеме я не выбираю самые длинные модели, вместо этого мне нужно принимать во внимание каждую модель, чтобы найти консенсус. Я надеюсь, что объяснил достаточно ясно для моей конкретной проблемы.

Спасибо!

Источник

2 ответа

Вот мое решение с рандомизированным порядком ввода для повышения достоверности теста.

import re
import random

data_values = """TLYLQMNSLRAED
TLYLQMNSLRAEDT
  YLQMNSLRAED
  YLQMNSLRAEDT
  YLQMNSLRAEDTA
  YLQMNSLRAEDTAV
    APRLLIYGASS
    APRLLIYGASSR
    APRLLIYGASSRA
    APRLLIYGASSRAT
    APRLLIYGASSRATG
    APRLLIYGASSRATGIP
    APRLLIYGASSRATGIPD
  GQAPRLLIY
KPGQAPRLLIYGASSR
KPGQAPRLLIYGASSRAT
KPGQAPRLLIYGASSRATG
KPGQAPRLLIYGASSRATGIPD
       LLIYGASSRATG
       LLIYGASSRATGIPD
   QAPRLLIYGASSR"""

test_li1 = data_values.split()
#print(test_li1)

test_li2 = ["abcdefghi", "defghijklmn", "hijklmnopq", "mnopqrst", "pqrstuvwxyz"]

def aggregate_str(data_li):
    copy_data_li = data_li[:]
    while len(copy_data_li) > 0:
        remove_li = []
        len_remove_li = len(remove_li)
        longest_str = max(copy_data_li, key=len)
        copy_data_li.remove(longest_str)
        remove_li.append(longest_str)
        while len_remove_li != len(remove_li):
            len_remove_li = len(remove_li)
            for value in copy_data_li:
                value_pattern = "".join([x+"?" for x in value])
                longest_match = max(re.findall(value_pattern, longest_str), key=len)
                if longest_match in value:
                    longest_str_index = longest_str.index(longest_match)
                    value_index = value.index(longest_match)
                    if value_index > longest_str_index and longest_str_index > 0:
                        longest_str = value[:value_index] + longest_str
                        copy_data_li.remove(value)
                        remove_li.append(value)
                    elif value_index < longest_str_index and longest_str_index + len(longest_match) == len(longest_str):
                        longest_str += value[len(longest_str)-longest_str_index:]
                        copy_data_li.remove(value)
                        remove_li.append(value)
                    elif value in longest_str:
                        copy_data_li.remove(value)
                        remove_li.append(value)
        print(longest_str)
        print(remove_li)

random.shuffle(test_li1)
random.shuffle(test_li2)
aggregate_str(test_li1)
#aggregate_str(test_li2)

Вывод из print ().

KPGQAPRLLIYGASSRATGIPD
['KPGQAPRLLIYGASSRATGIPD', 'APRLLIYGASS', 'KPGQAPRLLIYGASSR', 'APRLLIYGASSRAT', 'APRLLIYGASSR', 'APRLLIYGASSRA', 'GQAPRLLIY', 'APRLLIYGASSRATGIPD', 'APRLLIYGASSRATG', 'QAPRLLIYGASSR', 'LLIYGASSRATG', 'KPGQAPRLLIYGASSRATG', 'KPGQAPRLLIYGASSRAT', 'LLIYGASSRATGIPD', 'APRLLIYGASSRATGIP']
TLYLQMNSLRAEDTAV
['YLQMNSLRAEDTAV', 'TLYLQMNSLRAED', 'TLYLQMNSLRAEDT', 'YLQMNSLRAED', 'YLQMNSLRAEDTA', 'YLQMNSLRAEDT']

Edit1 - краткое объяснение кода.

1.) Найти самую длинную строку в списке

2.) Просмотрите все оставшиеся строки и найдите максимально возможное совпадение.

3.) Убедитесь, что совпадение не является ложным срабатыванием. Основываясь на том, как я написал этот код, следует избегать сопряжения одиночных перекрытий на концах терминала.

4.) При необходимости добавьте совпадение к самой длинной строке.

5.) Если больше ничего нельзя добавить к самой длинной строке, повторите процесс (1-4) для следующей самой длинной оставшейся строки.

Edit2 - исправлено нежелательное поведение при обработке таких данных, как ["abcdefghijklmn", "ghijklmZopqrstuv"]

Источник 
def main():
#patterns = ["TLYLQMNSLRAED","TLYLQMNSLRAEDT","YLQMNSLRAED","YLQMNSLRAEDT","YLQMNSLRAEDTA","YLQMNSLRAEDTAV"]
patterns = ["APRLLIYGASS","APRLLIYGASSR","APRLLIYGASSRA","APRLLIYGASSRAT","APRLLIYGASSRATG","APRLLIYGASSRATGIP","APRLLIYGASSRATGIPD","GQAPRLLIY","KPGQAPRLLIYGASSR","KPGQAPRLLIYGASSRAT","KPGQAPRLLIYGASSRATG","KPGQAPRLLIYGASSRATGIPD","LLIYGASSRATG","LLIYGASSRATGIPD","QAPRLLIYGASSR"]

test = find_core(patterns)
test = find_pre_and_post(test, patterns)
#final = "YLQMNSLRAED"
final = "KPGQAPRLLIYGASSRATGIPD"

if test == final:
    print("worked:" + test)
else:
    print("fail:"+ test)

def find_pre_and_post(core, patterns):
    pre = ""
    post = ""
    for pattern in patterns:
        start_index = pattern.find(core)
        if len(pattern[0:start_index]) > len(pre):
            pre = pattern[0:start_index]
        if len(pattern[start_index+len(core):len(pattern)]) > len(post):
            post = pattern[start_index+len(core):len(pattern)]

    return pre+core+post

def find_core(patterns):
    test = ""
    for i in range(len(patterns)):
        for j in range(2,len(patterns[i])):
            patterncount = 0
            for pattern in patterns:
                if patterns[i][0:j] in pattern:
                    patterncount += 1
            if patterncount == len(patterns):
                test = patterns[i][0:j]
    return test

main()

  1. Итак, что я делаю в первую очередь, это найти основное ядро ​​в find_core Функция начинается со строки длины два, так как одного символа недостаточно для первой строки. Затем я сравниваю эту подстроку и вижу, находится ли она во ВСЕХ строках как определение "ядра".

  2. Затем я нахожу индексы подстроки в каждой строке, чтобы потом найти подстроки pre и post, добавленные в ядро. Я отслеживаю эти длины и обновляю их, если одна длина больше другой. У меня не было времени, чтобы изучить крайние случаи, так что вот мой первый выстрел

Источник
Другие вопросы по тегам