Практикум 9. Выравнивание последовательностей

Задание 1. Подсчёт инделей

В этом задании нужно было написать программу на Python, которая принимает имя файла с выравниванием (результат работы needle или water) и считает количество инделей в каждой из двух последовательностей, а также их общее число.

Важно, что под инделем понимается не каждый отдельный символ -, а один непрерывный участок пропусков. Например, --- считается как один индель. Именно это определяет логику решения.

Решение состоит из нескольких шагов. Сначала программа получает имя файла из аргументов командной строки с помощью argv. Это сделано в соответствии с условием задания, где сказано, что программа не должна использовать input().

Далее файл открывается и считывается построчно. Поскольку формат вывода EMBOSS содержит служебные строки и разбитые на блоки фрагменты последовательностей, данные нельзя обрабатывать напрямую — сначала их нужно разобрать.

На первом проходе по файлу извлекаются имена последовательностей. Они находятся в строках, которые начинаются с # 1: и # 2:. С помощью метода split() берётся часть строки после двоеточия — это и есть нужное имя.

После этого выполняется второй проход по файлу. На этом этапе собираются сами последовательности. Строки, начинающиеся с #, пропускаются, так как они являются служебными. Остальные строки разбиваются на части, и если формат строки подходит, из неё извлекается фрагмент последовательности.

Поскольку выравнивание записано блоками, последовательность собирается постепенно: все фрагменты с одинаковым именем объединяются в одну строку.

Подсчёт инделей

Для подсчёта инделей используется отдельная функция count_indels. Она проходит по последовательности посимвольно и отслеживает, находимся ли мы внутри участка пропусков.

Если встречается символ - и до этого не было пропуска, значит начинается новый индель, и счётчик увеличивается. Если дефисы идут подряд, они считаются как один индель, потому что программа запоминает, что уже находится внутри разрыва.

Когда встречается обычный символ, это означает конец текущего инделя, и программа «сбрасывает» состояние. Благодаря этому корректно считаются именно непрерывные участки, а не отдельные дефисы.

Код программы

from sys import argv

def count_indels(seq):
    count = 0
    in_gap = False
    for ch in seq:
        if ch == '-':
            if not in_gap:
                count += 1
                in_gap = True
        else:
            in_gap = False
    return count

def main():
    filename = argv[1]
    f = open(filename)
    lines = f.readlines()
    f.close()

    name1 = ""
    name2 = ""

    for line in lines:
        if line.startswith("# 1:"):
            name1 = line.split(": ", 1)[1].strip()
        if line.startswith("# 2:"):
            name2 = line.split(": ", 1)[1].strip()

    seq1 = ""
    seq2 = ""

    for line in lines:
        if line.startswith("#"):
            continue
        parts = line.split()
        if len(parts) < 4:
            continue
        if parts[0] == name1:
            seq1 += parts[2]
        elif parts[0] == name2:
            seq2 += parts[2]

    indels1 = count_indels(seq1)
    indels2 = count_indels(seq2)

    print(name1, indels1)
    print(name2, indels2)
    print("Total", indels1 + indels2)

main()

Финальный этап

После того как обе последовательности собраны, для каждой из них вызывается функция подсчёта инделей. Затем результаты выводятся в нужном формате: сначала имя и число инделей для первой последовательности, потом для второй, и в конце — их сумма.

Задание 2. Глобальное парное выравнивание гомологичных белков

В рамках данного задания было проведено глобальное парное выравнивание трёх пар гомологичных белков из организмов Escherichia coli и Bacillus subtilis с использованием программы needle (алгоритм Нидлмана–Вунша). Выравнивание выполнялось с параметрами по умолчанию.

Для анализа были выбраны следующие белки:

• 6-phosphogluconate dehydrogenase (6PGD)
• 6-phosphogluconolactonase (6PGL)
• Acetyl-CoA carboxylase CT alpha subunit (ACCA)

Таблица 1. Характеристики глобального выравнивания

Глобальное выравнивание показало существенные различия в степени сходства между выбранными белками. Для 6PGD наблюдается высокая степень консервативности, тогда как для 6PGL сходство значительно ниже. Белок ACCA занимает промежуточное положение.

Задание 3. Локальное парное выравнивание гомологичных белков

На данном этапе было выполнено локальное выравнивание тех же пар белков с использованием программы water (алгоритм Смита–Ватермана). В отличие от глобального подхода, локальное выравнивание выделяет наиболее консервативные участки последовательностей.

Таблица 2. Характеристики локального выравнивания

Дополнительно были рассчитаны значения покрытия (coverage), отражающие долю последовательности, вовлечённую в выравнивание. Например, для белка 6PGL:

• Длина выравнивания: 265
• Длина последовательности E. coli: 331 → (265−1)/331 ≈ 79.8%
• Длина последовательности B. subtilis: 349 → (265−15)/349 ≈ 71.6%

Локальное выравнивание позволяет более точно выявить консервативные регионы, исключая участки с низким сходством.

Задание 4. Комментарии к выравниваниям

Было проведено сравнительное исследование результатов глобального (needle) и локального (water) выравнивания для трёх пар белков. Анализ включал оценку идентичности, сходства, числа разрывов, indels и покрытия.

1. 6PGD (6PGD_ECOLI vs 6PGD_BACSU)

