Предсказание генов прокариот

Информация о плазмиде

Моей задачей было предсказать гены в плазмиде CP015437. Она пренадлежит бактерии Anoxybacillus sp. B7M1 со следующем положением в систематике: Bacteria, Firmicutes, Bacilli, Bacillales, Bacillaceae, Anoxybacillus. Бактерии данного рода являются палочковидными и спорообразующими. Они встречаются в геотермальных источниках, навозе и на заводах по переработке молока.
Плазмида состоит из 54,976 пн и содержит 78 белок-кодирующих генов.

Сравнение предсказаний генов в базе данных GenBank и по данным Prodigal для плазмиды

Для начала я с помощью программы seqret пакета EMBOSS получила файл plasmid.fasta с последовательностью исследуемой плазмиды.

seqret embl:CP015437 plasmid.fasta

Затем я создала файл с особенностями plasmid.gff, в котором указаны предполагаемые гены.

seqret embl:CP015437 -feature plasmid.gff

Далее я запустила программу Prodigal с параметрами -i (указывает специфичный входной файл формата FASTA), -o (указывает специфичный выходной файл; по умолцанию stdout) и -f (выбор формата вывода: gbk, gff или sco; по умолчанию gbk). По совету я выбрала самый минималистичный формат sco, так как он удобен для дальнейшей работы.

Я получила файл plasmid.pro на выходе.

Для сравнения данных о предполагаемых генах я перевела файлы plasmid.gff и plasmid.pro к более удобному виду. Я воспользовалась программой gff.py, которая перевела файл plasmid.gff в gff.txt при вводе в командной строке:

python gff.py plasmid.gff gff.txt

Затем я воспользовалась программой pro.py, которая перевела файл plasmid.pro в pro.txt при вводе в командной строке:

python pro.py plasmid.pro pro.txt

Последним шагом было написание программы comparison.py, сравнивающей файлы gff.txt и pro.txt. Ниже представлен результат её работы.

Программа раборает по следующему алгоритму:

На вход она принимает два файла. Первым указывается файл gff.txt, с которым сравнивается полученный с помощью Prodigal файл pro.txt. Сначала строчки файлов переводятся в соответствующие массивы (gff_cds и pro_cds). Размер массива gff_cds (переменная all) будет равен числу генов в плазмиде. Далее он будет использоваться для подсчёта процентов. Затем я ввела и обнулила счётчики nc, n и c для посчёта числа полных совпадений, совпадений по N-концу и совпадений по C-концу соответственно. Следующим шагом я прохожу алгоритмом по каждой строчке из файла gff.txt (по каждому элементу из массива gff_cds) и разбиваю каждую строчку на отдельные "слова". Для конкретной строчки из gff_cds, если третье "слово" является символом "+", что обозначает положение гена на прямой цепи, то первое "слово" указывает на координаты начала рамки считывания (N-конец), а второе "слово" указывает на координаты её конца (С-конец). Если же запись относится к обратной цепи, то первое "слово" укажет С-конец, а второе "слово" - N-конец. Таким образом, я в переменные gff_n и gff_c записываю координаты начала и конца рамки считывания соответственно. Так, имея данные о конкретной строчке из файла gff.txt, я начинаю просматривать все строчки из файла pro.txt и аналогично записываю в переменные pro_n и pro_c координаты начала и конца рамки считывания с учётом направленности гена в цепи. Далее, если в файле pro.txt найдётся такая строчка, для которой хотя бы один из концов совпадёт с таковым в фиксированной строчке из файла gff.txt, то в случае совпадения обоих координат, прибавляется 1 к счётчику nc, иначе прибавляется 1 к счётчику n или c в зависимости от того, какие именно концы совпали. По завершении обоих циклов я имею:

  • Число полных совпадений (nc)
  • Число совпадений только по C-концу или число несовпадений по N-концу (c)
  • Число совпадений только по N-концу или число несовпадений по C-концу (n)
  • Число ошибочно предсказанных генов или число несовпадений по обоим концам (all - nc - c - n)

Для удобства я оформила результаты сравнения в таблице 1.

Таблица 1. Результаты сравнения предсказаний генов в базе данных GenBank и по данным Prodigal для плазмиды CP015437
Полностью совпавшие гены Ошибка только по N-концу Ошибка только по С-концу Ошибочно предсказанные гены
Число 76 1 0 1
Процент 97.4% 1.3% 0.0% 1.3%

Анализ несовпадений

Программа Prodigal сработала достаточно хорошо. Так что я имею всего 2 несовпадения. Рассмотрим сперва несовпадение только по N-концу. В файле gff.txt указано "28468 29013 +", а предсказано программой "28495 29013 +".

Последовательность ДНК с 28468 по 29013:

atgaaagacaggaaaaggggtattgaaatgaggaagtttatatcaacagctaatcggggc
atgggcctcgaaaattacattaacgccgcaaacacggcgtataaattaaaaaaaattgca
ttgattgacaaaaaacctgttcctatcaagtttagaactgataagaggacagggcgtatt
aaggaagcttggtttgaaaaacaaagtacggttgactactacggcattgtaaacggcgtc
cccgtcgcttttgagtgcaaaagcacacgagaagagactagattcccgctcaacaacatt
gaggagcatcaatacgagtatttgagggagtttcatttgcaaaaaggtgtcagcttcttg
cttgttgaattcgccgcacatcaagagatatacatactaaagtttgagcaactgcaacaa
tggtggttcgaatcacgccgtggtggccgcaagagcattccttataaatggttccaggaa
aattgtaagaagtgtggtcctggacgcggtctttctgttgactatctagcggcggcggga
ttatga

