Как разбить файл fasta

Этот код предназначен для извлечения и разделения последовательностей из файла fasta.

outfile=open('outf','w')
for line in open('input'):
      if line[0]==">":
         outfile.write('\n')
     else:
    outfile.write(line.strip())

   outfile.close()
 all_codons=[]
 for line in open('outf', 'r'):
     seq=line.strip()
     codons = [seq[i:i+3] for i in xrange(0, len(seq), 3) if len(seq[i:i+3])==3]
     all_codons.append(codons)

Затем из разделенных последовательностей я хочу взять три последовательности, длина которых составляет 9 (9 оснований), например:

    CGTAACAAG 
    AATCCGGAG 
    CCGCCTCGG

Я разбиваю первую последовательность на 3 подпоследовательности по три основания, так что из одной последовательности получаю 3 подпоследовательности, то же самое делаю для двух других последовательностей.

Так:

    CGT    AAC     AAG 
    AAT    CCG     GAG 
    CCG    CCT     CGG

пример:

identical_segment('CGT')

Я хочу применить эту функцию к каждой подпоследовательности трех последовательностей, а затем применить то же самое ко всем файлам fasta. Итак, цель состоит в том, чтобы получить матрицу, например, я беру первую подпоследовательность «CGT» и применяю функцию идентичный_сегмент(), она возвращает 28, то же самое для остальных 8 подпоследовательностей. Итак, я получаю матрицу (3,3):

28         2             3
4          23            35
23         4             27

Что я могу сделать?


python split fasta
person m28    schedule 21.07.2013    source источник
comment
вы не должны публиковать новый вопрос, чтобы уточнить старый! Просто нажмите кнопку «Изменить» в своем старом вопросе и добавьте туда информацию.   -  person seth    schedule 21.07.2013

Ответы (2)


arrow_upward
0
arrow_downward

В вашем коде есть некоторые проблемы.

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

def processLines(inputname):
    all_codons=[]
    for line in open(inputname):
        if line[0]==">":
            seq=line.strip()
            codons = [seq[i:i+3] for i in xrange(0, len(seq), 3) if
                      len(seq[i:i+3])==3]
        all_codons.append(codons)
    return all_codons

Кроме того, каждый вызов идентичных_сегментов будет генерировать словарь, который вы используете в качестве сопоставления от str к scores. Это может стать дорогим, когда число вызовов масштабируется. Чтобы этого избежать, можно попробовать два способа:

code={"a":0,"c":1,"g":2,"t":3} 
def identical_segment(input_string):
   .... # what you have written

или создайте класс, экземпляр которого содержит словарь.

Чтобы обработать несколько файлов, выполните:

output = [processLines(filename) for filename in filenames]
# filenames is an iterable

или если вы хотите сопоставить входное имя с выходным:

outputDict = {filename: processLines(filename) for 
              filename in filenames}

В конце концов, вызывайте свою анализирующую функцию для каждого выхода и записывайте их в выходной файл.

Подводя итог тому, что вы должны подобрать в этом посте:

  1. Выходные файлы могут быть не лучшим вариантом, так как файловый ввод-вывод стоит дорого. Если вы записываете его в какой-то файл, значит, вам придется его снова читать, что вдвойне дорого.

  2. Один и тот же объект не должен создаваться снова и снова. Проверьте свой код, чтобы этого не произошло.

  3. Разделите основную задачу на несколько небольших задач, а затем придумайте простой и интуитивно понятный способ начать выполнение каждой задачи. В этом примере у нас есть processfiles->analyse->output_result

  4. Понимание — это полезный способ повторения вещей в Python, и он более удобочитаем. Вы можете выполнить поиск по списку и словарю, чтобы узнать больше.

Попробуйте что-нибудь сами. Я буду более чем счастлив прочитать ваш улучшенный код здесь.

person Mai    schedule 21.07.2013

arrow_upward
0
arrow_downward

Попробуйте использовать BioPython, чтобы извлечь последовательности нуклеотидов из файла fasta. Используя этот пакет,

from Bio import AlignIO

for record in AlignIO.parse('filename.fasta', 'fasta'):
    print record.id, record.seq

# or store in a new file
seqs = []
for record in AlignIO.parse('filename.fasta', 'fasta'):
    seqs.append(record.seq + '\n')

with open(outfile, 'w') as out:
    out.writelines(seqs)
person wflynny    schedule 23.09.2013
comment
Это не работает, когда каждая последовательность имеет разную длину, например ValueError: Sequences must all be the same length - person Shadi; 25.05.2020