Для данной пары наблюдается высокая степень сходства по всей длине последовательностей. Идентичность составляет около 70%, сходство превышает 83%. Покрытие практически полное (≈99.8%), а количество разрывов минимально.

Это свидетельствует о том, что белки являются гомологичными по всей длине и выполняют одну и ту же функцию. Локальное выравнивание практически не отличается от глобального, поскольку вся последовательность уже хорошо согласуется.

Вывод: глобальное выравнивание полностью достаточно для анализа высококонсервативных белков.

2. 6PGL (6PGL_ECOLI vs 6PGL_BACSU)

В данном случае наблюдается значительно более низкая степень сходства: идентичность составляет 25–30%, сходство — около 42–49%. Покрытие неполное, а количество gaps и indels существенно выше.

Это указывает на наличие лишь частичной гомологии. Белки содержат отдельные консервативные участки, но не совпадают по всей длине.

Локальное выравнивание оказывается более информативным, поскольку позволяет выделить функционально значимые регионы. В глобальном выравнивании такие участки «размываются» из-за необходимости выравнивать всю последовательность.

Также наблюдается ситуация, когда некоторые совпадающие позиции в локальном выравнивании отсутствуют в глобальном. Это связано с различиями в алгоритмах: глобальный метод стремится оптимизировать выравнивание по всей длине, а локальный — только в пределах наиболее схожих областей.

Вывод: для белков с частичной гомологией локальное выравнивание значительно более информативно.

3. ACCA (ACCA_ECOLI vs ACCA_BACSU)

Для этой пары характерен средний уровень сходства: идентичность составляет около 51–54%, сходство — 66–70%. Покрытие остаётся высоким (>95%), а количество разрывов невелико.

Это позволяет говорить о гомологии по всей длине, но с меньшей степенью консервативности, чем в случае 6PGD.

Локальное выравнивание немного улучшает показатели, исключая слабосходные участки, однако принципиально новой информации не даёт.

Вывод: при умеренной гомологии оба метода дают сопоставимые результаты.

Общий вывод

Проведённый анализ показал, что выбор метода выравнивания зависит от степени сходства последовательностей:

• Высокая гомология (6PGD) → достаточно глобального выравнивания
• Частичная гомология (6PGL) → предпочтительно локальное выравнивание
• Средняя гомология (ACCA) → оба метода дают схожие результаты

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

Задание 5. Результат применения программ выравнивания к неродственным белкам

1. Постановка задачи и выбор белков

Целью данного задания было проверить, как программы глобального и локального выравнивания работают на паре белков, которые не являются близкими гомологами и имеют разные функции. Для анализа была выбрана случайная пара белков из Escherichia coli с разными мнемониками функций:

• RPOA_ECOLI — DNA-directed RNA polymerase subunit alpha;
• THIO_ECOLI — Thioredoxin 1.

Выбранные белки относятся к разным функциональным классам. RPOA_ECOLI является субъединицей ДНК-зависимой РНК-полимеразы и участвует в транскрипции, тогда как THIO_ECOLI является тиоредоксином и участвует в окислительно-восстановительных процессах. Поэтому данную пару можно использовать как пример неродственных или, по крайней мере, неочевидно родственных белков.

Длины исходных последовательностей:

• RPOA_ECOLI — 329 аминокислотных остатков;
• THIO_ECOLI — 109 аминокислотных остатков.

Для этой пары были выполнены глобальное выравнивание программой needle и локальное выравнивание программой water с параметрами по умолчанию.

2. Глобальное выравнивание программой needle

Глобальное выравнивание было выполнено программой needle, реализующей алгоритм Нидлмана–Вунша. Особенность глобального выравнивания состоит в том, что программа пытается сопоставить последовательности по всей длине, даже если они имеют разную длину и слабое сходство.

Таблица 3. Характеристики глобального выравнивания неродственных белков

По данным выдачи needle, длина глобального выравнивания составила 349 позиций. Из них идентичными оказались только 26 позиций, то есть 7.4%, а сходными — 42 позиции, то есть 12.0%. Количество позиций с гэпами составило 260, что соответствует 74.5% длины выравнивания.

Подсчёт инделей с помощью программы, написанной в задании 1, дал следующий результат:

RPOA_ECOLI 2
THIO_ECOLI 5
Total 7

Таким образом, в глобальном выравнивании было найдено 7 инделей. При этом важно различать количество позиций с гэпами и количество инделей: 260 позиций с гэпами не означают 260 отдельных событий вставки или делеции, потому что несколько подряд идущих дефисов считаются одним инделем.

3. Интерпретация глобального выравнивания

Глобальное выравнивание этой пары белков получилось слабым и плохо интерпретируемым с биологической точки зрения. Процент идентичности равен только 7.4%, что значительно ниже уровня, при котором можно уверенно говорить о сходстве белков по последовательности.

Особенно показательно большое количество гэпов: они занимают почти три четверти длины выравнивания. Это связано с тем, что программа needle обязана выровнять обе последовательности целиком. Поскольку RPOA_ECOLI примерно в три раза длиннее THIO_ECOLI, алгоритм вынужден вводить протяжённые разрывы, чтобы формально сопоставить последовательности по всей длине.

