Задания по предсказанию генов прокариот


Мне досталась плазмида с идентификатором CP010831. Она принадлежит бактерии Shigella sonnei штамма FORC_011 (рис. 1). Длина плазмиды - 83905 пар нуклеотидов, 85 кодирующих регионов, все из которых кодируют белки.


Рис. 1
Shigella sonnei

1. С помощью команд пакета EMBOSS
seqret embl:CP010831 -feature cp010831.gff
seqret embl:CP010831 cp010831.fasta
из базы данных embl были получены fasta-файл последовательности плазмиды (cp010831.fasta) и gff-файл (cp010831.gff), содержащий информацию об особенностях, благодаря параметру -feature.

2. Информация о координатах генов (начало, конец, ориентация) была извлечена из файла cp010831.gff с помощью команды
grep 'CDS' cp010831.gff | awk '{print $4,$5,$7}' >> info_genbank.txt
был получен файл info_genbank.txt. Команда grep находит в файле строки с паттерном CDS, затем awk находит 4-й, 5-й и 7-й столбцы gff-файла, где содержатся соответствующие координаты, и записывает в файл info_genbank.txt.

3. Для предсказания генов с помощью Prodigal, программа была запущена из командной строки таким образом:
prodigal.windows.exe -i cp010831.fasta -f sco -o prodigal -s prodical_genes
Параметр -i - входной файл, -o - выходной файл, -s - записывает потенциальные гены в соответствующий файл, -f - выбирает формат выходных файлов.
Получены файлы prodigal и prodigal_genes.
Prodical выдала, что, так как длина моей плазмиды всего 83905 пар нуклеотидов (а для оптимальной работы программы неоходимо хотя бы 100000), для получения лучших результатов стоит запустить ее с мета-опцией. А кто я такая, чтобы с ней спорить? Дальше будет видно, что так действительно результаты предсказания гораздо лучше.
В результате работы программы
prodigal.windows.exe -i cp010831.fasta -f sco -o prodigal -s prodical_genes -p meta
получены файлы prodigal_meta и prodical_genes_meta.

4. Для удобства дальнейших манипуляций, к файлам prodigal_meta и prodigal была пременена команда
grep '>' prodigal_meta |sed "s/_/ /g" | awk '{print $2,$3,$4}' >> info_prodigal.txt и grep '>' prodigal |sed "s/_/ /g" | awk '{print $2,$3,$4}' >> info1.txt
соответственно. Команда grep ищет строки, где содержится паттерн ">". Затем в найденных строках команда sed заменяет "_" на пробел: s (substitute) заменяет первое встреченное соответствие первому шаблону, на второй шаблон; g (global) гарантирует выполнение команды всеми найденными соответствиями внутри каждой из строк. Команда awk, как описано выше, находит столбцы 2, 3 и 4, которые затем записываются в файл info1.txt.

5. Далее с помощью скрипта были получены таблица 1 и 2, диаграмма 1 и 2. Причем, таблица и диаграмма 1 представляют оценку предсказания Prodigal, запущенной с мета-опцией, как программа и рекомендовала. Таблица 2 и диаграмма 2 - для Prodigal с параметрами по умолчанию. Очевидно, что с мета-опцией результат редсказания в разы лучше.
Скрипт вызывается так: python kek.py <файл со списком реальных генов> <файл со списком предсказанных генов> <название выходного файла без расширения>. Пояснения, касающиеся его работы, представлены файле.

Таблица №1
результат предсказаний Prodigal (meta)
Количество генов Процент
Абсолютно верно предсказанных 64 75,29
Только N-конец не совпадает 6 7,06
Только C-конец не совпадает 0 0
Оба конца не совпадают 15 17,65

Диаг. 1 Результат предсказания Prodigal (meta)

Таблица №1
Результат редсказаний Prodigal
Количество генов Процент
Абсолютно верно предсказанных 0 0
Только N-конец не совадает 0 0
Только C-конец не совадает 37 43,53
Оба конца не совпадают 48 56,47

Диаг. 2 Результат предсказаний Prodigal (auto)

6. Чтобы выяснить причину несовпадения в конкретных случаях, я посмотрела соответствующее место в геномном браузере.
Prodigal сочла сочетания нуклеотидов, характерные для старт-кодонов, началами гена, в то время как эти кодоны либо кодируют метионин в середине белковой цепи (рис. 1 и 2), либо, наоборот, нашла "лишний" старт-кодон (рис. 3 и 4).
На рисунках, полученных из геномного браузера видно, что белок и идентификатором ALZ58987.1 (26610..26894) заметно сократился из-за того, что программа определила его начало с другово места (с 26718 нуклеотида рямой цепи; на рис. 2 и 3 обозначен синим).
Ген белка с идентификатором ALZ59018.1 (51780..52682), напротив, согласно Prodigal, начат раньше (с 52706 нуклеотида; на рис. 4 и 5 обознаен зеленым).
На рисунках

Рис. 2 Ген белка с идентификатором ALZ58987.1

Рис. 3 Увеличенная картинка

Рис. 4 Ген белка с идентификатором ALZ59018.1

Рис. 5 Увеличенная картинка