Я выбрала мнемонику RNC. Ей соответсвуют белки RNase III (рибонуклеазы 3). Эта рибонуклеаза расщепляет двуцепочечную РНК, участвует в обработке предшественников рРНК, вирусных транскриптов, некоторых мРНК и как минимум одной тРНК(минорная форма метианиновой тРНК).
В файле bacteria-sw.fasta белков с такой мнемоникой нашлость 493 (grep "|RNC_" /P/y24/term4/bacteria-sw.fasta | wc -l)
Для анализа этой мнемоники были выбраны белки:
Получили последовательности этих белков:
seqret -filter @id.list rnc.fasta
В id.list лежат id этих белков в виде sw:RNC_* (вместо * дописываем соответсвующие названия).
После чего выравниванием полученные последовательности в Jalview (ссылка на fasta выравнивания)
Возьмем паттерн D-x-K-[TS]-x-L-Q-E-x-[LTAI]-[QA]. В выравнивании от имеет координаты 191-200, а в последовательности белка RNC_STRCO: 184-194.
C ним находится 389 белков, из которых с нужной нам мнемоникой тоже 389. То есть не нашлось 493 - 389 = 104 белка, но ложноположительных находок нет.
Поиск проводился с помощью программы fuzzpro:
fuzzpro -sequence /P/y24/term4/bacteria-sw.fasta -pattern 'D-x-K-[TS]-x-L-Q-E-x-x-[QA]' -outfile rnc_fuzzpro_2.txt
grep "Sequence:" rnc_fuzzpro_2.txt | wc -l
grep "Sequence: RNC_" rnc_fuzzpro_2.txt | wc -l
Попробуем ослабить паттерн. Укоротим его: D-x-K-[TS]-x-L-Q-E
Проведя аналогичный анализ получаем, что было найдено 433 находки, из которых 423 соответсвуют мнемоники. То есть 10 находок найдены без нашей мнемоники, а с мнемоникой не нашлось 493-423=70.
Этот паттерн лучше подходит, чем предыдущий.
Попробуем еще ослабить паттерн. Добавим x вместо [TS]: D-x-K-x-x-L-Q-E
Проведем анализ, и получаем, что было найдено 512 находок, из которых 423 соответсвуют мнемоники. То есть 512-423=89 находок найдены без нашей мнемоники, а с мнемоникой не нашлось по по-прежнему 70.
Получается, что это ослабление паттерна, было бессмысленно, и предыдущий (D-x-K-[TS]-x-L-Q-E) подходит больше всего.
Для тех же белков с мнемоникой RNC, которые использовались при построении паттерна, был выполнен поиск мотивов программой MEME. Последовательности были заданы как аминокислотные; использовалась модель oops, то есть по одному вхождению мотива на последовательность. Минимальная длина мотива составляла 8 аминокислот, максимальная — 15, число мотивов — до трёх.
Команда запуска MEME:
meme rnc.fasta -protein -mod oops -minw 8 -maxw 15 -nmotifs 3 -oc meme
MEME нашла три мотива. Все находки с низким p-value. Так как использовалась модель oops, каждый мотив искался в каждой входной последовательности. Найденные мотивы можно рассматривать как консервативные участки белков RNase III. Результаты MEME
Затем найденные мотивы были использованы для поиска по файлу /P/y24/term4/bacteria-sw.fasta программой MAST.
Команда запуска MAST:
mast meme/meme.html /P/y24/term4/bacteria-sw.fasta -oc mast
Результаты MAST
Было найдено 522 белка, из них 493 имели выбранную мнемонику RNC_*. Остальные 29 находок были отнесены к ложноположительным. Однако, они имеют сильно большее E-value, чем белки с нашей мнемоникой, то есть можно было бы поставить фильтрацию по значению E-value, и ложноположительных находок стало бы сильно меньше. В то же время, все 493 белка с нашей мнемоникой были найдены! Результаты гораздо лучше, чем при ручном подборе мотива.
Последовательность Шайна-Дальгарно — это участок бактериальной мРНК, расположенный перед старт-кодоном CDS. Он участвует в посадке рибосомы на мРНК за счёт комплементарного взаимодействия с 3'-концом 16S рРНК. Как консенсус последовательности был выбран мотив AGGAGG.
Поиск мотива проводился программой fuzznuc в геноме Klebsiella quasipneumoniae. Длина генома: 5323927 нуклеотидов.
Команда для поиска на прямой цепи (807 находок):
fuzznuc -sequence ~/term1/genome/GCF_016415705.1_ASM1641570v1_genomic.fna -pattern 'AGGAGG' -outfile sd_forward
Команда для поиска с учётом комплементарной цепи(1644 на обеих цепях):
fuzznuc -sequence ~/term1/genome/GCF_016415705.1_ASM1641570v1_genomic.fna -pattern 'AGGAGG' -complement -outfile sd_both
На комплементраной цепи соответсвенно 1644-807=837 находок.
Для оценки ожидаемого числа случайных находок были рассчитаны частоты нуклеотидов в геноме.
Вероятность случайного появления мотива AGGAGG на прямой цепи была рассчитана как:
P(AGGAGG) = P(A)² * P(G)⁴
Ожидаемое число находок на прямой цепи вычислялось как (L - 6 + 1) * P(AGGAGG), где L — длина генома.
Для комплементарной цепи учитывалось, что мотиву AGGAGG соответствует обратнокомплементарная последовательность, поэтому ожидаемое число находок на комплементарной цепи рассчитывалось с использованием частот C и T.
| Цепь | Наблюдаемое число | Ожидаемое число | Z-score |
|---|---|---|---|
| Прямая цепь | 807 | 1663.25 | -20.995 |
| Комплементарная цепь | 837 | 1688.92 | -20.730 |
Число случайных вхождений короткого мотива в геноме можно приближённо описывать распределением Пуассона (геном большой, вероятность встретить мотив маленькая). Z-статистика показывает, на сколько стандартных отклонений наблюдаемое число находок отличается от пуассоновского ожидания.
Поэтому для оценки достоверности отличия использовался Z-score, где O — наблюдаемое число находок, а E — ожидаемое:
Z = (O - E) / sqrt(E)
В обоих случаях модуль Z-score намного больше 1.96, а значит на уровне значимости 5% мы отвергаем гипотезу. Следовательно, число найденных вхождений мотива AGGAGG достоверно меньше ожидаемого по случайной модели, учитывающей только частоты отдельных нуклеотидов.
Далее было проверено, какая часть найденных мотивов располагается в правильной позиции относительно старт-кодона CDS. Поскольку последовательность Шайна-Дальгарно обычно находится примерно за 8 нуклеотидов перед старт-кодоном, правильной считалась находка, расположенная на той же цепи, что и CDS, на расстоянии 5–12 нуклеотидов перед старт-кодоном.
Для CDS на прямой цепи использовался критерий:
5 ≤ CDS_start - SD_end ≤ 12
Для CDS на комплементарной цепи использовался критерий:
5 ≤ SD_start - CDS_end ≤ 12
Проверка была выполнена для всех найденных вхождений мотива с использованием геномной таблицы локальных особенностей.
Всего было проанализировано 1644 находки. Из них 168 находок располагались в правильной позиции относительно старт-кодона какого-либо CDS.
что составляет 10.22%.
Скрипт с анализом расположения относительно старт-кодонов CDS
Таким образом, только небольшая часть всех геномных находок мотива AGGAGG располагается в положении, характерном для последовательности Шайна-Дальгарно. Остальные находки, вероятно, являются случайными совпадениями, поэтому для распознавания потенциальной последовательности Шайна-Дальгарно важна не только сама последовательность мотива, но и её положение относительно старт-кодона CDS.