Следовательно, результат глобального выравнивания в данном случае не следует рассматривать как подтверждение гомологии. Он скорее показывает ограничение глобального алгоритма: если применить его к белкам разной длины и разной функции, программа всё равно построит выравнивание, но такое выравнивание может быть формальным и биологически неубедительным.

4. Локальное выравнивание программой water

Локальное выравнивание было выполнено программой water, реализующей алгоритм Смита–Ватермана. В отличие от глобального выравнивания, этот алгоритм не пытается сопоставить последовательности полностью, а ищет лучший локальный участок сходства.

Таблица 4. Характеристики локального выравнивания неродственных белков

По данным выдачи water, длина локального выравнивания составила 124 позиции. В этом участке было найдено 26 идентичных позиций (21.0%) и 42 сходные позиции (33.9%). Количество позиций с гэпами составило 43, то есть 34.7% длины локального выравнивания.

Подсчёт инделей для локального выравнивания дал следующий результат:

RPOA_ECOLI 2
THIO_ECOLI 3
Total 5

Таким образом, в локальном выравнивании было найдено 5 инделей. Это меньше, чем в глобальном выравнивании, поскольку программа water не выравнивает последовательности целиком, а выбирает только наиболее подходящий участок.

5. Расчёт покрытия локальным выравниванием

В локальное выравнивание вошли следующие участки последовательностей:

• RPOA_ECOLI: позиции 130–233;
• THIO_ECOLI: позиции 5–81.

Длины покрытых участков были рассчитаны следующим образом:

• для RPOA_ECOLI: 233 − 130 + 1 = 104 аминокислотных остатка;
• для THIO_ECOLI: 81 − 5 + 1 = 77 аминокислотных остатков.

Покрытие последовательностей локальным выравниванием:

• Coverage 1 = 104 / 329 × 100% = 31.6%;
• Coverage 2 = 77 / 109 × 100% = 70.6%.

Из этих значений видно, что локальное выравнивание покрывает только около трети длинного белка RPOA_ECOLI, но большую часть короткого белка THIO_ECOLI. Это связано с существенным различием длин исходных последовательностей.

6. Сравнение глобального и локального выравнивания

Сравнение результатов показывает, что локальное выравнивание даёт более высокие значения identity и similarity, чем глобальное. Для глобального выравнивания identity составляет 7.4%, а для локального — 21.0%. Similarity увеличивается с 12.0% до 33.9%.

Однако это повышение не означает, что белки становятся гомологичными. Локальный алгоритм специально выбирает наиболее похожий фрагмент и отбрасывает остальные части последовательностей. Поэтому при сравнении неродственных белков он может находить участки случайного сходства. Такие участки выглядят лучше по численным характеристикам, чем глобальное выравнивание, но сами по себе не являются достаточным доказательством общего происхождения белков.

Глобальное выравнивание в данном случае оказалось малоинформативным, потому что оно содержит очень много гэпов и имеет крайне низкую идентичность. Локальное выравнивание оказалось более компактным и позволило выделить участок с относительно большим числом совпадений, но этот участок занимает только часть последовательностей и не позволяет уверенно говорить о родстве белков.

7. Вывод

В задании была проанализирована пара белков RPOA_ECOLI и THIO_ECOLI, имеющих разные функции и разные мнемоники. Глобальное выравнивание программой needle показало очень низкую идентичность (7.4%) и чрезвычайно высокую долю гэпов (74.5%), что указывает на отсутствие убедительного сходства по всей длине.

Локальное выравнивание программой water выявило ограниченный участок с более высокими значениями identity (21.0%) и similarity (33.9%), однако покрытие длинного белка RPOA_ECOLI составило только 31.6%. Поэтому данный результат следует интерпретировать осторожно: он показывает наличие локального сходства, но не доказывает гомологию белков.

Таким образом, на примере неродственных белков видно, что глобальное выравнивание может создавать формальное сопоставление с большим количеством гэпов, а локальное выравнивание может находить отдельные похожие участки. При этом для вывода о гомологии необходимо учитывать не только проценты identity и similarity, но также длину выравнивания, покрытие последовательностей, число гэпов и биологическую функцию сравниваемых белков.

6. Множественное выравнивание белков и импорт в Jalview

1. Выбор мнемоники и характеристика белка

Для выполнения данного задания была выбрана мнемоника функции 6PGD. Согласно базе данных UniProt/Swiss-Prot, для белка 6PGD_ECOLI рекомендуемое полное имя — 6-phosphogluconate dehydrogenase, decarboxylating.

Данный белок является ферментом, участвующим в окислительной фазе пентозофосфатного пути — важного метаболического процесса, обеспечивающего клетки восстановительным эквивалентом (NADPH) и предшественниками для синтеза нуклеотидов. Функциональная значимость этого фермента предполагает наличие консервативных участков, связанных с каталитической активностью.

Поиск в базе Swiss-Prot показал, что для данной мнемоники было найдено 53 белка, что указывает на широкую представленность данной функции в различных организмах.

2. Выбор последовательностей

Для проведения множественного выравнивания было выбрано 7 белков, включающих обязательные последовательности из организмов ECOLI и BACSU, а также дополнительные белки с той же мнемоникой функции:

