Сообщение об ошибке вызвано действиями ранее в коде, а именно:
entry = iterator.__next__
Что должно было быть:
entry = iterator.__next__()
Или в Python 3:
entry = next(iterator)
Моя переработка вашего кода:
from Bio import SeqIO
def batch_iterator(iterator, batch_size):
entry = True # Make sure we loop once
while entry:
batch = []
while len(batch) < batch_size:
try:
entry = next(iterator)
except StopIteration:
entry = False
if not entry:
# End of file
break
batch.append(entry)
if batch:
yield batch
record_iter = SeqIO.parse('/Users/nermin/mainfolder/MTB_NITR203.fasta', 'fasta')
for i, batch in enumerate(batch_iterator(record_iter, 1000), start=1):
filename = 'group_{}.fasta'.format(i)
count = SeqIO.write(batch, filename, 'fasta')
print('Wrote {} records to {}'.format(count, filename))
Мне нужен простой скрипт, чтобы разделить массивную последовательность ДНК на меньшие .fasta для поиска BLAST.
Здесь у меня есть вопрос по вашему вопросу. Из вашего описания я ожидал, что будет входной файл вроде:
>ZB3243.4 Platypus mRNA for dDFD-w, complete cds.
TTTTAATTTTGCTTTCAATATGACGGCTGTCAATGTTGCCCTGATTCGTGATACCAAGTG
GCTGACTTTAGAAGTCTGTAGAGAATTTCAGAGAGGAACTTGCTCTCGAGCTGATGCAGA
TTGCAAGTTTGCCCATCCACCAAGAGTTTGCCATGTGGAAAATGGTCGTGTGGTGGCCTG
... extremely long sequence
GCTAGGACAGGGAACAGGGAAGCACTTACAATTATTCCTTGATTTATTCAAAAGAACTGG
GAAAGATGGTTGTAGTTGTCTTTAGCTTCGGTTCAACTGAGTTTCGTTTTGTTAAACAGT
TCAGACCCTCTCACATCATAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
Но то, что обрабатывает ваш код, на самом деле таково:
>ZB3243.4 Platypus mRNA for dDFD-w, complete cds.
TTTTAATTTTGCTTTCAATATGACGGCTGTCAATGTTGCCCTGATTCGTGATACCAAGTG
GCTGACTTTAGAAGTCTGTAGAGAATTTCAGAGAGGAACTTGCTCTCGAGCTGATGCAGA
TTGCAAGTTTGCCCATCCACCAAGAGTTTGCCATGTGGAAAATGGTCGTGTGGTGGCCTG
>VF42354.1 Rhino Ig active H-chain V-region, subgroup VH-II
TTCATGCAAATATGCTCTCTTTCTTTAGAATATTTCTGTAGGTTTCTTGGGACTGACATT
TAAAACGCCTCACTTTTGAATGTGCACAAAACCTGCTCATTAACATGCATGTGTATAATT
TGTACCTGCAGATCTGATGTTGCATAATACAATCAAATTACTAGATTTTTTAAAGAGAGA
>GS45345.54 Aardvark binding protein.
AACAACACCTGCCACCAGCGTTCCGTTCGCTGCACCAACTACAGGCAATCAGCTGAAATT
CTGAACAGCAGAGTTATGGAGTATCAGAATCTTTCCATGGAAACCTCCATATGGCCTTTC
TATATATATTCTCGTATGTCTTATTCTACCAACACAACAATAAGCGTGTTGCAGTCAATG
... extremely large number of sequences
>GR343245.2 Eggplant subgroup VH-II, mRNA.
TTGCCGCTATGCTCACCCTACTGATGCTTCCATGATTGAAGCGAGTGATAATACTGTGAC
AATCTGCATGGATTACATCAAAGGTCGATGCTCGCGGGAGAAATGCAAGTACTTTCATCC
TCCTGCACACTTGCAAGCCAGACTCAAGGCAGCTCATCATCAGATGAACCATTCAGCTGC
>FG345252.3 Bedbug binding protein 4
TTTTGATTCTCTAAAGGGTCGGTGTACCCGAGAGAACTGCAAGTACCTTCACCCTCCTCC
ACACTTAAAAACGCAGCTGGAGATTAATGGGCGGAACAATCTGATTCAACAGAAGACTGC
CGCAGCCATGTTCGCCCAGCAGATGCAGCTTATGCTCCAAAACGCTCAAATGTCATCACT
>MD2435324.5 Mantis subgroup VH-II, mRNA.
TAAAAATATGCTAATTACAAGTTATAAATCAAACGGAGAGATGGGGGCATGGAGATAGTT
TTTACGTACTGGAGGAAAGTGTGTAAAACCATGGCAATGTCACCTTTTACACAAATGCCA
TTTTCCAAATGCAAATGGCTCATGCTCTTTAGACTACTCTTTGAATAACAAGTAAGATGC
Вы обрабатываете огромное количество последовательностей или одну большую последовательность?
Входная фаста - это всего лишь 1 массивная последовательность с примерно 4 миллионами базовых пар (точно так же, как в вашем первом примере), но из-за ее размера мы не можем выполнять быстрый поиск. Единственное, что нам нужно, это разделить его на более мелкие файлы fasta, желательно менее 1 мб каждый
Мы установили, что ваш код не предназначен для того, чтобы делать то, что вы хотите. Итак, давайте атакуем вашу проблему, просматривая данные одной последовательности, разбивая ее на отдельные файлы FASTA:
from Bio import SeqIO
from Bio.SeqRecord import SeqRecord
FILENAME = '/Users/nermin/mainfolder/MTB_NITR203.fasta'
BPS_PER_FILE = 500000 # base pairs per FASTA file
# There should be one and only one record, the entire genome:
large_record = SeqIO.read(FILENAME, 'fasta')
large_sequence = large_record.seq
large_description = large_record.description
limit = len(large_sequence)
for i, start in enumerate(range(0, limit, BPS_PER_FILE), start=1):
stop = min(start + BPS_PER_FILE, limit)
small_sequence = large_sequence[start:stop]
small_description = large_description + "; base pairs {} to {}".format(start, stop)
small_record = SeqRecord(small_sequence, id=large_record.id, description=small_description)
filename = 'group_{}.fasta'.format(i)
count = SeqIO.write(small_record, filename, "fasta")
assert (count == 1), "Incorrect record count!"
print('Wrote {} base pairs to {}'.format(stop - start, filename))
Каждый выходной файл представляет собой отдельную последовательность, как и исходный, но меньшего размера, с измененным описанием, которое теперь включает диапазон базовых пар в каждом конкретном файле.
person
cdlane
schedule
20.06.2017
