Standalone BLAST

Поиск в геноме участков, кодирующих белки, похожие на заданный

Команды:
  • makeblastdb -in bl_genome.fasta -dbtype nucl, где bl_genome.fasta - файл с геномом Baccilus licheniformis - делаем базу данных для бласта.
  • blastn -task blastn -query seq.txt -db bl_genome.fasta -evalue 0.001 -out blast.out, где seq.txt - файл с последовательностью OxdC.

    • Поиск гомологов белка OXDC_BACSU в геноме Baccilus licheniformis
    Число находок с E-value<0,001 2*
    E-value лучшей находки 0.0
    Название последовательности с лучшей находкой AE017333 Bacillus licheniformis DSM 13, complete genome.
    Координаты лучшей находки (от-до) 3420589 - 3421728
    Доля последовательности OxdC, вошедшая в выравнивание с лучшей находкой 380/385=0,99%
    *Как известно из прошлых практикумов, у OxdC есть близкий паралог - OxdD.

    Поиск гомологов некодирующих последовательностей программой BLASTN

    Нужно, используя файл с последовательностями всех тРНК из генома B. subtilis, узнать, сколько гомологов находится в геноме Baccilus licheniformis для этих последовательностей.
    Команда: blastn -task blastn -query trna_bacsu.fasta -db bl_genome.fasta -out trnabl.out -evalue 0.01 -outfmt 7, где trna_bacsu.fasta - файл с последовательностями тРНК, -evalue [цифра] - порог e-value, a -outfmt [цифра] - вид вывода результатов, в данном случае таблица
    Далее с помощью команд bash было необходимо подсчитать количество результатов для каждой тРНК. Для этого используется следующий скрипт. С помощью FAR переведём его в формат для Unix и активируем командой chmod +x. Запуск скрипта через Putty производим так: ./scr.sh, если находимся в той директории, где находится и файл скрипта. Таблицу с результатами см. в следующем задании.

    Поиск гомологов при изменённых параметрах программы BLASTN

    При попытках изменить reward и penalty - параметры подсчёта веса выравнивания ("вознаграждение" за совпадение и штраф за миссметч), программа не работает и пишет следущее: 
    10 and 6 are supported existence and extension values
8 and 6 are supported existence and extension values
25 and 10 are supported existence and extension values
Any values more stringent than 25 and 10 are supported
    Выберем в качестве штрафов за открытие и продление гэпа 26 и 11.
    Команда: blastn -task blastn -query trna_bacsu.fasta -db bl_genome.fasta -out trnabl.out -evalue 0.01 -outfmt 7 -reward 5 -penalty -4 -gapopen 26 -gapextend 11
    Теперь изменим и длину слова:
    blastn -task blastn -query trna_bacsu.fasta -db bl_genome.fasta -out trnabl.out -evalue 0.01 -outfmt 7 -reward 5 -penalty -4 -gapopen 26 -gapextend 11 -word_size 4
    Результаты (число найденных тРНК каждого типа) представлены в этой ТАБЛИЦЕ.

    Анализ полученных результатов

    При изменении штрафов результаты изменяются незначительно, а вот в сочетании с маленькой длиной слова число результатов резко подскакивает, возможно, в связи с нахождением гомологов среди тРНК других типов (не к той же аминокислоте) или просто случайных последовательностей, всё-таки длина любой тРНК не слишком большая. Рассмотрим несоответствие: сериновая тРНК, которая находится только при уменьшении длины слова: Выравнивание needle (см. ниже) довольно плохое - имеет маленькую длину: 
    Aligned_sequences: 2
# 1: AE017333
# 2: BSn5_t20970
# Matrix: EDNAFULL
# Gap_penalty: 10.0
# Extend_penalty: 0.5
# Length: 93
# Identity:      27/93 (29.0%)
# Similarity:    27/93 (29.0%)
# Gaps:          62/93 (66.7%)
# Score: 119.0

AE017333           0 --------------------------------------------------      0
                                                                       
BSn5_t20970        1 ggaggaatacccaagtccggctgaagggatcggtcttgaaaaccgacagg     50

AE017333           1 ------------cgcgagggttcaaatccctccttctccgcca     31
                                 ||||.||||||.||||||||.|.||||||||
BSn5_t20970       51 ggtgtcaaagcccgcgggggttcgaatccctcttcctccgcca     93
    Как и ожидалось, никакой тРНК в геноме в этом месте не нашлось, т.е. это просто случайный кусок генома.