Заменить последовательность FASTA в соответствии с ее идентификатором пробелами - PullRequest
1 голос
/ 22 апреля 2020

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

 >T1        
AQSFDRATVFEARRAGYQRESRVALGKSTGVLEWHVFHVWAPRETTILVETLSQIECSGKGIADRRQENPLTI------ATAI--TSQLLELVDIVIMDFKAITQFFL     
>T484      
AQSFDRATVFEARRAGYQREARVALGKSTGKLEWQVFHVWAPRETTILVETLSQLECSGKGIADRRQENPLKI------ATAI--TSQLLELVDIVIMDFKAITQFFL     
>T582      
AQSFDRATVFEKRRAGYQREARVALGKSTGKLEWQVFHVWAPRETTILVETLSQLECSGKGIADRRQENPLKI------ATAI--TSQLLELVDIVIMDFKAITQFFL     
>T1424     
AQSFDRATVFEKRRAGYQLEARVALGKSTGKLEWQVFHVWAPRETTILVETLSQLECSGKGIANRRQENPLKI------ATAI--TSQLLELVDILIMDFKAITQFFL     
>T1552     
AQSFDRATIFEKRRAGYQIEARVALGKSTGKLEWQVYHAWAPRETTILVETLSQLENAGKGVANRRHENPLKI------ATAI--TSQLLELVDILMMDFKAITQFFL     
.
.
.

I wi sh, чтобы случайным образом нарисовать последовательность f из N (N = 2000 последовательностей).

Для Например, если f = 2, я случайным образом рисую 2 последовательности из 2000 последовательностей.

f=2 

l=[] 

for i in range(f):
    x=randint(1, N)
    l.append(x)

В моем списке l я бы, например, [291, 566]. Затем я хочу нарисовать 291-ю и 566-ю последовательности:

> T1424
AQSFDRATVFEKRRAGYQLEARVALGKSTGKLEWQVFHVWAPRETTILVETLSQLECSGKGIANRRQENPLKI ------ ATAI - TSQLLELVDILIMDFKAITQFFL
> T1552
AQSFDRATIFEKRRAGYQIEARVALGKSTGKLEWQVYHAWAPRETTILVETLSQLENAGKGVANRRHENPLKI ------ ATAI - TSQLLELVDILMMDFKAITQFFL

Я хочу заменить эти последовательности на пробелы ("-") длиной 55:

> T1424
-----------------------------------------------------------------------------------------------------------------------------
> T1552
-----------------------------------------------------------------------------------------------------------------------------

Я пытался с этим кодом:

from Bio import SeqIO 
from random import *

records = list(SeqIO.parse("fichier1.txt", "fasta")) 
#print(records[0].id)  # first record
#print(records[0].seq)

N=len(records) #2000 

f=2 

l=[] 

for i in range(f):
    x=randint(1, N)
    l.append(x)

d={} 
for i2 in l:
    print(records[i2].id,records[i2].seq)
    d[str(records[i2].id)]=str(records[i2].seq)
    with open("fichier1.txt") as fichier, open("newfile.txt", "w") as newfile:
        texte = fichier.read()

        new_text += texte.replace(str(records[i2].seq), "---------------------------------------")
        newfile.write(new_text)

print(d)

Что не работает, потому что иногда может быть одна и та же последовательность в файле, но с другим идентификатором.

Из идентификатора и соответствующей последовательности, Я хотел бы изменить последовательность, чтобы ввести пробелы.

1 Ответ

2 голосов
/ 22 апреля 2020

Я могу решить вашу проблему кодирования (хотя я до сих пор не убежден, что это лучший способ решить вашу научную проблему c - то есть отсутствие определенных рибсомальных белков). В любом случае, здесь есть решение, вместо того чтобы полагаться на уникальность последовательностей, я просто отслеживал индексы последовательностей:

import random

from Bio import SeqIO
from Bio.Seq import Seq

num_to_sample = 2

# WARNING: puts all sequences in memory
lst = list(SeqIO.parse("fichier1.txt", "fasta"))

sample = set(random.sample(range(len(lst)), num_to_sample))

for idx, record in enumerate(lst):
    if idx in sample:
        record.seq = Seq(len(record.seq) * "-")
    print(record.format("fasta"), end='')

Пример вывода:

>T1
AQSFDRATVFEARRAGYQRESRVALGKSTGVLEWHVFHVWAPRETTILVETLSQIECSGK
GIADRRQENPLTI------ATAI--TSQLLELVDIVIMDFKAITQFFL
>T484
AQSFDRATVFEARRAGYQREARVALGKSTGKLEWQVFHVWAPRETTILVETLSQLECSGK
GIADRRQENPLKI------ATAI--TSQLLELVDIVIMDFKAITQFFL
>T582
------------------------------------------------------------
------------------------------------------------
>T1424
------------------------------------------------------------
------------------------------------------------
>T1552
AQSFDRATIFEKRRAGYQIEARVALGKSTGKLEWQVYHAWAPRETTILVETLSQLENAGK
GVANRRHENPLKI------ATAI--TSQLLELVDILMMDFKAITQFFL
...