Моделирование и реконструкция эволюции гена.

Заданное мне дерево (модель эволюции гена KAD_ECOLI) в скобочной модели*

(А:120,((((B:10,C:10):50,D:65):5,E:70):5,F:70):30);

*Расстояния даны как число мутаций на 100 нуклеотидных остатков.

• Изображение дерева, описанного этой формулой:

Листья и длины ветвей подписаны на рисунке. Кроме того, на рисунке добавлены внутренние узлы H,G,T,P, не описанные формулой. Эти же обозначения использовались при написании скрипта.
Данное дерево - укоренённое неультраметрическое.

• Описание ветвей дерева как разбиения множества листьев

Считаем дерево бескорневым:

A	B	C	D	E	F
-	-	-	*	-	-
-	*	*	-	-	-
-	*	*	*	-	-
-	-	-	-	*	-
*	-	-	-	-	*
-	-	-	-	-	*
*	-	-	-	-	-
-	*	-	-	-	-
-	-	*	-	-	-

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

A	B	C	D	E	F
-	*	*	-	-	-
-	*	*	*	-	-
*	-	-	-	-	*

•  Получение искуственных  мутантных последовательностей

Считаем, что в корне находится последовательность гена белка KAD_ECOLI, полученная при выполнении второго упражнения занятия "Банк EMBL" третьего семестра. Длина гена - 645 нуклеотидов. Поэтому формула для вычисления необходимого числа мутаций имеет вид:

645 * R/ 100

 Где R - это длина ветви (эволюционное расстояние, измеряющееся в числе мутаций, приходящемся на сто нуклеотидов).
Для того, чтобы воспользоваться программой msbar пакета EMBOSS, был написан скрипт, строки которого имеют вид:

msbar infile outfile -point 4 -count n -auto

Где infile - имя файла с последовательностью,в которую надо внести мутации, а outfile - имя файла, в котором будет находится измененная последовательность, n - число мутаций, которое нужно внести. Опция "-point" позволяет выбрать вид точечных мутаций, которые нужно "применить" к данной последовательности. Можно выбрать следующие значения: 0 - никакие; 1 - любые из перечисленных ниже; 2 - вставки; 3 - делеции; 4 - замены; 5 - дупликации; 6 - перемещения.
Файл сохранён как текст UNIX.

• Реконструкция дерева алгоритмами UPGMA, Neighbor-joining и максимального правдоподобия.  

Для этого был создан файл в fasta-формате, содержащий последовательности, соответствующие листьям дерева в виде "выравнивания". Далее, при помощи полученного файла, были применены следующие алгоритмы построения деревьев:

1. Алгоритм максимального правдоподобия

fdnaml offsprings.fasta -ttratio 1 -auto

Подробнее об опциях программы fdnaml написано здесь. Опция "-ttratio" определяет коэффициент отношения числа транзиций к трансверсиям. В данном случае, число транзиций принимается равным числу трансверсий, то есть мы считаем эти события равновероятными. Результатом работы команды стали два файла: с расширением .treefile и с расширением .fdnaml. Файл с разрешением .fdnaml содержит "текстово-графическое" изображение дерева:

+------------------F +--4 +--3 +-----------------E | | | +--------------D | | +--B 2---------1 | +-C | +-----------------------------------A

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

2. Алгоритм Neighbor-joining

Сначала были посчитаны попарные расстояния между последовательностями программой fdnadist:

fdnadist offsprings.fasta -ttratio 1 -auto

Получился файл с расширением .fdnadist, который нужно подать на вход программе fneighbor:

fneighbor offsprings.fdnadist -auto

Результат работы программы:

+--B            +---------1   !         +-C            !   !  +----------------D            2--3   !  ! +------------------E            !  +-4   !    +------------------F            !   +-------------------------------------A

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

3. Алгоритм UPGMA

Использовалась команда:

fneighbor offsprings.fdnadist -treetype u -auto

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

+------------------------------------------------------A           !   !                                                  +---B           !                        +-------------------------1 --5                    +---2                         +---C           !                    !   !   !                 +--3   +-----------------------------D           !                 !  !   +-----------------4  +---------------------------------E                             !                     +------------------------------------F

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