• 6PGD_ECOLI
• 6PGD_BACSU
• 6PGD_STAAC
• 6PGD_SCHPO
• 6PGD_RAOTE
• 6PGD_HAEIN
• 6PGD_BACLI

Таким образом, были выполнены требования задания: использованы две обязательные последовательности и добавлены пять дополнительных гомологичных белков.

3. Методика проведения множественного выравнивания

Работа была выполнена в несколько последовательных этапов:

• Сформирован список идентификаторов белков Swiss-Prot;
• Получен FASTA-файл с аминокислотными последовательностями;
• Проведено множественное выравнивание с помощью программы muscle;
• Результат сохранён в формате FASTA;
• Файл импортирован в программу Jalview для визуального анализа.

В Jalview была применена раскраска по percent identity, что позволило визуально выделить консервативные и вариабельные области.

Гиперссылка на файл проекта Jalview:

Скачать файл

4. Анализ качества множественного выравнивания

Полученное множественное выравнивание включает 7 последовательностей белка 6-phosphogluconate dehydrogenase. Длина выравнивания составляет 495 позиций. При этом фактические длины белков без учёта гэпов различаются: для большинства последовательностей они находятся примерно в пределах 467–492 аминокислотных остатков, тогда как последовательность 6PGD_RAOTE короче и содержит 445 аминокислотных остатков, поскольку в файле она представлена как фрагмент.

Поэтому выравнивание нельзя описывать как полностью равномерное по всей длине. Основная часть последовательностей действительно хорошо сопоставляется, однако на концах выравнивания заметны участки с пропусками. Особенно это видно для 6PGD_RAOTE: у этой последовательности отсутствует начальный участок выравнивания примерно в позициях 1–15, а также большая концевая область примерно в позициях 464–495. Это согласуется с тем, что данная запись является фрагментом белка.

Кроме того, отдельные гэпы встречаются и внутри выравнивания. Они расположены, например, в областях около позиций 49, 169, 250 и 317. Эти пропуски не нарушают общую структуру выравнивания, но показывают, что последовательности не являются полностью одинаковыми по длине и содержат локальные вставки или делеции.

Таким образом, выравнивание имеет неоднородную структуру: центральная часть белков выровнена хорошо и содержит много консервативных позиций, тогда как начальные и концевые области более вариабельны и хуже сопоставимы между всеми последовательностями.

Гомология выбранных белков

Несмотря на наличие гэпов и различия в длине, выбранные белки можно считать гомологичными. Этот вывод основан не только на одинаковой мнемонике 6PGD, но и на структуре самого множественного выравнивания. В выравнивании присутствует протяжённое консервативное ядро, в котором аминокислотные остатки хорошо совпадают у большинства или у всех последовательностей.

Наиболее убедительно гомология подтверждается центральными участками выравнивания, где наблюдаются длинные блоки консервативных позиций. Например, хорошо сохраняются области примерно в позициях 182–202, 266–272, 286–293, 340–344, 370–374 и 449–460 выравнивания. В этих фрагментах многие остатки совпадают у большинства белков, а некоторые позиции являются полностью идентичными для всех семи последовательностей.

При этом степень сходства между последовательностями не одинакова. Например, 6PGD_ECOLI особенно близок к 6PGD_RAOTE, а также хорошо сопоставляется с 6PGD_BACSU и 6PGD_STAAC. Последовательности из более удалённых организмов, например 6PGD_SCHPO и 6PGD_BACLI, сохраняют общую структуру выравнивания, но имеют больше замен.

Структура выравнивания

По выравниванию можно выделить несколько характерных областей:

• Начальная область, примерно позиции 1–15: вариабельный участок. У последовательности 6PGD_RAOTE здесь стоят гэпы, так как она представлена фрагментом. У остальных белков начало присутствует, но также содержит несколько пропусков и замен.
• Область примерно 16–170: большая часть последовательностей хорошо сопоставляется, но встречаются отдельные локальные различия. В районе позиции 49 у нескольких белков присутствует короткий гэп, что указывает на небольшую вставку или делецию в этой части белка.
• Центральная область примерно 180–205: один из наиболее консервативных фрагментов выравнивания. В этом участке хорошо сохраняется мотив, включающий последовательность вида GAGH/VKMVHNGIEY. Наличие такого консервативного блока указывает на функциональную значимость данной области.
• Средняя часть примерно 245–322: сохраняется общий порядок аминокислотных остатков, но появляются локальные гэпы и замены. Это показывает, что участок остаётся гомологичным, но менее консервативным, чем центральное ядро.
• Область примерно 440–463: хорошо выраженный консервативный участок. Здесь у большинства последовательностей сохраняется фрагмент, близкий к LPANLIQAQRDYFGAHT. Этот участок можно считать одним из наиболее надёжно выровненных в данном множественном выравнивании.
• Концевая область примерно 464–495: наиболее вариабельный участок. У 6PGD_RAOTE эта область полностью отсутствует, а у ряда других последовательностей имеются длинные гэпы. Поэтому этот конец выравнивания нельзя считать одинаково хорошо представленным у всех белков.

5. Биологическая интерпретация

Полученные результаты соответствуют биологической роли белка 6PGD. 6-phosphogluconate dehydrogenase является ферментом, участвующим в пентозофосфатном пути, поэтому для него ожидаемо сохранение участков, важных для структуры белка и каталитической активности.

