Извлечь последовательность fasta, фланкирующую аминокислоту

Я пытаюсь придумать сценарий python для извлечения последовательности из 12 аминокислот, фланкирующей данную аминокислоту (по 6 в каждом направлении) последовательности fasta.

Вход

У меня есть 2 входа: файл фаста и фрейм данных панды.

Файл fasta выглядит так:

> sp|P00001| some text here 1
MKLLILTCLVAVALARPKHPIKKVSPTFDTNMVGKHQGLPQEVLNENLLRFFVAPFPEVFGKEKVSLDAGPGMCSRNE
>sp|P00002| some text here 2
MSSGNAKIGHPAPNFKATAVMPDGQFKDISLSDYKGKYVVFFFYPLDFTFVCPTGLGRSSYRATSCLPALCLP
>sp|P00003| some text here 3
MSVLDSGNFSWKMTEACMKVKIPLVKKKSLRQNLIENGKLKEFMRTHKYNLGSKYIREAATLVSEQPLQN

Вот мой второй ввод, фрейм данных pandas (2 столбца «ProteinID» и «Phosphopeptide»)

ProteinID  -- Phosphopeptide
P00001   --   KVSPT*FDTNMVGK
P00001  --    SLDAGPGMCS*R
P00003   --   LDS*GNFSWKMTEACMK

Цель

Что мне нужно сделать, так это следующее. Для каждого «фосфопептида» мне нужно найти белок (ProteinID) в заголовке файла fasta (начинается с «>»). Затем мне нужно извлечь 6 аминокислот, предшествующих и следующих за аминокислотой, отмеченной звездочкой.

Вывод

Мой вывод - это новый столбец, записанный во фрейм данных, и он будет выглядеть так:

ProteinID  -- Phosphopeptide  --   NewColumn
P00001   --   KVSPT*FDTNMVGK  --   IKKVSPTFDTNMV 
P00001   --   SLDAGPGMCS*R    --   AGPGMCSRNE
P00003   --   LDS*GNFSWKMTEACMK -- MSVLDSGNFSWK

Обратите внимание, что последние 2 ряда содержат пептиды в конце или в начале соответствующих белков, поэтому у нас нет 12 аминокислот для извлечения в этих случаях.

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


person user3386437    schedule 06.03.2014    source источник
comment
Немного сбивает с толку, можете ли вы написать так: это ввод, а это ожидаемый результат. Так что любой, кто не знает структуру белка, легко поймет. для меня все струны.   -  person James    schedule 06.03.2014
comment
где заголовок файла fasta (начинается с «› »)?   -  person James    schedule 06.03.2014
comment
извините, я не понял, что фаста была деформирована. Хорошо, теперь это исправлено. Благодарность   -  person user3386437    schedule 07.03.2014
comment
Ознакомьтесь с справкой по редактированию, чтобы узнать, как правильно форматировать текст. Что касается самого вопроса: покажите нам, что вы пробовали и каковы ваши конкретные проблемы.   -  person Konrad Rudolph    schedule 07.03.2014


Ответы (2)


Вот функция, которая извлекает соответствующую подстроку:

def flank(seq, pp):
    # 1: find the position of the AA preceding the '*' marker in the
    # phosphopeptide
    marked_pos = pp.find('*') - 1
    if (marked_pos < 0):
        raise ValueError("invalid phosphopeptide string")

    # 2: find the phosphopeptide (without '*') in the sequence
    pp_pos = s.find(pp.replace('*', ''))
    if pp_pos == -1:
        raise ValueError("phosphopeptide not found in the sequence")

    # avoid a negative starting index
    start = max(0, pp_pos + marked_pos - 6)

    # 3: use slicing to produce the result
    return seq[start : pp_pos + marked_pos + 7]

Пример:

seq = "MKLLILTCLVAVALARPKHPIKKVSPTFDTNMVGKHQGLPQEVLNENLLRFFVAPFPEVFGKEKVSLDAGPGMCSRNE"
pp = "KVSPT*FDTNMVGK"
print(flank(seq, pp))

печатает:

IKKVSPTFDTNMV
person Brave Sir Robin    schedule 06.03.2014
comment
Большое спасибо вам обоим, rmartinjak и sjcipher за вашу помощь. Ура, Алекс - person user3386437; 08.03.2014

Привет, проверьте это: мой файл fasta называется txt:

Фрагмент:

#!/usr/bin/python
import re

protein_dict = [
    ('P00001', 'KVSPT*FDTNMVGK'),
    ('P00001', 'SLDAGPGMCS*R'),
    ('P00003', 'LDS*GNFSWKMTEACMK')
    ]

protein_id = None

def prepare_structure_from_fasta(file):
    fasta_structure = dict()
    with open(file, 'r') as fh:
        for line in fh:
            if '>' in line:
                protein_id = line.split('|')[1]
            else:
                if not protein_id:
                    raise Exception("Wrong fasta file structure")
                fasta_structure[protein_id] = line.strip()
    return fasta_structure


def match(pattern, string):
    matc = re.search(pattern, string)
    if matc:
        return matc.groups()[0]
    return None

fasta_struct = prepare_structure_from_fasta('txt')
final_struct = []

for pro_d in protein_dict:

    pro_id = pro_d[0]
    pep_id = pro_d[1]
    first, second = pep_id.split('*')

    if len(first) <= 6:
        f_count = 7 - len(first)
    else:
        first = first[len(first) - 7:]
        f_count = 0
    if len(second) <= 6:
        s_count = 7 - len(second)
    else:
        second = second[0:6]
        s_count = 0

    _regex = '([A-Z]{0,%d}%s%s[A-Z]{0,%d})' % (f_count,first,second,s_count)
    final_struct.append((pro_id, pep_id, match(_regex, fasta_struct[pro_id])))

for pro in final_struct:
    print pro

Вывод:

('P00001', 'KVSPT*FDTNMVGK', 'IKKVSPTFDTNMV')
('P00001', 'SLDAGPGMCS*R', 'AGPGMCSRNE')
('P00003', 'LDS*GNFSWKMTEACMK', 'MSVLDSGNFSWK')
person James    schedule 07.03.2014