Цель: взять указанный фрагмент генома длиной 7 тысяч нуклеотидов и найти в нем предполагаемых гомологов белков сенной палочки.
Для начала получим фрагмент: seqret -sask, по запросу программы укажем запись embl:AEWH01000012, начало 126 001, конец 133 000, оставим кусочек генома в файле myfragment.fasta.
Далее с помощью seqret sw:*_BACSU -out bacsu.fasta получим полный (почти) протеом сенной палочки (файл весит полтора мегабайта и был удален, как только перестал быть нужен). С помощью makeblastdb -in bacsu.fasta -out bacsu получим индексные файлы для поиска (здесь также не приводятся).
Трансляции рамок считывания (между старт- и стоп-кодонами) нам выдаст getorf -sequence myfragment.fasta -table 11 -minsize 240 -find 1 -out orfs.fasta. Отсюда информацию в нужном виде легко получить такой командой: grep ">" orfs.fasta | sed -e "s/>//" -e "s/\[//" -e "s/\]//" -e "s/- //" -e "s/ AE.*//" -e "s/AE.*_//" -e "s/(REVERSE SENSE)/reverse/" > orfs_processed.txt. ,
Запустив blast: blastp -db bacsu -query orfs.fasta -evalue .001 -outfmt 6 > blast_6.txt, прогоним выходной файл через такое: for i in {1..14}; do grep -Pc "AEWH01000012_$i\t" blast_6.txt; done и получим количество находок BLAST для каждого белка; всё вместе скомпонуем в orfs.xls.
Сделаем таблицу про рамки, для которых нашлось хоть одно совпадение:
| # | Начало | Конец | Направление | Найдено | ID ближайшего | E-value |
| 2 | 1565 | 2065 | прямое | 1 | YKHA_BACSU | 1e-29 |
| 4 | 6760 | 6999 | прямое | 2 | CWLJ_BACSU | 7e-18 |
| 7 | 6089 | 5475 | обратное | 1 | UBIE_BACSU | 2e-05 |
| 9 | 5240 | 4512 | обратное | 1 | YGAJ_BACSU | 4e-82 |
| 11 | 4426 | 3185 | обратное | 6 | YFIS_BACSU | 7e-14 |
| 12 | 2868 | 2341 | обратное | 1 | YDFB_BACSU | 1e-18 |
| 14 | 1109 | 219 | обратное | 4 | CZCD_BACSU | 3e-51 |
Схематическое представление положения найденных открытых рамок считывания, одновременно представляющее собой карту генов:
Теперь нам нужно найти, где находятся гены, соответствующие найденным самым близким белкам, в геноме сенной палочки. Для этого возьмем запись с полным ее геномом. И, в общем, на этом придется и закончить, потому что гены расшвыряло по геному так, что какой-то гомологии в расположении не получится найти совсем никак: никакие два из них не находятся даже ближе нескольких десятков тысяч нуклеотидов друг от друга, расположение на комплементарных цепях тоже беспорядочно.