• Анализ полученных данных

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

Алгоритм	Рисунок дерева в бескорневом виде	Описание ветвей дерева как разбиения множества листьев
Исходное дерево		ABCDEF -**--- -***-- *----*
Максимального правдоподобия		ABCDEF -**--- ---*** ----**
Neighbor-joining		ABCDEF -**--- ---*** ----**
UPGMA		ABCDEF -**--- -***-- *----*

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

Дополнительное задание

Сравнение разных способов оценки эволюционных расстояний между генами.

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

• Определить "истинные" эволюционные расстояния между конечными мутантами в cозданной модели эволюции (т.е. расстояния в числе замен между листьями заданного дерева)

Для выполнения этого задания была создана рабочая книга Excel "Dist.xls", содержащая 2 листа c названиями "All_data" и "Comparison".
На стр. "All_data" размещена таблица с "истинными" попарными расстояниями (число замен на 100 нуклеотидов):

узел	A	B	C	D	E	F	P	T	G	H
A	0	220	220	225	225	220	210	160	155	150
B		0	20	125	135	140	10	60	65	70
C			0	125	135	140	10	60	65	70
D				0	140	145	115	65	70	75
E					0	145	125	75	70	75
F						0	130	80	75	70
P							0	50	55	60
T								0	5	10
G									0	5
H										0

 Названия узлов соответствуют обозначениям первого рисунка.

• Определить попарные эволюционные расстояния между всеми последовательностями (включая исходную) c помощью программы distmat пакета EMBOSS.

На вход программе подается множественное выравнивание, но т.к. в исходной модели были только замены, то просто файл с последовательностями и соответствует биологически значимому выравниванию.
Методы оценки эволюционных расстояний программой distmat:
         0 : Uncorrected
         1 : Jukes-Cantor
         2 : Kimura
         3 : Tamura
         4 : Tajima-Nei
         5 : Jin-Nei Gamma
Выбирая соответствующие пункты меню программы, получены 3 матрицы попарных расстояний:

1. матрица попарных различий (D) (uncorrected distances);

distmat -sequence offsprings.fasta -outfile offsprings.distmat -nucmethod 0

Получившаяся матрица попарного различия: среднее число несовпадающих нуклеотидов на 100 позиций:

1 2 3 4 5 6 0.00 65.89 65.27 67.44 68.06 67.44 A 1 0.00 12.40 53.18 56.59 58.76 B 2 0.00 53.18 54.73 58.91 C 3 0.00 56.59 57.52 D 4 0.00 58.29 E 5 0.00 F 6

Матрица перенесена на лист "All_Data".

2. матрица попарных расстояний, вычисленных по формуле Джукса - Кантора (JC)

distmat -sequence offsprings.fasta -outfile offsprings.distmat -nucmethod 1

Матрица попарных эволюционных расстояний, вычисленных по методу Джукса – Кантора:

1 2 3 4 5 6 0.00 158.12 153.18 172.11 178.54 172.11 A 1 0.00 13.56 92.59 105.34 114.75 B 2 0.00 92.59 98.12 115.47 C 3 0.00 105.34 109.23 D 4 0.00 112.63 E 5 0.00 F 6

Матрица перенесена на лист "All_Data".

3. матрица попарных расстояний, вычисленных по модели М.Кимура.

distmat -sequence offsprings.fasta -outfile offsprings.distmat -nucmethod 2

Матрица попарных эволюционных расстояний, вычисленных по методу М.Кимура:

1 2 3 4 5 6 0.00 161.10 156.23 172.41 178.55 172.94 A 1 0.00 13.56 92.60 105.52 114.75 B 2 0.00 92.68 98.20 115.51 C 3 0.00 106.77 109.25 D 4 0.00 112.78 E 5 0.00 F 6

Основным отличием модели М.Кимура от модели Джукса - Кантора является то, что транзиция и трансверсия разновероятны.
Матрица перенесена на лист "All_Data".

• Сравнение методов оценки эволюционных расстояний, предложеннных программой distmat

Была составлена сравнительная таблица, отсортированная по убыванию "истинных расстояний":

	"истинные" попарные расстояния	Матрица попарного различия	метод Джукса – Кантора	метод М.Кимура
