Извлечь последовательность ДНК по функции из файла genbank
У меня есть файлы genbank, состоящие из нескольких аннотированных контигов. То, что я хотел бы сделать, это разделить это на базу данных, в которой есть отдельные записи генов, содержащие каждую особенность CDS, а также ее последовательность ДНК и аминокислотную последовательность. Пока это работает просто отлично:
for record in SeqIO.parse(open_file, 'gb'):
for feature in record.features:
if feature.type == 'CDS':
gene_record = {
'locus_tag':feature.qualifiers['locus_tag'][0],
'translation':feature.qualifiers['translation'][0],
}
Проблема в том, что я получаю последовательность ДНК. Файл genbank отформатирован так:
FEATURES Location/Qualifiers
source 1..29869
/organism="Arthrobacter"
/mol_type="genomic DNA"
/strain="strain_name"
gene complement(4..462)
/locus_tag="ArthroDRAFT_00001"
CDS complement(4..462)
/locus_tag="ArthroDRAFT_00001"
/product="hypothetical protein"
/translation="LSTGKELLNYQSALNDIHDEFSRAQQSDAGVSHLSVAKITEKLS
YLKATALQMDDLFSVLRKQGVSLRSTGLADWASVPTIQDDKEEGKTEPSLAKKEISSR
TTSKPNKIDFPKFEYPDHGQPTNKIRVGTILDTFSESAFSYEWINVALQD"
gene complement(1126..1842)
/locus_tag="ArthroDRAFT_00002"
CDS complement(1126..1842)
/locus_tag="ArthroDRAFT_00002"
/product="hypothetical protein"
/translation="VPRAFIYGSCVGGDTANVFPSDWDRPTYVARQSIISAAFGPTSV
EGDIELTSAFQRSMLEGDIEATAFPRLRQELPTHDVLILDIVDERLGVYELAPGKYLT
RSMELISSKLIGKQPVTPRLIEFGSDEHYGLWTRSVDMLVDVVKHGGIPVFALLPPWS
EKSIQGEDLTWHSVSVDLMNNKYARYNEYLVQSEFTVVTVPDEEALGDAEHKWGLAPF
HYTESVYESLRDQILVGVSS"
...etc
ORIGIN
1 ccctcaatcc tgaagagcca cattgatcca ttcgtatgag aatgcagatt ccgagaatgt
61 atcaagaatt gttcctacac gtattttgtt tgtcggctgg ccgtggtccg gatactcaaa
121 ttttggaaag tcaatcttgt ...
Итак, у меня есть местоположение для каждой функции, но мне нужно как-то разобрать последовательность ДНК, чтобы получить правильную последовательность ДНК, и я вроде как в растерянности, как это сделать. feature.location
дает мне то, что кажется полезным выводом:
[3:462](-)
[1125:1842](-)
[2159:3755](-)
[5190:5532](-)
[6226:6493](+)
Но я знаю, что использовать это можно только с помощью регулярных выражений, а затем record.seq[start_num:finish_num]
, А потом, если он на - прядь, беги через что бы то ни было reverse_complement
,
Это кажется слишком запутанным, и я подумал, что у биопиона должен быть более эффективный способ, но я не могу его найти. Я заметил, что feature.location
уже кажется, что 1 вычтен из него для учета индекса 0 в python, так что это число должно быть пригодным для использования... верно?
РЕДАКТИРОВАТЬ: SeqFeature.location
кажется, что это то, что я хочу использовать, но не могу понять, как: http://biopython.org/DIST/docs/api/Bio.SeqFeature.FeatureLocation-class.html
EDIT2: похоже, я могу сделать
record.seq[feature.location.start:feature.location.end]
Это еще не работает для негативных вещей, хотя...
2 ответа
Я использую "извлечение" метода определения местоположения объекта.
dna_seq = feature.location.extract (record)
Это отлично работает.
Хорошо, не берите в голову - больше прибегая к помощи и экспериментирования привели меня к лучшему решению. Не уверен, что это лучшее решение, я оставлю вопрос без ответа, если у кого-нибудь появятся лучшие идеи. Что я сделал, чтобы получить последовательность ДНК:
if feature.location.strand == 1:
dna_seq = record.seq[
feature.location.start:feature.location.end
]
elif feature.location.strand == -1:
dna_seq = record.seq[
feature.location.start:feature.location.end
].reverse_complement()
Не самое красивое решение, но, похоже, работает...