АвторТема: Ручное построение дерева (wertner)  (Прочитано 12672 раз)

0 Пользователей и 1 Гость просматривают эту тему.

Оффлайн CenturionАвтор темы

  • 100% Earth (Solar System) genofond
  • Администратор
  • *****
  • Сообщений: 9548
  • Страна: ru
  • Рейтинг +571/-2
Ручное построение дерева (wertner)
« : 17 Август 2009, 15:59:47 »
Пост Вадима (wertner) от 10.11.2008 (с Родства)

---------------------------------------------------------------------------------------

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

Свои формулы буду демонстрировать на МакДональдах R1a, потомках Сомерленда (http://dna-project.clan-donald-usa.org/DNAresults.htm). В один пост не управлюсь, поэтому буду рассказывать поэтапно.

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

Первой парой будут W.O. McDaniel(&2QF9S) и D. McDaniel(&642PV), которые отличаются друг от друга на мутацию в маркере CDY b. Второй парой будут W.H. McDaniel (&SVKFL) и B.G. McDaniel (&PDDRM), которые отличаются друг от друга на мутацию в маркере DYS 537.
Скорости мутаций CDY b и DYS 537 по Чандлеру равняются 0,03531 и 0,00057 мутаций на поколение и отличаются в 60 раз.

Скорей всего, гаплотип общего предка первой пары совпадает с одним из их двух гаплотипов и соответственно, отличается от другого на одну мутацию.
Рассмотрим такой набор гипотез: гипотеза 1 – это предположение, что их общие предок жил поколение назад, т.е. был их отцом, гипотеза 2 – то, что общие предок жил два поколения назад и так далее, в общем случае, гипотезы N, что ближайший общий предок жил N поколений назад.

Для каждой гипотезы посчитаем вероятность ее осуществления, пока не используя факт, что одна из них точно осуществилась. Т.е. априорную вероятность.
Для гипотезы 1 (ближайший общий предок - отец): вероятность того, что за одно поколение не произошло мутаций равняется произведению 67 вероятностей не-мутации (т.е. не-мутации) в каждом маркере. Не буду приводить весь ряд, но вот его начало и конец [Вероятность отсутствия мутаций в одном поколении] = (1-0,00076)*(1-0,00311)*…*(1-0,00087) = 0,797746589733803. Вероятность того, что произошла мутация в маркере CDY b и только в нем равна произведению вероятностей мутации в маркере CDY b (0,03531) и вероятностей не-мутации в 66 остальных маркерах: (1-0,00076)*(1-0,00311)*…* 0,03531*…*(1-0,00087) = 0,029199465.
Вероятность того, что у отца двух сыновей у одного сына не было мутаций, а у другого сына была мутация в маркере CDY b равна произведению этих вероятностей: 0,797746589733803 * 0,029199465 = 0,023293774.

Для гипотезы 2 (ближайший общий предок - дед): вероятность того, что за два поколения не произошло мутаций равна квадрату вероятности отсутствия мутаций в одном поколении: 0,797746589733803*0,797746589733803. Что же касается мутации в CDY b, то она могла произойти два поколения назад или одно поколение назад – тогда у общего предка (деда) были два сына с совпадающим с ним гаплотипом, а у одного из внуков уже произошла мутация. Т.е. генеалогические линии уже разошлись, но мутация наступила позже. Вероятность того, что за два поколения произошла мутация в CDY b будет произведением вероятности того, что в одном поколении (от деда к отцу) была мутация в маркере CDY b и только в нем, вероятности того, что в одном поколении (от отца к сыну) не было мутаций, т.е. 0,029199465*0,797746589733803 и количества таких вариантов = 2 (мутация была в первом поколении или во втором). Общая вероятность гипотезы: 0,797746589733803*0,797746589733803*0,029199465*0,797746589733803*2 = 0,029648298.

Аналогично подсчитываем вероятность для гипотезы 3 (ближайший общий предок - прадед): вероятность того, что за три поколения не произошло мутаций равна кубу вероятности отсутствия мутаций в одном поколении: 0,797746589733803*0,797746589733803*0,797746589733803. Что же касается мутации в CDY b, то она могла произойти три, два, одно поколение назад. Вероятность того, что за три поколения назад произошла мутация в CDY b будет произведением вероятности того, что в одном поколении (от прадеда к деду) была мутация в маркере CDY b и только в этом маркере, вероятности того, что в двух поколениях не было мутаций количества таких вариантов (3), т.е. 0,029199465*0,797746589733803*0,797746589733803*3. Общая вероятность гипотезы: 0,797746589733803*0,797746589733803*0,797746589733803*0,029199465*0,797746589733
803*0,797746589733803*3 = 0,028302248.

Аналогичны вычисления для остальных гипотез. Обращаю внимание, что при переходе к каждой последующей только добавляется умножение на квадрат вероятности не-мутирования гаплотипа (0,797746589733803*0,797746589733803) и меняется коэффициент равный количеству возможных моментов мутаций на одной из ветвей (1,2,3 и т.д. для соответствующих гипотез). Эти изменения не зависят от скорости мутации маркера CDY b.

Известно, что одна из этих гипотез воплотилась. Значит, для расчета вероятности гипотезы с учетом этого факта (т.е. апостериорной вероятности) по теореме Байеса (http://ru.wikipedia.org/wiki/%D0%A2%D0%B5%D0%BE%D1%80%D0%B5%D0%BC%D0%B0_%D0%91%D0%B0%D0%B9%D0%B5%D1%81%D0%B0) надо поделить априорную вероятность гипотезы на сумму априорных вероятностей всех гипотез.
Я не стал выводить/искать формулу суммы такого ряда, а просто взял сумму вероятностей гипотез от 1 до 300 (априорная вероятность гипотезы, что общий предок жил 300 лет назад равнялась 8,19622E-58).
Сумма всего вероятностей этих 300 гипотез равняется 0,176194034.

Теперь делим априорные вероятности на это число, и получаем соответственно нужные нам апостериорной вероятности:
Гипотеза 1: 0,023293774 / 0,176194034 = 0,132205235 = 13,22%
Гипотеза 2: 0,029648298 / 0,176194034 = 0,168270723 = 16,83%
Гипотеза 3: 0,028302248 / 0,176194034 = 0,160631137 = 16,06%
Гипотеза 4: 0,024015387 / 0,176194034 = 0,136300793 = 13,63%
Гипотеза 5: 0,019104229 / 0,176194034 = 0,108427216 = 10,84%
Гипотеза 6: 0,014589509 / 0,176194034 = 0,082803647 = 8,28%
Гипотеза 7: 0,010832217 / 0,176194034 = 0,061478911 = 6,15%
Гипотеза 8: 0,007878422 / 0,176194034 = 0,044714464 = 4,47%
Гипотеза 9: 0,005640553 / 0,176194034 = 0,032013301 = 3,20%
Гипотеза 10: 0,003988495 / 0,176194034 = 0,022636948 = 2,26%
Гипотеза 11: 0,002792104 / 0,176194034 = 0,015846759 = 1,58%
Гипотеза 12: 0,00193843 / 0,176194034 = 0,011001678 = 1,10%

Мат. ожидание - 4,50, медиана между 3 и 4 поколениями.
Теперь проведем такие же вычисления для второй пары, с мутацией в маркере DYS 537. Скорость мутации маркера DYS 537 по Чандлеру равняется 0,00057.
Значит вероятность того, что в одном поколении мутация произошла в маркере DYS 537 и только в нем равняется (1-0,00076)*(1-0,00311)*…* 0,00057 *…*(1-0,00087) = 0,000454975.
Соответственно, априорная вероятность гипотезы 1 вычисляется как 0,797746589733803 * 0,000454975 = 0,000362955. Априорная вероятность гипотезы 2 совершенно также вычисляется как 0,797746589733803 * 0,797746589733803 * 0,797746589733803 * 0,000454975 * 2 = 0,000461968. Ну и априорная вероятность гипотезы 3 равняется 0,797746589733803 * 0,797746589733803 * 0,797746589733803 * 0,797746589733803 * 0,797746589733803 * 0,000454975 * 3 = 0,000440995.

Сумма вероятностей всех 300 независимых априорных гипотез равна 0,002745388.

Соответственно, вычисляем апостериорные вероятности для второй пары:
Гипотеза 1: 0,000362955 / 0,002745388 = 0,132205235 = 13,22%
Гипотеза 2: 0,000461968 / 0,002745388 = 0,168270723 = 16,83%
Гипотеза 3: 0,000440995 / 0,002745388 = 0,160631137 = 16,06%
Гипотеза 4: 0,000374199 / 0,002745388 = 0,136300793 = 13,63%
Гипотеза 5: 0,000297675 / 0,002745388 = 0,108427216 = 10,84%
Гипотеза 6: 0,000227328 / 0,002745388 = 0,082803647 = 8,28%
Гипотеза 7: 0,000168783 / 0,002745388 = 0,061478911 = 6,15%
Гипотеза 8: 0,000122759 / 0,002745388 = 0,044714464 = 4,47%
Гипотеза 9: 8,78889E-05 / 0,002745388 = 0,032013301 = 3,20%
Гипотеза 10: 6,21472E-05 / 0,002745388 = 0,022636948 = 2,26%
Гипотеза 11: 4,35055E-05 / 0,002745388 = 0,015846759 = 1,58%
Гипотеза 12: 3,02039E-05 / 0,002745388 = 0,011001678 = 1,10%

Мат. ожидание - 4,50 поколений, медиана между 3 и 4 поколениями.

Итого имеем, что в обоих случаях вероятность того, что ближайший общий предок был одно поколение назад – 13,22%, два поколения назад – 16,83%, три поколения назад – 16,06%, 4 поколения назад – 13,63% и так далее.

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

Оффлайн CenturionАвтор темы

  • 100% Earth (Solar System) genofond
  • Администратор
  • *****
  • Сообщений: 9548
  • Страна: ru
  • Рейтинг +571/-2
Re: Ручное построение дерева (wertner)
« Ответ #1 : 17 Август 2009, 16:00:26 »
Продолжение...

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

Например, оценим возраст ближайшего общего предка для пар с двумя мутациями между гаплотипами. Причем отдельно для пары, где ближайший общий предок отличается на одну мутацию от каждого из потомков, и для пары, где гаплотип ближайшего общего предка совпадает с гаплотипом одного из потомков и отличается на две мутации от гаплотипа другого потомка. В построенном мной дереве МакДональдов (его представлю в конце своего рассказа) первой из этих пар являются D.B. MacDonald (&LT9RK) и R.R. MacDonald (&IEBH9) – их ближайший общий предок отличается от D.B. MacDonald'а на одну мутацию в маркере GATA H4 и от R.R. MacDonald’а на одну мутацию в маркере CDY b. Для примера второй пары подходят L.W. MacDonald (&R0PJB) и R.A. MacDonald (&EDC4L) – гаплотип их ближайшего общего предка совпадает (по версии моего дерева) с гаплотипом L.W. MacDonald'а и отличается на единичные мутации в маркерах DYS 449 и DYS 576 от гаплотипа R.A. MacDonald’а.

Подсчитаем для первой пары априорную вероятность гипотезы 1 (что ближайший общий предок - отец). Вероятность того, что произошла одна мутация в маркере GATA H4 (скорость мутации по Чандлеру 0,00208) и только в нем равна (1-0,00076)*(1-0,00311)*…* 0,00208 *…*(1-0,00087) = 0,001662771. Вероятность мутации за одно поколение в маркере CDY b и только в нем уже вычислена нами в предыдущем примере и равна 0,029199465. Их произведение и будет вероятностью гипотезы 1. Итого: 0,001662771*0,029199465 = 0,000048552.

При расчете гипотезы 2 подсчитаем вероятность того, что за два поколения произошла одна мутация в маркере GATA H4 и только в нем, равна произведению вероятности не-мутации гаплотипа, вероятность мутации за одно поколение только в маркере GATA H4 и количеству вариантов момента мутации (одно или два поколения назад) равным 2. Итого: 0,797746589733803*0,001662771*2 = 0,002652941. Вероятность мутации за два поколения в маркере CDY b и только в нем уже вычислена нами в предыдущем примере и равна 0,046587548. Итого: 0,002652941*0,046587548 = 0,000123594.

Для гипотезы 3 аналогично получаем априорную вероятность равную 0,000176974.

Рассчитав априорные вероятности для 300 гипотез и поделив отдельно априорные вероятности гипотез на их сумму (0,001652814) получаем апостеорные вероятности:
Гипотеза 1: 0,000048552 / 0,001652814 = 0,029375388 = 2,94%
Гипотеза 2: 0,000123594 / 0,001652814 = 0,074777943 = 7,48%
Гипотеза 3: 0,000176974 / 0,001652814 = 0,107074473 = 10,71%
Гипотеза 4: 0,000200224 / 0,001652814 = 0,121141607 = 12,11%
Гипотеза 5: 0,000199098 / 0,001652814 = 0,120460114 = 12,05%
Гипотеза 6: 0,000182457 / 0,001652814 = 0,11039151 = 11,04%
Гипотеза 7: 0,000158046 / 0,001652814 = 0,095622296 = 9,56%
Гипотеза 8: 0,00013137 / 0,001652814 = 0,079482766 = 7,95%
Гипотеза 9: 0,000105811 / 0,001652814 = 0,064018859 = 6,40%
Гипотеза 10: 8,31336E-05 / 0,001652814 = 0,050298244 = 5,03%
Гипотеза 11: 6,40165E-05 / 0,001652814 = 0,038731838 = 3,87%
Гипотеза 12: 4,84841E-05 / 0,001652814 = 0,02933426 = 2,93%

Мат. ожидание - 6,64 поколений, медиана между 5 и 6 поколениями.

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

Для вычисления априорной вероятности гипотезы 1 нам уже известна вероятность первой ветви – вероятность того, что мутация не произошла в одном поколении, 0,797746589733803. Вероятность того, что мутации произошли в маркерах DYS 449 (скорость 0,00838) и DYS 576 (скорость 0,01022) и только в них равна (1-0,00076)*(1-0,00311)*…*0,00838*… *0,01022*…*(1-0,00087) = 0,0000696107. Соответственно, априорная вероятность гипотезы 1 это 0,797746589733803*0,0000696107 = 0,0000555317.

Вычисление вероятности гипотезы 2 этой пары отличается от вычисления вероятностей гипотез 2 других пар только количеством возможных вариантов момента мутации: каждая из двух мутаций может произойти в первом или втором поколении, а значит количество вариантов равно квадрату двойки – четырем. Итого: 0,797746589733803*0,797746589733803*0,0000555317*0,797746589733803*4 = 0,000141361.

Для гипотезы 3 количество вариантов равняется 3*3=9, а значит вся вероятность вычисляется как 0,797746589733803*0,797746589733803*0,797746589733803*0,0000555317*0,79774658973
3803*0,797746589733803*9 = 0,000202415.

После вычисления вероятностей для 300 гипотез вычислим апостеорные вероятности:
Гипотеза 1: 0,0000555317 / 0,001890415 = 0,029375388 = 2,94%
Гипотеза 2: 0,000141361 / 0,001890415 = 0,074777943 = 7,48%
Гипотеза 3: 0,000202415 / 0,001890415 = 0,107074473 = 10,71%
Гипотеза 4: 0,000229008 / 0,001890415 = 0,121141607 = 12,11%
Гипотеза 5: 0,00022772 / 0,001890415 = 0,120460114 = 12,05%
Гипотеза 6: 0,000208686 / 0,001890415 = 0,11039151 = 11,04%
Гипотеза 7: 0,000180766 / 0,001890415 = 0,095622296 = 9,56%
Гипотеза 8: 0,000150255 / 0,001890415 = 0,079482766 = 7,95%
Гипотеза 9: 0,000121022 / 0,001890415 = 0,064018859 = 6,40%
Гипотеза 10: 9,50846E-05 / 0,001890415 = 0,050298244 = 5,03%
Гипотеза 11: 7,32193E-05 / 0,001890415 = 0,038731838 = 3,87%
Гипотеза 12: 5,54539E-05 / 0,001890415 = 0,02933426 = 2,93%

Мат. ожидание - 6,64 поколений, медиана между 5 и 6 поколениями.

Распределение совершенно совпадает с распределением вероятности для первой пары. Значит, разнесение количества мутаций по ветвям от одного предка влияния не оказывает.

Теперь давайте, выясним какое влияние оказывает на распределение вероятности времени ближайшего предка информация о промежуточных предках.
Возьмем гипотетическую тройку гаплотипов и предположим, что каждый имеет отличия от гаплотипа ближайшего общего предка на одинарные мутации в двух маркерах. Причем, все мутировавшие маркеры отличны от друг друга. Например, первый гаплотип отличается от гаплотипа ближайшего общего предка трех человек на мутации в маркерах
DYS 537 и GATA H4, второй гаплотип на мутации в маркерах DYS 442 и CDY b, третий на мутации в маркерах DYS 449 и DYS 576.
Для сравнения с этой тройкой рассмотрим ситуацию, с таким же количеством мутаций, но когда двое из трех людей имеют промежуточного общего предка, который не является предком третьего человека и отличается по значением гаплотипу как от гаплотипов тройки, так и от гаплотипа ближайшего общего предка. Пусть промежуточного предка отличает от общего предка мутация в маркере DYS 537, а потомков промежуточного предка отличают от него мутации в маркерах GATA H4 и CDY b (у каждого по одной мутации в одном из этих маркеров). Гаплотип оставшегося, третьего человека пусть отличается от гаплотипа ближайшего общего предка на мутации в маркерах DYS 449 и DYS 576.

Расчет распределения вероятности времени ближайшего общего предка для первой тройки не отличается от предыдущих расчетов.
Вероятность мутации в маркерах DYS 537 (0,00057) и GATA H4 (0,00208) и только в них за одно поколение составит (1-0,00076)*(1-0,00311)*…* 0,00208 * … * 0,00057 *…* (1-0,00087) = 9,4832E-07, для маркеров DYS 442 (0,00324) и CDY b (0,03531) - 9,49138E-05, для маркеров DYS 449 (0,00838) и DYS 576 (0,01022) - 6,96107E-05.
Проведя вычисления, пример уже которых неоднократно описывался выше получим следующее распределение априорной вероятности:
Гипотеза 1: 6,26557E-15 / 1,35074E-10 = 4,63863E-05 = 0,005%
Гипотеза 2: 2,0358E-13 / 1,35074E-10 = 0,001507179 = 0,15%
Гипотеза 3: 1,17727E-12 / 1,35074E-10 = 0,008715798 = 0,87%
Гипотеза 4: 3,35819E-12 / 1,35074E-10 = 0,024861909 = 2,49%
Гипотеза 5: 6,50369E-12 / 1,35074E-10 = 0,048149238 = 4,81%
Гипотеза 6: 9,85921E-12 / 1,35074E-10 = 0,072991411 = 7,30%
Гипотеза 7: 1,26217E-11 / 1,35074E-10 = 0,093443128 = 9,34%
Гипотеза 8: 1,42779E-11 / 1,35074E-10 = 0,105704608 = 10,57%
Гипотеза 9: 1,46952E-11 / 1,35074E-10 = 0,108793744 = 10,88%
Гипотеза 10: 1,40383E-11 / 1,35074E-10 = 0,103930672 = 10, 93%
Гипотеза 11: 1,2626E-11 / 1,35074E-10 = 0,093474837 = 9,35%
Гипотеза 12: 1,08042E-11 / 1,35074E-10 = 0,079987284 = 8,00%

Мат. ожидание - 10,33 поколений, медиана между 9 и 10 поколениями.

Расчет распределения возраста ближайшего общего предка для второй тройки несколько сложнее.

Т.к. мы предположили, что есть промежуточный предок у двух гаплотипов из трех, то значит возраст ближайшего общего предка для всей тройки не менее двух поколений (сначала должна произойти мутация в маркере DYS 537, а уж затем отдельно в маркерах GATA H4 и CDY cool.gif. Т.е. вероятность гипотезы 1 равняется нулю.

Если справедлива гипотеза 2 (общий предок всей тройки - дед), то мутация DYS 537 произошла за одно поколение (вероятность этого вычислена в предыдущих примерах, 0,000454975), мутации в маркерах GATA H4 и CDY b произошли тоже за одно поколение (0,001662771 и 0,029199465, соответственно), а мутация в маркерах DYS 449 и DYS 576 произошли за два поколения (0,000222127). Отсюда априорная вероятность гипотезы 2 равняется 0,000454975 * 0,001662771 * 0,029199465 * 0,000222127 = 4,90677E-12.

Для гипотезы 3 надо рассмотреть два варианта: промежуточный предок жил два поколения назад и что промежуточный предок жил поколение назад. В первом случае мутация в маркере DYS 537 произошла за одно поколение (0,000454975), а мутации в маркерах GATA H4 и CDY b за два (0,002652941 и 0,046587548, соответственно). Вовтором случае мутация в маркер DYS 537 произошла за два поколения (0,000725909), мутации в маркерах GATA H4 и CDY b за одно (0,001662771 и 0,029199465, соответственно). В обоих случаях мутации в маркерах DYS 449 и DYS 576 произошли за три поколения (0,000398702). Значит, общая априорная вероятность гипотезы равняется (0,000454975*0,002652941*0,046587548+0,000725909*0,001662771*0,029199465)* 0,000398702 = 3,64719E-11.

Принцип вычисления априорной гипотезы 4 схож. Надо рассмотреть случаи когда мутация в маркере DYS 537 произошла за одно, два, три поколения и мутации в маркерах GATA H4 и CDY b произошли за три, два, одно поколение соответственно. Во всех случаях мутации в маркерах DYS 449 и DYS 576 произошли за четыре поколения.
Итого: (0,000454975*0,003174561*0,055747586+0,000725909*0,002652941*0,046587548+0,00086
8638*0,001662771*0,029199465)* 0,000565446 = 1,20107E-10.

Вычислив, таким образом, априорные вероятности 300 гипотез получим распределение апостериорной вероятности:
Гипотеза 1: 0 = 0%
Гипотеза 2: 0,000436438 = 0,04%
Гипотеза 3: 0,003244032 = 0,32%
Гипотеза 4: 0,010683036 = 1,07%
Гипотеза 5: 0,023307236 = 2,33%
Гипотеза 6: 0,039506759 = 3,95%
Гипотеза 7: 0,056432219 = 5,64%
Гипотеза 8: 0,071228319 = 7,12%
Гипотеза 9: 0,081876433 = 8,19%
Гипотеза 10: 0,08748282 = 8,75%
Гипотеза 11: 0,08815340 = 8,82%
Гипотеза 12: 0,08467113 = 8,47%

Мат. ожидание - 12,6, медиана между 11 и 12.

Как мы видим, информация о промежуточных предках, о структуре выборке или как говорит Анатолий Алексеевич – «выделение молодых предков» увеличила возраст общего предка примерно на 20% (в этом конкретном случае).


Оффлайн CenturionАвтор темы

  • 100% Earth (Solar System) genofond
  • Администратор
  • *****
  • Сообщений: 9548
  • Страна: ru
  • Рейтинг +571/-2
Re: Ручное построение дерева (wertner)
« Ответ #2 : 17 Август 2009, 16:01:00 »
Продолжение...

Структура ветвей во второй тройке предыдущего примера (наличие промежуточного общего предка для двух человек из трех) хорошо подходит и для демонстрации еще одного вычисления.
Рассчитав распределение возраста ближайшего общего предка мы получили дополнительную информацию для корректировки распределения возраста промежуточного общего предка. Ведь, действительно, если бы третий человек отличался от ближайшего общего предка всей тройки на 10 мутаций, а не на 2, то это увеличило бы не только предполагаемый возраст ближайшего общего предка всей тройки, но и позволило бы предположить, что промежуточный предок двух человек все же жил раньше, несмотря на маленькое количество мутаций. И наоборот, если бы третий человек отличался от ближайшего общего предка всей тройки на 0 мутаций, а не на 2, то это уменьшило бы не только предполагаемый возраст ближайшего общего предка всей тройки, да и позволило бы предположить, что промежуточный предок двух человек все же жил позже, несмотря на наличие мутаций.
Теперь приведу пример таких вычислений. Вычислим вероятности того, что промежуточный предок жил одно поколение назад, два поколения назад и т.д., называя эти предположения соответственно гипотеза 1, гипотеза 2 и т.д.
Для вычисления вероятности гипотезы 1 выделим из вероятности каждого предположения о возрасте ближайшего общего предка всей тройки соответствующую составляющую гипотезы 1.
В случае, когда общий предок жил три поколения назад (вероятность 0,04%) промежуточный предок обязательно жил одно поколение назад. Значит, составляющая равняется 100% и в общую копилку вероятности гипотезы 1 кладем 0,04%*100% = 0,04%.
В случае, когда общий предок жил два поколения назад вероятность 0,32% складывается из составляющей гипотезы 1 и составляющей гипотезы 2. Сравнивая их вклад (0,002652941*0,046587548 и 0,000725909*0,001662771*0,029199465) получаем вклад 61% и 39%, соответственно. Тогда вклад этого случая составляет 0,32%*69%=0,20%.
В случае, когда общий предок жил четыре поколения назад (вероятность 1,07%) из всех составляющих (0,000454975*0,003174561*0,055747586+0,000725909*0,002652941*0,046587548+0,00086
8638*0,001662771*0,029199465) выделяем составляющую гипотезы 1: 0,000454975*0,003174561*0,055747586 и вычисляем ее долю: 38%. Соответственно, вклад этого случая будет 1,07%*38%=0,41%.
Вычислив подобных образом вклад для всех оставшихся 296 случаев получим итоговую вероятность гипотезы 1 равную 3,2%.

Рассчитав вероятности всех 299 гипотез получим распределение вероятности возраста промежуточного предка с учетом распределения возраста общего предка всей тройки:
Гипотеза 1: 0,032915974 = 3,29%
Гипотеза 2: 0,08695746 = 8,87%
Гипотеза 3: 0,126135397 = 12,61%
Гипотеза 4: 0,141492072 = 14,15%
Гипотеза 5: 0,136913423 = 13,69%
Гипотеза 6: 0,120135091 = 12,01%
Гипотеза 7: 0,098249674 = 9,82%
Гипотеза 8: 0,07616951 = 7,62%
Гипотеза 9: 0,056612887 = 5,66%
Гипотеза 10: 0,040661165 = 4,07%
Гипотеза 11: 0,028386025 = 2,84%
Гипотеза 12: 0,019346758 = 1,93%

Мат. ожидание - 5,84 поколений, медиана между 4 и 5 поколениями.

Сравнивая значения вероятностей или мат.ожидание распределения видим, что учет количества мутаций от общего предка до третьего человека оказало влияние на предполагаемы возраст двух более близких родственников и уменьшило предполагаемый возраст примерно на 10%. Примечательный результат не самой величиной влияния (в других случаях он может достигать и 30-40%) но самим фактом того, что структура и возраст одних ветвей оказывает влияние на возраст других.

Во всех предыдущих вычислениях у нас выходило, что для вычисления возраста предка не важны скорости мутации отдельных маркеров, а важна только скорость мутации всего гаплотипа. Что же, скорости маркеров совсем не важны? Конечно, важны. Но для другой задачи – выбора более правильного дерева.
Допустим, у нас есть 4 человека, у которых значения аллелей 65 маркеров равны, но в маркерах CDY b и DYS 537 их значения составляют классический квадрат:
У человека «А» CDY b = 37, DYS 537 = 11, у «Б» CDY b = 38, DYS 537 = 11, у «В» CDY b = 37, DYS 537 = 12, у «Г» CDY b = 38, DYS 537 = 12. Из анализа соседних гаплотипов к этой группе известно, что предковый гаплотип для этой группы совпадает с гаплотипом А. Но как построить генеалогическое древо? Ясно, что Б ближе к А, чем к В. И В ближе к А, чем к Б. И можно считать, что это разные генеалогические ветви, ведущие свое происхождение от общего предка, имеющего гаплотип А. Куда же отнести Г? Вспоминая скорости мутации маркеров CDY b и DYS 537 (0,03531 и 0,00057 по Чандлеру) приходим к выводу, что Г скорее ближе к В, т.к. скорость мутации CDY b выше и вероятней дублирование мутации этого маркера, чем маркера DYS 537. Что ж, в этом случае ответ довольно очевиден. Но как быть, когда мы решаем об отнесении к определенным ветвям не одного гаплотипа, а целой группы, где на чашах весов лежит выбор между дублирующими и возвратными мутациями нескольких маркеров и к тому же надо учесть тот факт, что в одном случае от близкого вроде предка вдруг произошло десятки мутаций. Или наоборот, что от далекого предка вдруг не произошло ни одной мутации за сотню поколений.
Пример такого вычисления я приведу для этого простого случая с упомянутыми гаплотипами А, Б, В, Г.

Мы сравниваем два дерева. В одном из них общий предок, имеющий гаплотип, совпадающий с гаплотипом А (назовем этого общего предка А') имеет потомков трех ветвей. На первой ветви находится собственно А, на второй ветви Б, на третьей ветви находится В и Г, причем их промежуточный общий предок имеет гаплотип, совпадающий с гаплотипом В – назовем его В'. Используя аналогичные обозначения опишем второе дерево: общий предок А' имеет потомков трех ветвей. На первой ветви находится собственно А, на второй ветви находится Б и Г с промежуточным предком Б', на третьей ветви находится В.
Первым шагом вычислим возраста А' и В' для первого дерева и А' и Б' для второго. Примеры подобных вычислений я уже приводил выше, поэтому сразу перейдем к результату. Для первого дерева мат. ожидание распределения возраста А' = 5,62 поколения, мат. ожидание распределения возраста В' = 2,78. Для второго дерева тоже мат. ожидание распределения возраста А' = 5,62 поколения, мат. ожидание распределения возраста Б' = 2,78, что неудивительно, т.к. деревья имеют зеркальную структуру и одинаковое количество мутаций на ветвях.
Чтобы вычислить вероятность того, что от В' произошла бы мутация в маркере CDY b и в других маркерах не произошло мутации и появился бы на свет человек с гаплотипом Г, нам пришлось бы каждую вероятность из распределения вероятности возраста В' умножить на вероятность того, что за это количество поколений произошла мутация в маркере CDY b и только в нем, а после сложить получившиеся произведения. При переходе к более далеким предкам нам пришлось рассматривать уже количество вариантов равное количество возможных поколений (я использовал в ходе рассуждений 300) в степени количества отрезков на дереве (в первом дереве примера их 5: А'А, А'Б, А'В', В'В, В'Г). При больших деревьях такие вычисления очень долги и я сознательно пошел на огрубление оценки вероятности дерева. Я предлагаю вычислять не все варианты, а только вероятность произошедших мутаций за округленную разницу между мат. ожиданием распределений вероятностей возраста предков. В рассматриваемом примере разница между мат. ожиданием возраст В' и возрастом Г (равном 0) составит 2,76. Округляем до 3. Вероятность того, что за 3 поколения произошла одна мутация в маркере CDY b и только в нем нами уже давно вычислена в начале всех примеров вычислений и составляет 0,029199465*0,797746589733803*0,797746589733803*3 = 0,055747586. Т.е. вероятность отрезка В'Г равна 0,055747586. Вероятность отрезка В'В
равна 0,797746589733803*0,797746589733803*0,797746589733803 = 0,507685628. Длину отрезка A'В' принимаем равной разнице мат. ожиданием возрастов этих предков и она равняется 5,62 – 2,76 = 2,86. Т.е. примерно 3 поколения. Соответственно, вероятность этого отрезка примерно равна вероятности того, что за три поколения произошла мутация в маркере DYS 537 и только в нем. Это 0,797746589733803 * 0,797746589733803 * 0,000454975 * 3 = 0,000868638. Вероятность отрезка А'Б равна вероятности того, что за 6 поколений (округлили 5,62) произошла мутация в маркере CDY b и только в нем, т.е. 0,056604496. Вероятность отрезка А'А равна вероятности того, что за 6 поколений не произошло ни одной мутации, т.е. 0,257744697. Перемножаем эти вероятности и получаем вероятность всего дерева: 0,055747586*0,507685628*0,000868638*0,056604496*0,257744697 = 3,58674E-07.
Для второго дерева проводим аналогичные вычисления:
Вероятность отрезка Б'Г - это вероятность мутации в маркере DYS 537 и только в нем за 3 поколения, т.е. 0,000868638 (сравните с аналогичной вероятностью отрезка В'Г первого дерева равной 0,055747586). Вероятность отрезка Б'Б – это вероятность не-мутации гаплотипа за три поколения = 0,507685628. Вероятность отрезка А'Б' – это вероятность мутации в маркере CDY b и только в нем за 3 поколения = 0,055747586. Вероятность отрезка А'В примерно равна вероятности мутации в маркере DYS 537 и только в нем за 6 поколений, т.е. 0,00088199. Вероятность отрезка А'А – это также как и в первом дереве вероятность не-мутации гаплотипа за 6 поколений - 0,257744697. Итого: 0,000868638*0,507685628*0,055747586*0,00088199*0,257744697 = 5,58872E-09
Сравнивая эти два дерева, получаем, что вероятность первого дерева в 64 раза выше, чем вероятность второго. На примере таких простых деревьев я показал способ оценки вероятностей дерева. И если в этом случае он не дал нам никакой новой информации, то в случае деревьев нескольких десятков, сотен гаплотипов такая оценка оказывает важнейшую помощь в принятии решений об отнесении гаплотипа или группы близких гаплотипов на определенную ветвь. Даже при моем небольшом опыте было множество случаев, когда без такой оценки нельзя было выбрать лучший вариант. Также было несколько случаев, когда оказывалось, что увеличение количества мутаций в дереве все же делало дерево более вероятным из-за того, что количество мутаций по разным ветвям становилось примерно равным и дерево принимало более сбалансированный, гармоничный вид.

Оффлайн Grigoriev

  • Администратор
  • *****
  • Сообщений: 3098
  • Страна: ru
  • Рейтинг +187/-3
    • Молекулярная генеалогия
  • Y-ДНК: O3a3c
  • мтДНК: J1c
Re: Ручное построение дерева (wertner)
« Ответ #3 : 17 Август 2009, 17:32:35 »
Замечательно!

Вадим, Роман, давайте оформим это в качестве статьи (я не параноик) - чтобы не терялось на форуме, и можно в случае чего написать "предквалификационная" статья, т.е. ожидающая последующей доработки и корректировки.

Оффлайн CenturionАвтор темы

  • 100% Earth (Solar System) genofond
  • Администратор
  • *****
  • Сообщений: 9548
  • Страна: ru
  • Рейтинг +571/-2
Re: Ручное построение дерева (wertner)
« Ответ #4 : 17 Август 2009, 18:20:38 »
Да, у Вадима прямо материал для статьи. Я бы на его месте все это отредактировал с картинками и статья готова.

Оффлайн wertner

  • ...
  • Сообщений: 1332
  • Страна: ru
  • Рейтинг +321/-0
    • YFull
  • Y-ДНК: E-V13->E-S2972->E-Z16661
  • мтДНК: U4a (xU4a3)
Re: Ручное построение дерева (wertner)
« Ответ #5 : 18 Август 2009, 01:28:42 »
Centurion, спасибо.

Да, у Вадима прямо материал для статьи. Я бы на его месте все это отредактировал с картинками и статья готова.
Черновик уже два месяца лежит. Уже было с переводчиком стал договариваться. Недельку бы еще над черновиком посидеть, да только пока душа не лежит. Зрею. Тренируюсь на работе: уже 20 небольших вики-статей на внутрикорпоративном портале написал  ;D

Оффлайн CenturionАвтор темы

  • 100% Earth (Solar System) genofond
  • Администратор
  • *****
  • Сообщений: 9548
  • Страна: ru
  • Рейтинг +571/-2
Re: Ручное построение дерева (wertner)
« Ответ #6 : 18 Август 2009, 01:32:01 »
Цитировать
уже 20 небольших вики-статей на внутрикорпоративном портале написал
По нашей теме?

Оффлайн wertner

  • ...
  • Сообщений: 1332
  • Страна: ru
  • Рейтинг +321/-0
    • YFull
  • Y-ДНК: E-V13->E-S2972->E-Z16661
  • мтДНК: U4a (xU4a3)
Re: Ручное построение дерева (wertner)
« Ответ #7 : 18 Август 2009, 02:16:47 »
По нашей теме?
Не, по профессиональной :)

Оффлайн I2a1a

  • ...
  • Сообщений: 10364
  • Страна: ee
  • Рейтинг +761/-8
Re: Ручное построение дерева (wertner)
« Ответ #8 : 23 Август 2009, 00:08:37 »
Вадим, не могли бы Вы подсказать какую-нибудь относительно простую и надежную прогу для расчетов Байеза?


Оффлайн wertner

  • ...
  • Сообщений: 1332
  • Страна: ru
  • Рейтинг +321/-0
    • YFull
  • Y-ДНК: E-V13->E-S2972->E-Z16661
  • мтДНК: U4a (xU4a3)
Re: Ручное построение дерева (wertner)
« Ответ #9 : 23 Август 2009, 03:05:53 »
Вадим, не могли бы Вы подсказать какую-нибудь относительно простую и надежную прогу для расчетов Байеза?
Перед тем, как написать свою прогу пользовался Excel. Других вариантов, к сожалению, не знаю.

Оффлайн I2a1a

  • ...
  • Сообщений: 10364
  • Страна: ee
  • Рейтинг +761/-8
Re: Ручное построение дерева (wertner)
« Ответ #10 : 07 Сентябрь 2009, 00:25:00 »
MrBayes -Программа для анализа филогений с помощью байезовской инференции.  Понимает формат Nexus.
Я конвертнул поправленный Вадимом У. ych файл в формат nex (отдельное спасибо Саше Лифанову за конвертер), немножко поправил файл согласно примерам, и запустил марковскую цепочку по методу Монте-Карло. Работает нормально, но медленно - на 145 гаплотипов примерно 3500 циклов за 20 минут.

Вадим, что вы думаете по этому поводу?

http://mrbayes.csit.fsu.edu/wiki/index.php/Tutorial#When_to_Stop_the_Analysis

Оффлайн wertner

  • ...
  • Сообщений: 1332
  • Страна: ru
  • Рейтинг +321/-0
    • YFull
  • Y-ДНК: E-V13->E-S2972->E-Z16661
  • мтДНК: U4a (xU4a3)
Re: Ручное построение дерева (wertner)
« Ответ #11 : 07 Сентябрь 2009, 01:07:05 »
MrBayes -Программа для анализа филогений с помощью байезовской инференции.  Понимает формат Nexus.
Я конвертнул поправленный Вадимом У. ych файл в формат nex (отдельное спасибо Саше Лифанову за конвертер), немножко поправил файл согласно примерам, и запустил марковскую цепочку по методу Монте-Карло. Работает нормально, но медленно - на 145 гаплотипов примерно 3500 циклов за 20 минут.

Вадим, что вы думаете по этому поводу?

http://mrbayes.csit.fsu.edu/wiki/index.php/Tutorial#When_to_Stop_the_Analysis
Очень интересно! Буду разбираться.

Оффлайн I2a1a

  • ...
  • Сообщений: 10364
  • Страна: ee
  • Рейтинг +761/-8
Re: Ручное построение дерева (wertner)
« Ответ #12 : 07 Сентябрь 2009, 01:09:10 »
На выходе после "прожига" деревьев генерирует деревья в формате Newick. Good.

Цитировать
The program will output a cladogram with the posterior probabilities for each split and a phylogram with mean branch lengths. The trees will also be printed to a file that can be read by tree drawing programs such as TreeView, MacClade, and Mesquite.

Quick Start Version

There are four steps to a typical Bayesian phylogenetic analysis using MrBayes:
Read the Nexus data file
Set the evolutionary model
Run the analysis
Summarize the samples
In more detail, each of these steps is performed as described in the following paragraphs.
1. At the MrBayes > prompt, type execute primates.nex. This will bring the data into the program. The data file (primates.nex) must be in the same directory as the MrBayes program. If you run into problems, refer to section 2.2 for possible solutions. If you are running your own data file, beware that it may contain some MrBayes commands that can change the behavior of the program; delete those commands or put them in square brackets to follow this tutorial.
2. At the MrBayes > prompt, type lset nst=6 rates=invgamma. This sets the evolutionary model to the GTR model with gamma-distributed rate variation across sites and a proportion of invariable sites. If your data are not DNA or RNA, if you want to invoke a different model, or if you want to use non-default priors, refer to the rest of this manual, particularly sections 3 to 5 and the Appendix.
3.1. At the MrBayes > prompt, type mcmc ngen=10000 samplefreq=10. This will ensure that you get at least 1,000 samples from the posterior probability distribution. For larger data sets you probably want to run the analysis longer and sample less frequently (the default sampling frequency is every 100th generation). You can find the predicted remaining time to completion of the analysis in the last column printed to screen.
3.2. If the standard deviation of split frequencies is below 0.01 after 100,000 generations, stop the run by answering no when the program asks "Continue the analysis? (yes/no)". Otherwise, keep adding generations until the value falls below 0.01.
4.1. Summarize the parameter values by typing sump burnin=250 (or whatever value corresponds to 25 % of your samples). The program will output a table with summaries of the samples of the substitution model parameters, including the mean, mode, and 95 % credibility interval of each parameter. Make sure that the potential scale reduction factor (PSRF) is reasonably close to 1.0 for all parameters; if not, you need to run the analysis longer.
4.2. Summarize the trees by typing sumt burnin=250 (or whatever value corresponds to 25 % of your samples). The program will output a cladogram with the posterior probabilities for each split and a phylogram with mean branch lengths. The trees will also be printed to a file that can be read by tree drawing programs such as TreeView, MacClade, and Mesquite.
It does not have to be more complicated than this; however, as you get more proficient you will probably want to know more about what is happening behind the scenes. The rest of this section explains each of the steps in more detail and introduces you to all the implicit assumptions you are making and the machinery that MrBayes uses in order to perform your analysis.

Оффлайн I2a1a

  • ...
  • Сообщений: 10364
  • Страна: ee
  • Рейтинг +761/-8
Re: Ручное построение дерева (wertner)
« Ответ #13 : 07 Сентябрь 2009, 13:01:42 »
Результаты пилотного запуска

Оффлайн I2a1a

  • ...
  • Сообщений: 10364
  • Страна: ee
  • Рейтинг +761/-8
Re: Ручное построение дерева (wertner)
« Ответ #14 : 07 Сентябрь 2009, 13:25:09 »
И кладограмма одного из 2001 вероятностных "деревьев"  (расчитанно в MrBayes, визуализированно в Mesquite)

 

© 2007 Молекулярная Генеалогия (МолГен)

Внимание! Все сообщения отражают только мнения их авторов.
Все права на материалы принадлежат их авторам (владельцам) и сетевым изданиям, с которых они взяты.