AD	225	67,44	172,11	172,41
AE	225	68,06	178,54	178,55
AB	220	65,89	158,12	161,1
AC	220	65,27	153,18	156,23
AF	220	67,44	172,11	172,94
DF	145	57,52	109,23	109,25
EF	145	58,29	112,63	112,78
BF	140	58,76	114,75	114,75
CF	140	58,91	115,47	115,51
DE	140	56,59	105,34	106,77
BE	135	56,59	105,34	105,52
CE	135	54,73	98,12	98,2
BD	125	53,18	92,59	92,6
CD	125	53,18	92,59	92,68
BC	20	12,4	13,56	13,56

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

Из графика видно, что кривые, соответствующие методам Джукса-Кантора и М.Кимура практически полностью совпадают. Для того, чтобы различия между методами были более заметны, построен график, на котором по оси абсцисс отложен условный номер пары расстояний, а по оси ординат - нормированные значения чисел мутаций (см лист "Comparison_2"). По совпадению с кривой "истинных" значений можно судить о точности метода.
 
По графику можно заметить, что хоть кривые, отображающие методы Джукса-Кантора и М.Кимура, почти всё время совпадают, зелёная кривая всё же выше малиновой, а, значит, оценки эволюционных расстояний, полученные методом М.Кимура ближе к "истинным", чем, полученные методом Джукса-Кантора. Такой результат объясняется отличием этих моделей*: в модели Джукса - Кантора транзиции и трансверсии считаются равновероятными, а в модели М.Кимура - нет (транзиции более вероятны). В связи с этим получается, что метод оценки эволюционного расстояния между нуклеотидными последовательностями М.Кимура более верен, потому что в его основе лежит более точная модель. И хотя на данном примере различия между этими методами минимальны, при увеличении числа мутаций - эволюционного расстояния, погрешность метода Джукса - Кантора будет накапливаться и, в силу этого, данный метод будет неприменим.
Метод попарных различий показал наихудший результат - его кривая почти параллельна оси абсцисс и совпадает с кривой "истинных" значений почти только в 1 точке.
Интересно отметить, что с увеличением номера пары все три метода стремятся в одну точку и при этом почти совпадают с "истинным" значением. Это объясняется тем, что в силу сортировки по "истинным" расстояниям первые рассмотренные пары находятся далеко друг от друга (см. исходное дерево и сравнительную таблицу), а последующие пары рассматривались с меньшим расстоянием. Отсюда делаем вывод, что чем меньше эволюционное расстояние между потомками, тем более точны все 3 рассмотренных метода, однако при увеличении эволюционнного расстояния они сильно расходятся как друг от друга, так и от "истинной" картины. Иными словами, получили такую картинку, только отражённую относительно оси ординат:

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

Способ оценки эволюционных расстояний	Диапозон близких** к "истинному" значений	Диапозон максимально отличающихся*** от "истинного" значений
Метод попарных различий	до 20-25 замен на 100 нуклеотидов	от 125 замен на 100 нуклеотидов
Метод Джукса-Кантора	до 20-25 замен на 100 нуклеотидов	от 220 замен на 100 нуклеотидов
Метод М.Кимура	до 20-25 замен на 100 нуклеотидов	от 222 замен на 100 нуклеотидов

Вывод:
Значит, для очень "дальних родственников" лучше использовать метод М.Кимура, для более близких последовательностей - Джукса - Кантора, а для очень близких - матрицу попарных различий, поскольку этот метод не требует особых компьютерных ресурсов.


* Сходство моделей Джукса-Кантора и М.Кимура:

Модель Джукса-Кантора	Модель М.Кимура
Расстояние = -b ln (1-D/b) D - число мутаций b - константа. b= 3/4 для нуклеотидных последовательностей и 19/20 - для белковых.	Расстояние = -0.5 ln[ (1-2P-Q)*sqrt(1-2Q)] P = число транзиций/npos Q = число трансляций/npos npos - число вычисляемых позиций

Как видно из формул расчёта эволюционных расстояний, применяемых в рассмотренных моделях, общий вид формулы имеет вид:

---

Главная  Первый семестр  Второй семестр  Третий семестр