АвторТема: Влияние топологии дерева на сдвиг оценки ТМRCA по внутрикладовой АSD.  (Прочитано 12303 раз)

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

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

  • Главный модератор
  • *****
  • Сообщений: 1785
  • Рейтинг +390/-3
  • Omnia mea mecum porto
  • Y-ДНК: R1a-M198
  • мтДНК: U4a
Недавно, пытаясь (с подачи Вадима Веренича) "скрестить" ASD и ро-статистику в расчётах ВБОП по полученным в результате филогении деревьям, я "по-ходу" получил оценку для популяционного влияния дерева с заданной топологией на расчёт его возраста по ASD. Хотя полученная формула имеет лишь косвенное отношение к первоначальной затее, она может оказаться даже интересней в свете недавних разговоров о "популяционном смещении" в ASD-методах.
 
Итак, пусть задано произвольное дерево со следующими параметрами:
 
N - число современных гаплотипов на нём
Т - возраст дерева (число поколений от корня до вершины)
Li и Ki - cоответственно длина (в поколениях) и кратность i-того ребра дерева.
 
Здесь под "ребром" понимается участок дерева без ветвлений на нём, а под "кратностью" - число финальных потомков, к которым ведёт этот участок.
 
Для начала рассмотрим формулу для "смещенной" ASDc = SUMij(ai-aj)2/N2, переводя которую к гистограммному виду (по частоте аллелей в выборке) и рассматривая случай единичной мутации на любом i-том ребре с кратностью Кi получаем приращение дисперсии за счёт этой мутации равное 2*Ki(1-Ki/N)/N. Усредняя этот результат по всем рёбрам (с весом их длины), получаем, для средней  ASDс в результате ожидаемого числа мутаций на дереве со средней скоростью m мутаций/поколение ASDc = 2mSUMi(LiKi(1-Ki/N)/N), и, учитывая, что для любого дерева SUMi(LiKi/N) = T, имеем итоговую оценку смещения относительно истинного возраста Т:
 
Т - ASDС/2m = SUMi(LiKi2/N2)                            (1)
 
Однако, для "смещённой" дисперсии мы не можем чётко выделить именно "популяционное смещение" , так как чисто статистически эта оценка также смещена. Чтобы его всё же найти, рассмотрим формулу "несмещённой"  ASD:
 
ASDН = SUMij(ai-aj)2/(N(N-1)) = ASDС*N/(N-1)
 
и, домножая обе части уравнения (1) на N/(N-1), после нехитрых преобразований, получаем:
 
Т - ASDН/2m  = SUMi(LiKi(Ki-1))/(N(N-1))

что и является искомой формулой чисто "популяционного" смещения для внутрикладовой ASD. Интересны два факта, вытекающих из этой формулы:
 
1. Терминальные ветви дерева, кратность которых всегда равна 1 не влияют на смещение дисперсии, давая нулевой результат, а всё "зло" кроется исключительно в общих ветвях, имеющих двух или более современных потомков из выборки. Собственно, их наличие и является истинной причиной популяционного сдвига.
 
2. Смещение квадратично возрастает с увеличением кратности ветвей и линейно - с их длиной. То есть, наличие большего числа выживших потомков предка имеет гораздо больший эффект, чем его близость по возрасту к современникам, которая тоже увеличивает смещение - но в меньшей степени.
 
Надеюсь, что эта формула может помочь, если не в конкретных расчётах и прикидках, то в понимании сути "популяционного эффекта" - наверняка.
 
P.S. Eщё раз подчеркиваю, что величины ASDC и ASDH в приведённых здесь формулах следует понимать как матожидание (среднее) ASD, а не любое (выборочное) её значение.
« Последнее редактирование: 14 Декабрь 2010, 19:19:37 от Овод »

Оффлайн Nimissin

  • Сообщений: 2347
  • Рейтинг +724/-0
  • Y-ДНК: N-M178 L839+ P298+ M2019+ M2118+ M1991+ M1988+
  • мтДНК: C4b12a
Уважаемый Овод, полученная Вами формула мне очень нравится. Но пока чисто интуитивно. Возьму некоторое время для изучения.

Оффлайн shekhol

  • Сообщений: 752
  • Страна: ru
  • Рейтинг +127/-2
  • Y-ДНК: N1*(L1025+)
Спасибо, Овод!

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

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

  • Главный модератор
  • *****
  • Сообщений: 1785
  • Рейтинг +390/-3
  • Omnia mea mecum porto
  • Y-ДНК: R1a-M198
  • мтДНК: U4a
правильно ли я понял,
что, при разветвлённой топологии дерева, оценивая ВОП по выжившим потомкам мы можем получить несколько удревнённое ВОП.

Наоборот, омоложенное.

Оффлайн Nimissin

  • Сообщений: 2347
  • Рейтинг +724/-0
  • Y-ДНК: N-M178 L839+ P298+ M2019+ M2118+ M1991+ M1988+
  • мтДНК: C4b12a
