Создадим файл с записью D89965 банка EMBL:
entret embl: D89965
Назовём открытой рамкой считывания нуклеотидную последовательность, начинающуюся старт-кодоном и заканчивающуюся стоп-кодоном (-find 1). Запустим программу getorf так, чтобы она выдавала все возможные трансляции (по умолчанию) всех открытых рамок считывания длиной не менее 30 нуклеотидов (-minsize 30) при бактериальном коде (-table 11):
getorf -sequence d89965.entret -table 11 -minsize 30 -find 1
Получим файл с такими трансляциями: d89965.orf.
В поле CDS файла d89965.entret указаны координаты участка [163-435]. В результатах работы программы не нашлось такой рамки считывания, наиболее похожими на нужную оказались меньшая [176 - 316] и большая [19 - 432]. Это может быть связано с несоответствием кода и организма: в млекопитающих (крыса) способ кодировки может слегка отличаться от бактерий. Однако, как будет видно из документа Swoss-Prot, этот ген способен экспрессироваться кишечной палочкой.
Чтобы сравнить полученные рамки с записью Swiss-Prot получим нужный файл hslv_ecoli.entret:
entret sw:p0a7b8
И получим аминокислотную последовательность, указанную в swiss-prot.
С помощью программы TBLASTN обнаружим, что третья рамка [176-316] не соответствует данному белку, в отличие от пятой рамки [19-432], соответствующей идеально.
С помощью программы BLASTP найдём, что записи hslv_ecoli.entret банка Swiss-Prot соответствует только последняя рамка [375-1] с кодирующей последовательностью на комплементарной цепи.
Такие результаты вполне объяснимы, так как белок, синтезирующийся по цепи mRNA будет иметь "комплементарную" к ней последовательность и практически такую же длину. Однако для меня удивительно, что обнаруженные две последовательности, кодирующие соответственно белок и мРНК, находятся не строго друг над другом: [375-1] и [19-432]. Это может быть связано с особенностями транскрипции и трансляции в исследуемом организме.

В файле trna_ecoli.fasta содержатся все возможные тРНК бактерии E. coli. Попытаемся найти гомологов данных тРНК в геноме бактерии Xanthomonas campestris. Будем искать соответствия между данными тРНК и геномной ДНК бактерий программой BLASTN:
blastall -p blastn -d xc_ -i trna_ecoli.fasta -m8 -o trna_blast.fasta
В получившемся файле для каждой тРНК есть несколько находок. С помощью Excel и grep напишем скрипт: trna_linux.scr.
Промежуточные этапы для создания скрипта:
grep ">" trna_ecoli.fasta > names.fasta
(список всех тРНК)
Создадим файл Excel: trna.xlsx.
noreturn trna.scr trna_linux.scr
(Перевели файл в формат linux)
chmod +x trna_linux.scr
(сделали скрипт исполняемым)
./trna_linux.scr
(запустили скрипт)
В результате работы скрипта был получен файл kolvo.txt, содержащий только количества находок.
Аналогично напишем скрипт2 и получим данные с ограничением порога e-value 0.001.
Запишем результаты в таблицу.
Как и ожидалось, количество находок при ограничении e-value снизилось до 1-4, а для некоторых тРНК соответствий найдено не было.

Повторим предыдущий опыт с использованием программ megablast и discontigous megablast.
Этапы работы:
megablast -d xc_ -i trna_ecoli.fasta -m8 -o trna_megablast.fasta
Напишем скрипт: trna_megablast.scr
noreturn trna_megablast.scr trna_megablast_linux.scr
chmod +x trna_megablast_linux.scr
./trna_megablast_linux.scr
Аналогично выполним задачу для программы discontigous megablast, запустив первой команду:
megablast -d xc_ -i trna_ecoli.fasta -m8 -o trna_megablastd.fasta -t 18 -W 11 -N 0,
где -t 18 - длина последовательности в шаблоне, -W 11 - длина последовательности, по которой ведётся поиск, -N 0 - поиск по кодирующей последовательности (даёт те же результаты, что и при поиске по некодирующей последовательности).
Напишем скрипт trna_megablastd.scr.
Информацию запишем в уже упоминавшийся файл Excel.

В файле Excel, являющимся результатом предыдущего задания, найдем тРНК argX, для которой BLASTN обнаружил гомологов, а MEGABLAST - нет, что неудивительно, т.к. megablast ищет мало отличающиеся друг от друга последовательности (при значительных различиях выравнивание считается плохим).
Вырежем гомологичную последовательность из бактерии Xanthomonas campestris в отдельный файл:

seqret -sask
Reads and writes (returns) sequences
Input (gapped) sequence(s): xc_genome.fasta:AE012187
     Begin at position [start]: 1868
       End at position [end]: 1891
        Reverse strand [N]:
output sequence(s) [ae012187.fasta]: argx_xc.fasta

Длина: 77
Идентичность: 24/77 (31.2%)
Сходство:     24/77 (31.2%)
Гэпы:         53/77 (68.8%)
Вес:          120.0