Наиболее консервативные области множественного выравнивания, особенно центральные фрагменты и участок около позиций 440–463, вероятно соответствуют функционально важным частям белка. Такие области хуже переносят аминокислотные замены, поэтому в ходе эволюции они сохраняются у разных организмов.

Менее консервативные участки, наоборот, расположены преимущественно на концах выравнивания и в отдельных внутренних областях с локальными гэпами. Эти различия могут отражать эволюционные вставки и делеции, а также то, что одна из выбранных последовательностей, 6PGD_RAOTE, является неполной фрагментарной записью.

Таким образом, множественное выравнивание показывает не равномерное совпадение по всей длине, а сочетание хорошо сохраняемого гомологичного ядра и более изменчивых периферических областей.

6. Итоговый вывод

В задании было построено множественное выравнивание семи белков с мнемоникой 6PGD. Выравнивание подтверждает, что выбранные последовательности относятся к одной гомологичной группе, так как они имеют общую функцию и содержат протяжённые консервативные участки.

Однако анализ FASTA-выравнивания показывает, что сходство распределено неравномерно. Центральная часть белков выровнена лучше всего и содержит много консервативных позиций, тогда как начальные и особенно концевые области более вариабельны. Наиболее заметное отличие связано с последовательностью 6PGD_RAOTE, которая является фрагментом и поэтому содержит длинные концевые гэпы.

Следовательно, корректный вывод состоит в том, что белки гомологичны и имеют сохраняемое функциональное ядро, но выравнивание не является одинаково качественным по всей длине. Наиболее надёжными для биологической интерпретации являются центральные консервативные участки, тогда как концевые области следует оценивать осторожно.

7. Параметры программ needle и water

Целью данной работы было выяснить, какие параметры запрашивают программы needle и water при запуске без опции -auto, понять биологический и алгоритмический смысл этих параметров, а также проверить на конкретных примерах, как изменение настроек влияет на итоговое выравнивание последовательностей.

Особенно важным условием задания было показать, что изменение параметров может приводить не просто к изменению численного значения score, а именно к перестройке самого выравнивания, то есть к изменению того, какие аминокислотные остатки сопоставляются друг другу, где появляются гэпы и какие участки последовательностей считаются оптимально выровненными.

Исходные данные и выбор белковых пар

Для выполнения задания были использованы две пары белков, различающиеся по степени родства и поэтому позволяющие проследить влияние параметров в разных ситуациях.

1. Гомологичная пара для исследования программы needle

В качестве гомологичной пары была выбрана пара 6PGL_ECOLI — 6PGL_BACSU. Эта пара использовалась ранее в заданиях 2–3 и оказалась удобной для анализа, поскольку для неё были получены не максимально высокие показатели сходства и покрытия:

• % identity = 30.6%
• % similarity = 48.7%
• coverage 1 = 79.8%
• coverage 2 = 71.6%
• indels = 6

Такие характеристики показывают, что белки, с одной стороны, действительно являются гомологичными, но, с другой стороны, их выравнивание не является полностью тривиальным и однозначным. Это делает пару особенно подходящей для изучения влияния штрафов за гэпы и параметров, связанных с концевыми разрывами. Если бы белки были практически идентичны по всей длине, изменение параметров могло бы не привести к заметной перестройке выравнивания. В данном случае предполагалось, что эффект будет выражен заметнее.

2. Неродственная пара для исследования программы water

В качестве неродственной пары была выбрана пара RPOA_ECOLI — THIO_ECOLI:

• RPOA_ECOLI — DNA-directed RNA polymerase subunit alpha
• THIO_ECOLI — Thioredoxin 1

Данные белки различаются по функции, а также по общей организации последовательностей, поэтому они хорошо подходят для исследования поведения программы water, выполняющей локальное выравнивание. Для неродственных белков особенно важно то, что алгоритм ищет не сходство по всей длине, а лишь лучший локальный участок. Поэтому в такой задаче изменение параметров способно заметно повлиять на то, какой именно участок последовательностей будет выбран как оптимальный.

Какие параметры запрашивают программы needle и water без -auto

При запуске программ needle и water без опции -auto они переходят в интерактивный режим и последовательно запрашивают у пользователя те параметры, которые необходимы для построения выравнивания, если они заранее не были указаны в командной строке.

В ходе выполнения работы было установлено, что в базовом случае обе программы запрашивают следующие параметры:

• Gap opening penalty — штраф за открытие гэпа;
• Gap extension penalty — штраф за продление уже открытого гэпа;
• Output alignment — имя файла, в который будет записан результат выравнивания.

Таким образом, уже при обычном интерактивном запуске становится видно, какие характеристики алгоритм считает принципиально важными для своей работы: стоимость появления разрывов и способ сохранения результата.

Смысл исследуемых параметров

1. Параметр gapopen

Параметр gapopen задаёт штраф за открытие нового гэпа, то есть за появление нового разрыва в одной из выравниваемых последовательностей. Биологически это можно интерпретировать как цену события вставки или делеции. С алгоритмической точки зрения этот параметр определяет, насколько охотно программа будет начинать новый гэп.