Появились следующие вопросы:

1. Вызвает сомнение верность формулы
SUMi(LiKi/N) = T.
Похоже, она верна только для случая, когда все ветви растут непосредственно от общего предка. В иных случаях после суммирования должен получаться средневзвешенный по кратности ветвей Ki возраст ветвей <L>, меньший чем T.

2. Рассмотрим самый простой случай двух ветвей, равных по кратности и происходящих от самого общего предка. Тогда L1 = L2 = T, K1 = K2 = N/2. Предположим также, что N велико и вместо (N-1) можно использовать N. Подставим значения в формулу для смещенного ASDc:

Т - ASDс/2m = T (1/4 + 1/4) = T/2.

Тогда ASDс/2m = T/2.

Результат неверный, так как из начальных условий должно быть

ASDс/2m = T.

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

Оффлайн dima75

  • Сообщений: 558
  • Страна: il
  • Рейтинг +100/-1
    • Tartakovsky surname project
  • Y-ДНК: E1b1b1-PF1975 (E-Z830-A)
  • мтДНК: K1a1b1a
правильно ли я понял,
что, при разветвлённой топологии дерева, оценивая ВОП по выжившим потомкам мы можем получить несколько удревнённое ВОП.

Наоборот, омоложенное.
Звучит логично даже без формул: понятно, что без внесения поправок ВБОП будет омоложен в сторону самой большой ветви.

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

  • Главный модератор
  • *****
  • Сообщений: 1785
  • Рейтинг +390/-3
  • Omnia mea mecum porto
  • Y-ДНК: R1a-M198
  • мтДНК: U4a
1. Вызвает сомнение верность формулы
SUMi(LiKi/N) = T.
Похоже, она верна только для случая, когда все ветви растут непосредственно от общего предка. В иных случаях после суммирования должен получаться средневзвешенный по кратности ветвей Ki возраст ветвей <L>, меньший чем T.

Формула верна абсолютно для любого мыслимого дерева. Ещё раз подумайте, и если придумаете хотя бы одно дерево, для которого это соотношение нарушено, то я Вам выдам персональную нобелевскую премию в виде ящика коньяка. Или хотите равноправное пари?

2. Рассмотрим самый простой случай двух ветвей, равных по кратности и происходящих от самого общего предка. Тогда L1 = L2 = T, K1 = K2 = N/2. Предположим также, что N велико и вместо (N-1) можно использовать N. Подставим значения в формулу для смещенного ASDc:

Т - ASDс/2m = T (1/4 + 1/4) = T/2.

Тогда ASDс/2m = T/2.

Результат неверный, так как из начальных условий должно быть

ASDс/2m = T.

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

Результат верный. Для двух рёбер и гаплотипов может быть только два. И уж никак нельзя считать в этом случае, "что N велико". Вы употребили "смещённую" формулу и получили ожидаемое 50%-ное статистическое смещение, равное для двух гаплотипов T/2.
 
В таких случаях можно считать только по "несмещённой" формуле. Подставьте свои цифири в неё и убедитесь, что она даёт верный ответ, то есть ASDH /2m = Т.
 
Я вообще рекомендую народу считать по "несмещённой" формуле всегда, хоть она и не слушком удобна. И уж, во всяком случае, для N<100, чтобы ошибка не превысила 1%.

Оффлайн Yurgan

  • Кто везёт, тому везёт
  • ...
  • Сообщений: 8255
  • Страна: gt
  • Рейтинг +825/-7
  • Потомок кузнеца Ильмаринена
    • Сибирский родословец
правильно ли я понял,
что, при разветвлённой топологии дерева, оценивая ВОП по выжившим потомкам мы можем получить несколько удревнённое ВОП.

Наоборот, омоложенное.
Согласен.

Оффлайн Каржавин

  • ...
  • Сообщений: 1806
  • Рейтинг +144/-2
