Рисование нескольких последовательностей из одного файла на основе общих полей в другом файле - PullRequest
Новая
       73

Рисование нескольких последовательностей из одного файла на основе общих полей в другом файле

0 голосов
/ 11 декабря 2018

Я пытаюсь запустить скрипт Python для рисования последовательностей из отдельного файла (merged.fas), относительно списка (gene_fams_eggnog.txt), который я имею в качестве вывода из другой программы.

код выглядит следующим образом: 

from Bio import SeqIO
import os, sys, re
from collections import defaultdict
sequences = "merged.fas"
all_seqs = SeqIO.index(sequences, "fasta")
gene_fams = defaultdict(list)

gene_fams_file = open("gene_fams_eggnog.txt")
for line in gene_fams_file:
     fields = re.split("\t", line.rstrip())
     gene_fams[fields[0]].append[fields[1]]

for fam in gene_fams.keys():
    output_filename = str(fam) + ".fasta"
    outh = open(output_filename, "w")
    for id in gene_fams[fam]:
        if id in all_seqs:
            outh.write(">" + all_seqs[id].description + "\n" + str(all_seqs[id].seq) + "\n")
        else:
            print "Uh oh! Sequence with ID " + str(id) + " is not in the all_seqs file!"
            quit()
    outh.close()

Список выглядит следующим образом:

1 Saccharomycescerevisiae_DAA09367.1
1 bieneu_EED42827.1
1 Asp_XP_749186.1
1 Mag_XP_003717339.1
2 Mag_XP_003716586.1
2 Mag_XP_003709453.1
3 Asp_XP_749329.1

Поле 0 обозначает группировку, основанную на сходстве между последовательностями.Сценарий должен был взять все последовательности из merged.fas, которые соответствуют коду в поле 1, и записать их в файловую базу в поле 0.

Так что в случае части списка Iпоказали, что все последовательности, которые имеют 1 в поле 0 (Saccharomycescerevisiae_DAA09367.1, bieneu_EED42827.1, Asp_XP_749186.1, Mag_XP_003717339.1), были бы записаны в файл с именем 1.fasta.Это должно продолжаться с 2.fasta - сколько бы групп ни было.

Так что это сработало, однако оно не включает в себя все последовательности, которые есть в группе, оно будет включать только последнююуказан как часть этой группы.Используя мой пример выше, у меня был бы только файл (1.fasta) с одной последовательностью (Mag_XP_003717339.1) вместо всех четырех.

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

1 Ответ

0 голосов
/ 31 декабря 2018

Хотя я не определил причину проблемы, на которую вы жаловались, я удивлен, что ваш код вообще работает с этой ошибкой: 

gene_fams[fields[0]].append[fields[1]]

, то есть append[...] вместо append(...).Но, возможно, это также, «не там, в самом скрипте, который я запускаю».Я переписал ваш сценарий ниже, и он отлично работает для меня.Одной из проблем было использование вами имени переменной id, которая является встроенной в Python.Вы увидите, я иду в крайность, чтобы избежать таких ошибок: 

from Bio import SeqIO
from collections import defaultdict

SEQUENCE_FILE_NAME = "merged.fas"
FAMILY_FILE_NAME = "gene_families_eggnog.txt"

all_sequences = SeqIO.index(SEQUENCE_FILE_NAME, "fasta")
gene_families = defaultdict(list)

with open(FAMILY_FILE_NAME) as gene_families_file:
    for line in gene_families_file:
        family_id, gene_id = line.rstrip().split()
        gene_families[family_id].append(gene_id)

for family_id, gene_ids in gene_families.items():
    output_filename = family_id + ".fasta"

    with open(output_filename, "w") as output:
        for gene_id in gene_ids:
            assert gene_id in all_sequences, "Sequence {} is not in {}!".format(gene_id, SEQUENCE_FILE_NAME)

            output.write(all_sequences[gene_id].format("fasta"))
...