Если значение gapopen небольшое, алгоритм легче вводит новые гэпы. В этом случае выравнивание становится более гибким: программа чаще допускает вставки и удаления ради того, чтобы лучше сопоставить совпадающие или сходные аминокислоты в соседних участках. Если же значение gapopen велико, открытие каждого нового гэпа становится дорогим, и алгоритм начинает предпочитать несовпадения или менее выгодные сопоставления вместо создания дополнительных разрывов.

2. Параметр gapextend

Параметр gapextend задаёт штраф за продление уже открытого гэпа на один символ. Иными словами, после того как программа решила открыть гэп, она должна “платить” дополнительную цену за каждый следующий символ этого разрыва.

Небольшое значение gapextend делает длинные гэпы относительно выгодными: если разрыв уже открыт, его можно продолжать без очень больших дополнительных потерь. Напротив, высокое значение gapextend препятствует образованию длинных гэпов и делает алгоритм менее склонным к протяжённым вставкам и делециям.

Следовательно, параметры gapopen и gapextend работают совместно. Первый отвечает за вероятность появления новых гэпов, а второй — за их длину. Именно поэтому изменение этих параметров может менять не только общий балл выравнивания, но и его внутреннюю структуру.

3. Параметр outfile

Параметр outfile задаёт имя выходного файла, в который записывается результат выравнивания. Хотя он не влияет на сам алгоритм сопоставления последовательностей, этот параметр обязателен с практической точки зрения, поскольку без него невозможно сохранить результат и затем сравнить разные варианты выравнивания между собой.

4. Дополнительные параметры, исследованные в работе

Помимо параметров, непосредственно запрашиваемых в интерактивном режиме, в работе были исследованы и дополнительные настройки:

• для программы needle: -endweight, -endopen, -endextend;
• для программы water: -datafile, то есть выбор матрицы замен.

Эти параметры не случайны. Для глобального выравнивания принципиально важно, штрафуются ли концевые гэпы, а для локального выравнивания большое значение имеет система оценок замен аминокислот, поскольку именно она определяет, какие пары остатков алгоритм считает выгодными для сопоставления.

Краткое различие между программами needle и water

Программа needle реализует алгоритм глобального выравнивания. Это означает, что она стремится сопоставить последовательности по всей длине, от начала до конца, даже если отдельные участки совпадают слабо. Поэтому для needle особенно существенны штрафы за гэпы, включая концевые разрывы, поскольку такие штрафы напрямую влияют на то, как будут выровнены начало и конец последовательностей.

Программа water, напротив, выполняет локальное выравнивание. Она ищет не лучшее сопоставление по всей длине, а лучший локальный участок сходства. Поэтому для water изменение параметров может приводить к тому, что оптимальным окажется вообще другой фрагмент последовательностей. Это особенно важно при работе с неродственными или слабо сходными белками, где существует несколько возможных коротких участков умеренного сходства, и итоговый выбор сильно зависит от настроек алгоритма.

Ход работы и результаты

Эксперимент 1. Влияние gapopen и gapextend на результат needle для пары 6PGL_ECOLI / 6PGL_BACSU

В первом эксперименте исследовалось, как изменение штрафов за открытие и продление гэпов влияет на глобальное выравнивание гомологичной пары белков. Для этого были выполнены три запуска программы needle с разными значениями параметров:

needle sw:6pgl_ecoli sw:6pgl_bacsu 6pgl_needle_lowgap.fasta  -auto -aformat fasta -gapopen 3  -gapextend 0.1
needle sw:6pgl_ecoli sw:6pgl_bacsu 6pgl_needle_base.fasta    -auto -aformat fasta -gapopen 10 -gapextend 0.5
needle sw:6pgl_ecoli sw:6pgl_bacsu 6pgl_needle_highgap.fasta -auto -aformat fasta -gapopen 20 -gapextend 2.0

После получения результатов из файлов были удалены FASTA-заголовки, чтобы сравнивать только сами выровненные последовательности:

grep -v '^>' 6pgl_needle_base.fasta > a.txt
grep -v '^>' 6pgl_needle_lowgap.fasta > b.txt
grep -v '^>' 6pgl_needle_highgap.fasta > c.txt

diff -u a.txt b.txt
diff -u a.txt c.txt

Наблюдения

Сравнение показало, что изменение параметров gapopen и gapextend приводит к изменению самого выравнивания, а не только его итоговой численной оценки. Это означает, что при разных штрафах программа выбирает разные способы сопоставления аминокислотных позиций.

По результатам сравнения можно сделать несколько содержательных наблюдений. Во-первых, изменяется количество и расположение гэпов. Во-вторых, в некоторых участках один вариант выравнивания содержит более длинные разрывы, тогда как другой — несколько более коротких. В-третьих, происходят смещения в сопоставлении отдельных позиций, то есть аминокислоты, которые в одном варианте считались соответствующими друг другу, в другом варианте уже оказываются выровненными иначе.

Наиболее заметны такие различия в участках, где исходное сходство между белками не очень велико. Это закономерно: чем менее однозначен фрагмент, тем сильнее выбор оптимального решения зависит от параметров оценки. В хорошо консервативных участках выравнивание обычно остаётся стабильным, а в слабо сходных сегментах может заметно перестраиваться.

