Построение профиля подсемейства и проверка его работы.
Рис. 1. Укорененноное дерево. Выбранное подсемейство выделено красным.
Для создания профиля из выравнивания, полученнного в предыдущем практикуме в отдельный файл j.fasta были извлечены представители, составляющие подсемейство. Для построения и калибровки профиля применялся пакет HMMER, установленный на kodomo. Использовались следущие команды: Построение профиля:
hmm2build profile15 j.fastaКалибровка:
hmm2calibrate profile15Полученный файл: profile15 Далее из Uniprot я по частям скачал все белки, содержащие мой домен. Скачанные файлы были объединены в один файл db.fasta. После этого был произведен поиск по всем белкам UniProt, содержащим домен.
Использованная команда:
hmm2search profile11 db.fasta >> outfileПолученный файл: outf
Далее я вставил находки в файл и отметил среди находок представителей подсемейства, построил гистограмму весов находок и ROC-кривую. Для построения гистограммы использовался скрипт, который обработал файл и выдал гистогамму в виде строки, котрая графически была построена в Excel(см. итоговый файл)
Рис. 1. Гистограмма весов находок. Вес: 1- >500 2 -(501;400) 3 - (401;300) 4-(301;200) 5-(201;100) 6-(101;0) 7-(1;-100) 8-<-100
Для построения ROC-кривой, я на основание столбца Profile(1-если последовательность входит в подсемейство, 0 -если нет) рассчитал значения чувтвительности(SE) и величины 1- специфичность(1-SP). Далее я построил на основе этих данных ROC-кривую(представлена на рисунке 2).
Рис. 2. ROC-кривая.
В результате я получил итоговый файл.
На основании построенной ROC-кривой было выбрано пороговое значение E-value. Критерием выбора являлось максимальное значение разности [Чувствительность - (1-Специфичность)]. Полученный порог E-value - 0,019224. При таком пороге чувствительноть профиля составляет 1,0000000, а специфичность 0,98077. Это достаточно хорошие значения, поэтому построенный профиль вполне можно использовать для выделения подсемейства.