Блок 4: домены белка

Задание 1

В базе pfam.xfam.org было найдено семейство HBB (идентификатор в базе: PF06777). Белки семейства входят в качестве доменов в ДНК-репарирующие хеликазы эукариот. Хеликаза формирует тороидальную структуру для разделения цепей ДНК и поиска повреждений, а HBB-домен помогает зажимать ДНК. Если хеликаза обнаруживает повреждение, запускается механизм репарации. Для семейства HBB существует 17 доменных архитектур в 492 последовательностях у 450 особей. Для оценки были выбраны архитктуры:
  • DEAD_2, HBB, Helicase_C_2 (402 последовательности)(pdb со структурой)
  • DEAD_2 x 2, HBB, Helicase_C_2: (15 последовательностей).

Скачивать выравнивание из базы JalView отказался (при полной работоспособности других сетевых сервисов навроде muscle и даже после попытки скачать предлагаемый программой пример AC PF03760), поэтому файл с выравниванием всех последовательностей домена был получен непосредственно с pfam.org и покрашен в JalView (проект JalView).

Далее с помощью скрипта python swisspfam-to-xls.py -p PF06777 -o PF067777_all.xls -z была получена доменная структура и идентификаторы всех последовательностей домена. Идентификаторы (с удалением дубликатов на новом листе в excel) были вытащены в файл, с помощью которого из uniprot.org (Retrieve/id mapping, загрузить файл-список идентификаторов) были получены последовательности домена. Наверное, их можно было получить гораздо проще, просто убрав гэпы из выравнивания последовательностей в JalView и сохранив невыровненными, но в подсказках к заданию почему-то было сказано действовать другим путём. Возможно, смысл в том, что выравнивание, полученное с pfam, статично и не обновлялось: при запросе по всем 492-м AC, входящим в домен PF06777 по версии pfam, uniprot выдал только 480, остальные оказались удалены (с пометкой deleted в столбце protein names). Так что реально в домене теперь 480 последовательностей.

Файл с последовательностями домена из uniprot (осторожно, 2,5Мб текста)был подан на вход скрипту python uniprot-to-taxonomy.py -i PF06777_from_uniprot.txt -o PF06777_taxo.xls, на выходе получен список таксонов, входящих в домен. Список таксонов и доменная структура скопированы в основной файл. Далее на отдельный лист вынесен столбец с AC последовательностей и столбец с доменной структурой, столбцу с AC поставлены с помощью ВПР([номер ячейки с искомым];[таблица, в которой ищем]; [номер столбца таблицы (в которой ищем), откуда брать результат для найденного(в единственной строке таблицы) искомого]; ЛОЖЬ). Теперь можно отбирать представителей архитектуры визуально. Надо выбрать 2 комбинации таксон-подтаксон так, чтоб в каждой было не менее 5 представителей с каждой из выбранных доменных архитектур. Итого получится 4 группы:

  1. Eukaryota/Metazoa с (DEAD_2, HBB, Helicase_C_2) 12 последовательностей
  2. Eukaryota/Metazoa c (DEAD_2,DEAD_2, HBB, Helicase_C_2) 3 последовательностей
  3. Eukaryota/Fungi с (DEAD_2, HBB, Helicase_C_2) 11 последовательностей
  4. Eukaryota/Fungi c (DEAD_2,DEAD_2, HBB, Helicase_C_2)4 последовательностей

Также к выборке была добавлена последовательность XPD_THEAC, для которой есть 3D-структура. Идентификаторы были собраны в текстовый файл, переданный в retieve/ID mapping на uniprot с получением fasta-файла с 31ой отобранной последовательностью.

Из выравнивания была удалена последовательность J9DJ02 (Eukaryota/Fungi c (DEAD_2,DEAD_2, HBB, Helicase_C_2)), плохо выровненная с остальными. Последовательности SHEEP и THEKT, несмотря на длинные невыровненные участки, в целом совпадают по архитектуре, поэтому были оставлены.

Выравнивание показало, что второй домен DEAD есть в большинстве последовательностей, хотя в pfam его наличие почему-то не учитывается:

В файле выравнивания начиная с 80ой позиции можно видеть консервативные колонки второго DEAD-домена и у тех последовательностей, которые в pfam значатся с архитектурой DEAD_2, HBB, Helicase_C_2, то есть с одним DEAD-доменом. В то же время структуру DEAD_2, DEAD_2, HBB, Helicase_C_2 имело большинство последовательностей, значившиеся "deleted" на uniprot. Возможно, первый DEAD-домен не имеет функциональной нагрузки и потому не учитывается в семействе.

Задание 2


Для выполнения второго задания - построения дерева по выравниванию - само выравнивание было немного подправлено: сокращены описания последовательностей и помечены таксоны, а также убрана невыровненная последовательность. Не входящая в отобранные таксоны последовательность XPD_THEAC оставлена для укоренения дерева. Сокращения по таксонам:
  • Eu_Me = Eukaryota/Metazoa
  • Eu_Fu = Eukaryota/Fungi
  • Ar_Eu = Archaea/Euryarchaeota

Fasta-файл c последовательностями выровнен в MEGA, удалены невыровнявшиеся N- и C-концы. По Fasta-файлу с выравниванием там же и построено дерево по алгоритму максимального сходства. Алгоритмы Neighbor Joining и Minimal Evolution дают ту же топологию. Видно, что разница в архитектуре, найденная pfam, никак не отражается на положении в дереве. Однако домены внутри таксонов оказываются более схожими, чем между ними: млекопитающие оказались в одном районе дерева.

Дерево доступно в Newick-формате

Задание 3


Для построения профиля была выбрана клада дерева, отмеченная на рисунке:

Из выравнивания, использованного при построении дерева в задании 2, получен (копированием нужных участков) файл с 6-ю последовательностями выбранной клады, откуда в JalView были удалены пустые (где были только гэпы) столбцы. Fasta-файл c выравниванием был использован для построения профиля в HMM:
hmm2build hmm_built for_profile_cut_and_aligned.fasta, ( построенный профиль)
Затем профиль калибруется:
hmm2calibrate --histfile hmm_histo hmm_built, ( гистограмма калибровки)
После калибровки профиль используется для поиска последовательностей в файле со всеми последовательностями, содержащими семейство PF06777:
hmm2search hmm_built PF06777_from_uniprot.fasta > hmm_search, (результаты поиска)
Из результатов поиска извлекаются данные по находкам, вставляются в xlsx-файл и упорядочиваются по возрастанию e-value. Для построения ROC-кривой требуется рассчитать чувствительность (отношение числа последовательностей, сочтённых профилями, при том, что они являются профилями, к общему числу профилей, True Positive) и специфичность (отношение числа последовательностей, не сочтённых за профиль, при том, что они профилем и не являются, к общему числу не-профилей, False Negative). Чувствительность для ячейки K считаем как число профилей в ячейках K+1...N, делённое на общее число профилей. Специфичность для ячейки K считаем как число не-профилей в ячейках 1...K, делённое на общее число не-профилей. Потом берём значения (1-специфичность), то есть неспецифичность, в качестве аргумента (ось абсцисс) и строим по ним и значениям чувствительности график:

Для граничного e-value по графику была найдена точка выполаживания кривой, однако выяснилось, что это пороговое значение нулевое: точка лежит на последней профильной последовательности F6YA33_MOUSE. То есть из-за малого объёма выборки вероятность ошибки первого рода оказалась очень высокой. Таблица результатов:
Входит в семействоНе входит в семействоСумма
Выше порога62430
Ниже порога0450450
Сумма6474480