Наш метод дает несмещенную максимально правдоподобную оценку, обладающую минимумом дисперсии ошибки как для TMRCA, так и для истинного предкового галпотипа, и в классе линейных оценок я не вижу способа здесь что-либо усовершенствовать.
На каждой итерации вычисляйте возраст ветви и более старшую ветвь учитывайте с бОльшим весом, чем младшую
Ведь у старшей ветви было меньше шансов на мутацию и ее предковый должен быть ближе к общему предковому.
В настоящее время
Метод выборочных пар дает НЕСМЕЩЕННУЮ оценку TMRCA. Этот факт не только аналитически получен (допустим, вывод формул неверный), но и проверен статистическим моделированием, причем, алгоритм моделирования весьма сложный, и использует массу популяционных параметров. Количество генеалогических деревьев, которое моделировалось для проверки метода, исчислялось десятками тысяч. Поведение качества оценки было вполне адекватным задаваемым параметрам скоростей мутаций, глубине генеалогических деревьев и прочим параметрам. Производилась проверка алгоритма моделирования на всяких "сингулярных" случаях.
Ваше предложение, как я понимаю, касается процедуры вычисления предкового гаплотипа. В настоящее время такая процедура, в том виде, как она изложена в нашей статье, уже дает несмещенную оценку, и здесь ничего не улучшить. Следовательно, разговор может идти только об уменьшении дисперсии ошибки вычисления предкового гаплотипа. Но наша исходная процедура основана на вычислении выборочного среднего, которое является максимально правдоподобной оценкой с минимальной границей дисперсии. Как оценку улучшить в этом плане? Есть одно соображение, а именно, уменьшить дисперсию ошибки за счет привлечения большего количества пар галоптипов, которые можно составить по данной выборке.
Смысл предложенного нами
метода заключается в том, что мы каждый раз выделяем собственный "генеалогический" путь по ветвям, который соединяет начального предка с каким-либо финальным потомком, и на этом генеалогическом пути любая мутация, независимо от того, давно она была или нет, одинаково сдвигает значение аллеля в "плюс" или в "минус". А тот факт, что влияние древних мутаций отзывается в значениях аллелей не одного, а многих финальных потомков, и чем старше мутация, тем в бОльших потомках она отзовется, проявляется в нашем
Методе тем, что дисперсия ошибки хоть и уменьшается с увеличением количества пар, по которым мы вычисляем оценку TMRCA, но значительно медленнее, чем это могло бы быть в случае, если бы все финальные потомки развивались по собственным генеалогическим путям независимо друг от друга.
Я приведу простой пример с вычислением дисперсии на основе выборки коррелированных нормально распределенных измерений. Корреляция в данном случае - это аналог того, что значения аллеля у разных финальных потомков связаны через прошлые общие мутации. Очевидно, и многие об этом забывают, что формула
D = (1/N)*SUM (Xi-M)^2 = (2/N^2) SUM SUM (Xi-Xj)^2 (1)
справедлива только для случая некоррелированных измерений, а если есть корреляция, то формула значительно усложняется. У нас аналогичная ситуация. Измерениями являются значения аллелей в каком-либо маркере у разных финальных потомков. И эти измерения (значения аллелей) коррелированы, поскольку от начального предка до каждого из финальных потомков есть общие мутации. Следовательно, и формула вычисления дисперсии должна была быть иная, а не как все привыкли в квадратичном и кросс-дисперсионном методах.
Вернемся к нашей выборке нормально распределенных коррелированных случайных величин {X(1),X(2),...,X(n)}. А вот интересно, в каких условиях можно пользоваться формулой (1) для вычисления дисперсии, если измерения на самом деле коррелированы? Очевидно, что если корреляция экспоненциально затухает и стремится к нулю только на бесконечности, то нельзя. А вот если имеется такое "расстояние" r между измерениями при котором корреляция становится нулевой, то можно. Для этого мы сформируем не все возможные пары {X(i),X(j)}, как это имеет место в правой части формулы (1), а только пары следующего типа: {X(1),X(r)}, {X(2),X(r+1)},{X(3),X(r+2)},...,{X(n-r),X(n)}. Таким образом, мы можем воспользоваться правой частью формулы (1), абсолютно наплевав на конкретный вид корреляции (по аналогии, на конкретный вид генеалогического древа), который имеется на расстояниях, меньших величины r. Как видим, нам достаточно только знания радиуса корреляции r и более ничего, чтобы определенным образом отобрать измерения для вычисления дисперсии. С оценкой матожидания (по аналогии, предкового галоптипа) все несколько проще, поскольку среднее значение оценки матожидания по коррелированной выборке не зависит от корреляции, но от нее зависит разброс этой ошибки.
У нас ситуация весьма похожая. Возьмем наугад пару i-го и j-го финальных потомков (аллели у них X(i) И X(j) соответственно). Найдем каким-то образом (не будем уточнять, каким) глубину в поколениях T1, на которой находился их БЛИЖАЙШИЙ общий предок. Понятно, что от этого ближайшего и до начального предка (глубина T0) у них имеется ровно столько общих предков, сколько поколений отделяет ближайшего общего предка от начального предка (заметим, что "начальный предок" - это в некоторой степени условность, некая точка отсчета. Имеется в виду то, что начальный предок, для субклада R1a это совсем не начальный предок гаплогруппы R1, и тем более гаплогруппы R.
Точно так же и в нашей ситуации, мы сначала выясняем, какие именно потомки (гаплотипы) сходятся к общему предку (пусть его глубину T1 мы и не знаем пока), а затем, зная структуру дерева (хотя бы до субкладов по снипам), мы отбираем только такие пары, которые не имеют взаимной корреляции, то есть, не имеют общих частей генеалогических путей, соединяющих каждого из них с общим предком. Вот только в этом случае мы имеем право использовать формулу (1). Очевидно, что для той же пары, но по отношению к более "глубокому" (глубже, чем Т1) их общему предку корреляция появится, и чем глубже мы возьмем их очередного общего предка, тем больше будет эта корреляция.
Ваша идея, уважаемый Wertner, может быть применена, если мы хотим использовать все возможные пары гаплотипов, и тогда действительно мы должны как-то учесть те самые корреляции (общие участки генеалогических путей) в виде весовых коэффициентов (а корреляции в формуле вычисления дисперсии коррелированных величин действительно используются как весовые коэффициенты). Но лично я пока не вижу, как это математически написать. Возможно, Ваша идея рекуррентных соотношений с весовыми коэффициентами позволит привлечь к вычислению TMRCA и предкового галоптипа большее количество пар исходных гаплотипов, и, тем самым, несколько понизит разброс обеих оценок.