При низких штрафах (gapopen = 3, gapextend = 0.1) выравнивание становится более гибким. Программа легче вставляет гэпы и может активнее “подстраивать” одну последовательность под другую, чтобы добиться лучшего сопоставления локально похожих участков. Такой режим делает алгоритм более чувствительным к слабому локальному сходству, но одновременно может приводить к менее строгой и более “разреженной” структуре выравнивания.

При высоких штрафах (gapopen = 20, gapextend = 2.0) выравнивание становится более жёстким. Алгоритм старается избегать гэпов и поэтому чаще сохраняет несовпадения вместо того, чтобы открывать новые разрывы. В результате некоторые локально похожие участки оказываются сопоставлены менее точно, зато общая схема выравнивания становится более консервативной.

Вывод по эксперименту 1

Для пары 6PGL_ECOLI / 6PGL_BACSU параметры gapopen и gapextend существенно влияют на структуру глобального выравнивания. Следовательно, выбор штрафов за гэпы является не формальной технической деталью, а важной частью настройки алгоритма, определяющей, какие аминокислотные позиции будут считаться соответствующими друг другу.

Эксперимент 2. Влияние штрафов за концевые гэпы в needle

Во втором эксперименте для той же пары белков исследовалось влияние учёта концевых гэпов. Были выполнены два запуска программы needle:

needle sw:6pgl_ecoli sw:6pgl_bacsu 6pgl_needle_noend.fasta -auto -aformat fasta -gapopen 10 -gapextend 0.5
needle sw:6pgl_ecoli sw:6pgl_bacsu 6pgl_needle_end.fasta   -auto -aformat fasta -gapopen 10 -gapextend 0.5 -endweight Y -endopen 10 -endextend 0.5

Сравнение результатов проводилось следующим образом:

diff -u 6pgl_needle_noend.fasta 6pgl_needle_end.fasta

Наблюдения

Сравнение показало, что включение учёта концевых гэпов действительно влияет на структуру выравнивания. Изменения были заметны прежде всего в концевых областях последовательностей, а также в соседних с ними участках. Это означает, что алгоритм по-разному распределяет разрывы и несовпадения в зависимости от того, включаются ли концевые гэпы в итоговый штраф.

Такой результат является вполне ожидаемым с теоретической точки зрения. В глобальном выравнивании начало и конец последовательностей входят в общую задачу оптимизации. Если концевые гэпы не штрафуются, программа может относительно свободно “оставить” часть последовательности вне плотного сопоставления, не теряя за это баллы. Если же концевые гэпы начинают штрафоваться, алгоритм вынужден искать другое компромиссное решение, поскольку теперь такие свободные концы становятся менее выгодными.

Именно поэтому изменение параметров -endweight, -endopen и -endextend способно влиять не только на крайние позиции, но и на внутреннюю структуру выравнивания. В некоторых случаях перестройка в конце последовательности приводит к смещению сопоставления и в соседних участках, поскольку глобальное выравнивание оптимизируется как единая система.

Для выбранной пары этот эффект оказался хорошо заметным, что согласуется с ранее полученными значениями покрытия. Поскольку покрытие обеих последовательностей не является полным, концевые области не были полностью однозначны, и поэтому решение алгоритма оказалось чувствительным к тому, штрафуются ли такие участки.

Вывод по эксперименту 2

Параметры -endweight, -endopen и -endextend оказывают заметное влияние на результат работы программы needle. Для глобального выравнивания это принципиально важно, поскольку учёт или неучёт концевых гэпов способен менять не только общий балл, но и само сопоставление аминокислотных позиций.

Эксперимент 3. Влияние gapopen и gapextend на результат water для пары RPOA_ECOLI / THIO_ECOLI

В третьем эксперименте исследовалось влияние штрафов за гэпы на локальное выравнивание неродственных белков. Были выполнены следующие команды:

water sw:rpoa_ecoli sw:thio_ecoli rpoa_thio_water_base.fasta    -auto -aformat fasta -gapopen 10 -gapextend 0.5
water sw:rpoa_ecoli sw:thio_ecoli rpoa_thio_water_lowgap.fasta  -auto -aformat fasta -gapopen 2  -gapextend 0.1
water sw:rpoa_ecoli sw:thio_ecoli rpoa_thio_water_highgap.fasta -auto -aformat fasta -gapopen 20 -gapextend 2.0

Для сравнения был использован следующий вызов:

diff -u rpoa_thio_water_base.fasta rpoa_thio_water_lowgap.fasta

Наблюдения

Сравнение показало, что между базовым вариантом и вариантом с низкими штрафами за гэпы имеются выраженные различия. При уменьшении значений gapopen и gapextend изменяется не только общий характер выравнивания, но и сам локальный участок, который программа считает наилучшим.

По результатам сравнения можно сделать вывод, что меняется длина локально выровненного сегмента, расположение гэпов и набор аминокислотных фрагментов, включённых в выравнивание. Иначе говоря, алгоритм находит другое решение задачи. Это особенно важно для локального выравнивания, поскольку его цель состоит в поиске одного наилучшего сходного участка, а не в сопоставлении последовательностей по всей длине.