Предполагаемый кодируемый белок:

MKDRKRGIEMRKFISTANRGMGLENYINAANTAYKLKKIALIDKKPVPIKFRTDKRTGRI
KEAWFEKQSTVDYYGIVNGVPVAFECKSTREETRFPLNNIEEHQYEYLREFHLQKGVSFL
LVEFAAHQEIYILKFEQLQQWWFESRRGGRKSIPYKWFQENCKKCGPGRGLSVDYLAAAG
L

Последовательность ДНК с 28495 по 29013:

                           atgaggaagtttatatcaacagctaatcggggc
atgggcctcgaaaattacattaacgccgcaaacacggcgtataaattaaaaaaaattgca
ttgattgacaaaaaacctgttcctatcaagtttagaactgataagaggacagggcgtatt
aaggaagcttggtttgaaaaacaaagtacggttgactactacggcattgtaaacggcgtc
cccgtcgcttttgagtgcaaaagcacacgagaagagactagattcccgctcaacaacatt
gaggagcatcaatacgagtatttgagggagtttcatttgcaaaaaggtgtcagcttcttg
cttgttgaattcgccgcacatcaagagatatacatactaaagtttgagcaactgcaacaa
tggtggttcgaatcacgccgtggtggccgcaagagcattccttataaatggttccaggaa
aattgtaagaagtgtggtcctggacgcggtctttctgttgactatctagcggcggcggga
ttatga

Предполагаемый кодируемый белок:

         MRKFISTANRGMGLENYINAANTAYKLKKIALIDKKPVPIKFRTDKRTGRI
KEAWFEKQSTVDYYGIVNGVPVAFECKSTREETRFPLNNIEEHQYEYLREFHLQKGVSFL
LVEFAAHQEIYILKFEQLQQWWFESRRGGRKSIPYKWFQENCKKCGPGRGLSVDYLAAAG
L

Мне кажется странным, что программа предложила ген, кодирующий на 9 аминокислот меньше, чем на самом деле. Возможно, это объясняется частотой кодонов. По аннотации последовательность гена принадлежит рекомбинантному белку из U семьи. Интересно, что запустив в BLASTP по очереди оба варианта белковой последовательности, мы в обоих случаях получим схожие находки, большинство которых относится к белку Holliday junction resolvase RecU. Окружение данного гена с N-конца весьма свободно:

Что же касается полного несовпадения, то для него в файле gff.txt указано "54398 54709 -", а предсказано программой "53903 54976 +". Ниже представлено окружение гена в промежутке 53903 - 54976.

Если бы предложенный программой ген существовал, то он бы перекрывался с предыдущим геном "53214 53906 +" на 4 нуклеотида, что в принципе возможно у прокариот, но не желательно. Последовательность такого гена "53903 54976 +":

atgagaaatatttcaccgacgttgctacagtccatttcgagccagttccagctagcggct
gatgatcgggaaccgaatcaatatgttttagtgcgaaatcgtttttatggccgttcggat
gatccaggactgatggatgtagaatgggcgcggcttcctaatgtcatcagcattgaaatt
gatgaaacggtcgaaacaccggctaaatcgttcacaattgtggtggctaacaacgatggt
cgtttttcaccggactattatgcgggaaaatggccgcaacgatatctatttggcgatgaa
gtaccggagagcccctgggcgcgacaactattcccgaataccgaaattgcgattcatctc
ggatacggtgaggaagtcatcccttatttgcggggactcattgatgaggtgaccatcaac
tccgacgcgcaaacgctaacgattactgggcgttctatgtataaaaaagtattgacgacg
acaatcaagcctccatcaaacaagcagtactttctgaacgataaaactatgaacataggc
gaggcggtaaaaaaattggtcgaaggagcgggcgtgccatgcaccgctccttttattacg
gttccaggcagtaacccgccagaatcgtatgagcttggtgtggctaggggcgaaaggttt
cagacgtttgacgatgcagttcgcgatattgtagattcgttgttttttacgatggttgaa
gatgcacaggggaaggttactctcaagtcgattcctcggtattcccaagcggatgaccca
aaggtggtagtcgacgactttctccatctgaaagaactcgagtacaccatgagtgacaat
gacctatatgggaccatcatcattaaatccgggaaaaaggcggatattttccatagcggt
tttatcaccaatcaaattttgcgcgggcagagaagagaggcagaatacgaagtcccttgg
gcgaacacgtcatttaagcggaaattagcggcgcaagcctatttcacgcaaatgttacac
cgctggagaaaactatcgattgctatccccgctaacccagcgattgaattatgg

Если сделать поиск по выше написанной последовательности ДНК в TBLASTX, то программа не выдаст достоверные находки.

Запустим BLASTP с последовательностью белка из аннотации гена "53214 53906 +":

MEKVVDYHLWVIRLGIPRNRLESNLPLCIFNHRKKQRIYNIANCIVKRLKPFAPSHTKLI
RFWRVTAWNRNKRSGAWHARSFDQFFYRLAYVHSFIVQKVLLV

В качестве единственной находки мы получим гипотетический белок GFC29_3885 из исследуемой нами плазмиды.

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