Итак, пусть задано произвольное дерево со следующими параметрами:
N - число современных гаплотипов на нём
Т - возраст дерева (число поколений от корня до вершины)
Li и Ki - cоответственно длина (в поколениях) и кратность i-того ребра дерева.
Здесь под "ребром" понимается участок дерева без ветвлений на нём, а под "кратностью" - число финальных потомков, к которым ведёт этот участок.
Т - ASDН/2m  = SUMi(LiKi(Ki-1))/(N(N-1))
Правильно ли я понимаю, что если все ветки начинаются с первопредка и не ветвятся, то кратность К становится равной N, суммирование вырождается в один член L=T, так как Ki(Ki-1)/(N(N-1) =1. При этом получается, что
T - ASD/2m = T

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

  • Главный модератор
  • *****
  • Сообщений: 1785
  • Рейтинг +390/-3
  • Omnia mea mecum porto
  • Y-ДНК: R1a-M198
  • мтДНК: U4a
Итак, пусть задано произвольное дерево со следующими параметрами:
N - число современных гаплотипов на нём
Т - возраст дерева (число поколений от корня до вершины)
Li и Ki - cоответственно длина (в поколениях) и кратность i-того ребра дерева.
Здесь под "ребром" понимается участок дерева без ветвлений на нём, а под "кратностью" - число финальных потомков, к которым ведёт этот участок.
Т - ASDН/2m  = SUMi(LiKi(Ki-1))/(N(N-1))
Правильно ли я понимаю, что если все ветки начинаются с первопредка и не ветвятся, то кратность К становится равной N, суммирование вырождается в один член L=T, так как Ki(Ki-1)/(N(N-1) =1. При этом получается, что
T - ASD/2m = T

Нет, не правильно. Внимательно прочтите определение кратности ребра, приведённое в цитирумом Вами участке (выделено мной жирным шрифтом), и Вы поймёте, что "кратность" - это не количество рёбер равной длины, а количество современников, к которым ведёт данное ребро.
 
В Вашем примере каждое ребро ведёт к одному потомку  и кратность каждого равна 1, а длина - количеству поколений до предка. И, следовательно, правильный ответ:
 
Т-ASDН/2m = 0

Оффлайн Каржавин

  • ...
  • Сообщений: 1806
  • Рейтинг +144/-2
Как я понимаю, Li берутся из филогенетического дерева, в котором эти Li на самом деле не соответствуют истинным длинам (короче).
Я нечто похожее вводил в оценку, причем, данные брал из дерева, которое мне известно было полностью (коли моделировал, значит, полностью известна структура дерева). Действительно, смещение оценки ЗНАЧИТЕЛЬНО уменьшалось, но увеличился разброс оценки. А ведь это самое главное - разброс. Смещение кое-как можно и калибровочными апроксимационными кривыми компенсировать.
Поскольку предлагаемая Вами формула (как и мои экзерсисы на данную тему) эмпирическая и не выведена из соответствующей статистической модели, хоть и тенденции она отслеживает в верном направлении, то увеличение разброса оценки я Вам гарантирую.
Я штуки 4-5 аналогичных формул учета филогении пробовал. Вывод однозначен: нужно, чтобы формула вытекала хоть из какой-то правдоподобной статистической модели и из соответствующих вероятностных распределений.

Оффлайн VVR

  • ...
  • Сообщений: 2462
  • Страна: ua
  • Рейтинг +616/-0
  • Y-ДНК: o.R1a1a1b1a2a1a1a1e~-YP569,YP1260+;м.R1a1a1b1a1a1a2~-L260,YP1337+
  • мтДНК: K1c1h
Недавно, пытаясь (с подачи Вадима Веренича) "скрестить" ASD и ро-статистику в расчётах ВБОП по полученным в результате филогении деревьям,
Не помню о каком сообщении Вадима речь, но если вдруг Вы перепутали В.В. с VVR, а ASD с ММП (http://forum.molgen.org/index.php/topic,1682.msg67422.html#msg67422) , то попробую объяснить, что я имел ввиду, тем более Каржавин тоже коснулся вопроса топологии.
Я имел ввиду при построении дерева применять ММП для проверки(уточнения) филогении, а затем на основе максимально правдоподобного дерева рассчитывать ро-статистику.

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

  • Главный модератор
  • *****
  • Сообщений: 1785
  • Рейтинг +390/-3
  • Omnia mea mecum porto
  • Y-ДНК: R1a-M198
  • мтДНК: U4a
Нет, я имел в виду именно Вадима. Он просил меня здесь http://forum.molgen.org/index.php/topic,744.msg47029.html#msg47029 дополнить мой скрипт по ро-статистике расчётом ASD ещё весной. Но я что-то закопался и так его и не закончил (зато получил формулу сдвига).
 
Вы же говорили о ММП.

Оффлайн VVR

  • ...
  • Сообщений: 2462
  • Страна: ua
  • Рейтинг +616/-0
  • Y-ДНК: o.R1a1a1b1a2a1a1a1e~-YP569,YP1260+;м.R1a1a1b1a1a1a2~-L260,YP1337+
  • мтДНК: K1c1h

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

  • Главный модератор
  • *****
  • Сообщений: 1785
  • Рейтинг +390/-3
  • Omnia mea mecum porto
  • Y-ДНК: R1a-M198
  • мтДНК: U4a
Я имел ввиду при построении дерева применять ММП для проверки(уточнения) филогении, а затем на основе максимально правдоподобного дерева рассчитывать ро-статистику.

Ваша идея мне понятна. Она вполне здравая. Есть несколько коммерческих программ, которые во главу угла при филогении ставят не принцип парсимонии, а формулу Байеса или ММП. Но пока они не конкурентны из-за медленной работы.
 
Многие так и поступают, как Вы предлагаете: сначала находят все оптимальные (и субоптимальные) деревья парсимонией, а потом выбирают наилучшее по ММП или Байесу. Где-то в сети был даже обзор на эту тему.

 

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

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