В базовом варианте локальное выравнивание соответствует одному набору фрагментов, тогда как при более низких штрафах программа может позволить себе более свободно вводить гэпы и тем самым соединять сходные позиции в более протяжённый, но менее жёсткий по структуре участок. Для неродственных белков такая ситуация особенно показательна, поскольку здесь нет длинного устойчивого гомологичного фрагмента, который оставался бы неизменным при любых настройках.

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

Вывод по эксперименту 3

Для пары RPOA_ECOLI / THIO_ECOLI изменение параметров gapopen и gapextend приводит к изменению самого локального выравнивания. Следовательно, для программы water штрафы за гэпы влияют не только на численную оценку, но и на выбор конкретного фрагмента последовательностей, который будет признан оптимальным.

Эксперимент 4. Влияние матрицы замен на результат water

В четвёртом эксперименте для той же пары неродственных белков исследовалось влияние матрицы замен. Были выполнены локальные выравнивания с разными матрицами:

water sw:rpoa_ecoli sw:thio_ecoli rpoa_thio_b62.fasta -auto -aformat fasta -datafile EBLOSUM62 -gapopen 10 -gapextend 0.5
water sw:rpoa_ecoli sw:thio_ecoli rpoa_thio_b45.fasta -auto -aformat fasta -datafile EBLOSUM45 -gapopen 10 -gapextend 0.5
water sw:rpoa_ecoli sw:thio_ecoli rpoa_thio_b80.fasta -auto -aformat fasta -datafile EBLOSUM80 -gapopen 10 -gapextend 0.5

Сравнение проводилось так:

diff -u rpoa_thio_b62.fasta rpoa_thio_b45.fasta
diff -u rpoa_thio_b62.fasta rpoa_thio_b80.fasta

Наблюдения

Сравнение показало, что замена матрицы EBLOSUM62 на EBLOSUM45 и EBLOSUM80 приводит к изменению локального выравнивания. Различия касаются длины выровненного фрагмента, положения гэпов и состава тех аминокислотных позиций, которые алгоритм считает выгодным сопоставить.

Причина такого эффекта заключается в том, что матрица замен задаёт систему оценок совпадений и замен аминокислот. Разные матрицы по-разному “оценивают” одни и те же пары остатков. Следовательно, при изменении матрицы меняется и функция качества выравнивания, а значит, изменяется и расположение оптимума.

Для неродственной пары этот эффект проявляется особенно наглядно. Поскольку локальное сходство между такими белками само по себе слабое и неоднозначное, даже небольшое изменение системы оценок может сместить оптимальное решение и привести к выбору другого участка. Именно поэтому выбор матрицы замен нельзя считать второстепенной технической настройкой: он непосредственно определяет биологическую интерпретацию результата.

Вывод по эксперименту 4

Матрица замен является одним из ключевых параметров локального выравнивания. Для пары RPOA_ECOLI / THIO_ECOLI изменение матрицы приводит к изменению самого выровненного участка, а не только итогового score. Это подтверждает, что при анализе слабого или сомнительного сходства выбор матрицы имеет принципиальное значение.

Обсуждение результатов

Проведённые эксперименты показали, что параметры программ needle и water действительно влияют на результат выравнивания, причём это влияние не сводится к формальному изменению численного балла. Во всех исследованных случаях было показано, что изменение параметров может менять саму структуру выравнивания.

Что показали эксперименты с needle

На гомологичной паре 6PGL_ECOLI / 6PGL_BACSU было установлено, что:

• изменение штрафов за открытие и продление гэпа меняет структуру глобального выравнивания;
• учёт концевых гэпов также способен заметно повлиять на итоговое сопоставление позиций;
• наиболее чувствительными к изменениям параметров оказываются участки с неполным покрытием и умеренным сходством.

Это особенно важно для практической биоинформатики, потому что глобальное выравнивание нередко воспринимается как нечто однозначное. На самом деле даже для гомологичных белков итоговый результат зависит от выбранной модели штрафов, особенно если белки совпадают не по всей длине и содержат вставки, делеции или слабо консервативные области.

Что показали эксперименты с water

На неродственной паре RPOA_ECOLI / THIO_ECOLI было установлено, что:

• изменение gapopen и gapextend меняет найденный локальный участок;
• изменение матрицы замен также приводит к перестройке локального выравнивания;
• для неродственных белков локальное выравнивание особенно чувствительно к настройкам, поскольку устойчивого длинного гомологичного участка у них нет.

Таким образом, для локального выравнивания параметры оказываются в некотором смысле даже более критичными, чем для глобального. Они влияют не только на степень “красоты” выравнивания, но и на сам выбор того фрагмента, который будет интерпретирован как наиболее сходный.

Методический комментарий к способу сравнения результатов

Для сравнения результатов использовалась выдача в формате FASTA (-aformat fasta), после чего выравнивания сравнивались напрямую с помощью программы diff. Такой способ оказался удобным и методически оправданным, поскольку позволял оценивать не только численные характеристики, но и структуру самих выравниваний.

Это особенно важно в контексте данного задания. Если бы анализ ограничивался только сравнением score, нельзя было бы с уверенностью утверждать, что параметры действительно меняют биологически значимый результат. Использование diff позволило показать, что в разных вариантах изменяется расположение гэпов, длина выровненных участков и сопоставление конкретных аминокислотных позиций. Следовательно, главное условие задания было выполнено: различие результатов подтверждалось различием самих выравниваний.