Программа TNT.

[Valery]

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

[Alexander]

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

Да, при всех одинаковых (совпаденцах)  получился трезубец. ну хоть так! немного отлегло.

Но при трёх разных - нет. Вот 6-маркеров+3 признака снипа:

xread
9 4
Base     000000000
1     111111001
2     222222010
3     333333100
;
ccode +[0.9/1 0 /1 1 /1 2 /1 3 /1 4 /1 5 /99 6/99 7 /99 8 ;

Явно тут что-то не так?

[Овод]

Явно тут что-то не так?

Да. У Вас в трёх последних маркерах (снипах) остались попарные совпадения. Вот программа и предлагает варианты с объединением двух гаплотипов в одну ветку. То есть возможны три дерева.
 
Уберите эти три снипа - они здесь не нужны. И всё будет ОК.

[Alexander]

Да. У Вас в трёх последних маркерах (снипах) остались попарные совпадения. Вот программа и предлагает варианты с объединением двух гаплотипов в одну ветку. То есть возможны три дерева.
Уберите эти три снипа - они здесь не нужны. И всё будет ОК.

Сделел так - теперь 9 равновесных маркеров:
xread
9 4
Base     000000000
1     111111111
2     222222222
3     333333333
;
ccode +[0.9/1 0 /1 1 /1 2 /1 3 /1 4 /1 5 /1 6/1 7 /1 8 ;

Не троится.... Опять не догнал...

[I2a1a]

Не троится.... Опять не догнал...

У Вас есть право на звонок другу. Советую позвонить пользователю Valery, ибо он сейчас на проводе.

[Овод]

Не троится.... Опять не догнал...

Не удивительно. У Вас маркеры аддитивны - от предка прежде всего выделяется первый гаплотип - одна мутация во всех маркерах, затем - второй (ещё одна мутация)  и только затем - третьий.
 
Для того чтобы такие гаплотипы разошлись от корня, аллели должны быть неаддитивны.
 
Уберите объявление аддитивности.

[Valery]

Советую позвонить пользователю Valery, ибо он сейчас на проводе.

Ой не шутите так Вадим

 - Мадам, мадам, скорее! Ваш муж на проводе!
Прибегает - муж повесился (C)

[Alexander]

Для того чтобы такие гаплотипы разошлись от корня, аллели должны быть неаддитивны. 
Уберите объявление аддитивности.

Убрал, трезубец получил. Спасибо. Век живи, век учись - дураком помрёшь  .
Но "счастья" всё равно нет.
Как быть с разделением веток?
Это значит - нужно проследить, чтобы в трёх выборках (от разных веток) во всех гаплотипах не было попарно совпадающих значений аллелей - во всех маркерах (т.е. не было общих мутаций и ненулевых длин), и не было аддитивных маркеров?
Наверняка какие-то маркеры (их ведь 67) да совпадут...
И ещё про снипы.
если я задам три снипа, как 000,001, 010, то опять же появятся "общие мутации" в виде нулей и ветки не разделятся. Может быть в этом случае задать снипы, как 000, 222 и 555?

[Овод]

Отчего же, аддитивность - помеха только в Вашем конкретном случае. Совпадения маркеров - тоже.
 
У Вас же нет веток, а есть только три гаплотипа. Снипы же служат для объединения гаплотипов в ветви. А у Вас просто нечего объединять.
 
Словом, всё это - игрушки. Нужно работать с конкретными гаплотипами, а не кунстштюками, тогда и понимание реальных проблем появится. А не чисто умозрительных.

[Alexander]

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

Тем не менее, спасибо за терпение. И согласен - нужно соображать самому.

[Alexander]

Нижеизложенный поток мыслей – попытка осмысления принципа расчёта деревьев.
С одной стороны, написанное может помочь приступающим к изучению программы ТНТ.
С другой стороны, может быть найдутся знающие люди, которые смогут помочь объяснить неясные места…

Длина ветвей и дерева в ТНТ

Термины (из Вики):
Дерево состоит из вершин и рёбер.
Вершины филогенетического дерева делятся на три класса: листья, узлы и (максимум один) корень. Листья — это конечные вершины, то есть те, в которые входят ровно по одному ребру. Корень представляет общего предка всех рассматриваемых объектов. Ребра филогенетического дерева принято называть «ветвями».
Размер ствола соответствует количеству базовых гаплотипов. Длина ветвей – общему количеству мутаций на них. Длина (высота?) дерева – сумме длин ветвей.

Более сложное понятие – стоимость ветвей и дерева. В простых случаях стоимость равна длине, но в непростых – стоимость (и длина) как-то рассчитываются по загадочным формулам.
По сведениям из данной темы стоимость может быть равна взвешенному или невзвешенному количеству мутаций. Что это такое – нужно где-то почитать в литературе. В опциях ТНТ нашёл только опцию Branch lengths – в некоторых случаях там действительно количество мутаций, в других – что-то иное.

Чтобы яснее и нагляднее представить понятие о длине ветви и дерева в ТНТ, привожу несколько простых примеров с использованием искусственных 6-маркерных гаплотипов и опции ТНТ - Optimise-Branch lengths.

Пример  1

xread
6 4
Base     355555
1     123415
2     366666
3     955551     
;
ccode +0.5;  - аддитивность включена

A.tnt Matrix (6x4, 16 states).
1 trees retained.
 
Tree 0. min. branch lengths
  ┌─0 Base
├┤  ┌─10 3
  └0┤  ┌─5 2
      └0┴─12 1

Сравнив исходный файл и полученное аддитивное дерево с указанием длины ветвей, можно увидеть, что длина ветвей (тут - стоимость) равна суммарному числу мутаций от базового.
В первой ветке 2+3+2+1+4+0=12
Во второй ветке 0+1+1+1+1+1=5
В третьей ветке 6+0+0+0+0+4=10
А длина (стоимость) всего дерева в данном случае равна сумме длин ветвей=12+5+10=27 (мутаций)

Best score: 27 (Tree lengths = 27)

Пример 2 (теперь для этой же выборки отключим аддитивность)

xread
6 4
Base     355555
1     123415
2     366666
3     955551     
;                  - аддитивность отключена

1 trees retained
 
Tree 0. min. branch lengths
  ┌─0 Base
├┤  ┌─2 3
  └0┤  ┌─1 2
      └0┴─1 1

Best score (= Tree length): 12

Тут картина иная. Длина ветвей оценивается не через разности, а через число совпадений-несовпадений маркеров (по разницам):
В первой ветке 1 совпадение
Во второй ветке 0 совпадений

В третьей ветке почему-то стоит 2 – это как раз несовпадения (а совпадений – 4), но так уж, наверное, нужно (не знаю, почему).

Зато стоимость всего дерева оценивается по числу несовпадений = 5 (в первой ветке)+ 5 (во второй ветке)+2 (в третьей ветке)=12 разниц.
---------------------------------------------
Смотрим теперь, как соотношение маркеров влияет на коллапс ветвей (коллапс или «схлопывание» происходит при нулевой длине ребра), - так для примера №1 имеем коллапс-вид:

Tree 0. min. branch lengths
     ┌──0 Base
├──┤     ┌──10 3
     └──0┼──5 2
            └──12 1

Пример 3. Выборка посчти та же, только в первом гаплотипе последний маркер изменён – было 5 (совпадение с базовым), стало 4 (меньше базы)
   

xread
6 4
Base     355555
1     123414
2     366666
3     955551     
;
ccode +0.5;    - аддитивность включена     

Tree 0. min. branch lengths Best score (= Tree length): 27 = (5+9+12)
  ┌─0 Base
├┤  ┌─5 2
  └0┤ ┌─9 3
     └1┴─12 1

     Видим, что появилась 1 на ребре первого+третий гаплотипа. По взгляде на выборку понятно, что в последнем маркере исчезло како-либо совпадение с базой, и 2 разницы в минус (от базы) перевешивают одну разницу в +.
     Третий гаплотип изменён не был, поэтому общее число мутаций для него осталось 10, но 1 вылезла в суммарное ребро, поэтому на ребре от 3-го осталось 9 мутаций.
     В первом гаплотипе добавилась 1 мутация, стало в общем 13, на общем ребре 1, остаётся 12 мутаций.
     2-й гаплотип без изменений.
     Таким образом, видно, что если в выборке появляется локус, значение которого ни в одном гаплотипе не совпадает с базовым, либо по отношению к базе в плюс и в минус они неравновесны, то их разница выходит на общее ребро. Коллапса в этом случае никогда не будет – новая ветка пошла в рост с неким интервалом.

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


xread
6 4
Base     355555
1     123419
2     366666
3     955551     
;
ccode +0.5;

Tree 0. min. branch lengths
  ┌─0 Base
├┤  ┌─10 3
  └0┤ ┌─4 2
     └1┴─15 1

Best score 30
Tree lengths
       ┌──────┐
       │        +0│
┌──────┼──────┤
│     0│            30│
└──────┴──────┘

Пример 5 а) и в)– два аналогичных случая при отключённой аддитивности:
а)

xread
6 4
Base     355555
1     123414
2     366666
3     955551     
;

Tree 0. min. branch lengths №1 к базе – 1, №2 – 0, №3 – 1
                      №1 к №2 – 1, к №3 -1 (уже стало не понятно, что за длины?)
Best score (= Tree length): 13 = (6 несовпадений+5+2) – это как бы понятно
  ┌─0 Base
├┤  ┌─1 3
  └0┤ ┌─0 2
     └0┴─1 1
В)

xread
6 4
Base     355555
1     123419
2     366666
3     955551     
;

Tree 0. min. branch lengths (подсчёт аналогичен предыдущему и непонятен)
  ┌─0 Base
├┤  ┌─1 3
  └0┤ ┌─0 2
     └0┴─1 1
0.03 secs.
Best score (= Tree length): 13 (дерево аналогично предыдущему и понятно)

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

Попробуем обвесить пару ветвей.
Пример 6

xread
6 5
Base     355555
1     123419
2     366666
3     955551
4     355552   
;
ccode +0.5;

Best score: 30.  1 trees retained.
Tree 0. min. branch lengths
      ┌──0 Base
      │            ┌──0 4
├──┤     ┌──3┴──7 3
     └──0┤     ┌──4 2
            └──1┴──15 1
Тут как бы всё логично – счёт по мутациям.

Проверим, как влияет задание снипа. Дадим базе признак 000, а группирование в ветви переиначим №1-№3 (признак 001) и №2-№4 (признак 010), и признакам снипов зададим большой вес 99:

Пример 7
xread
9 5
Base     355555000
1     123419001
2     366666010
3     955551010
4     355552001   
;
ccode +0.5;
ccode /99 6.8;

Best score: 232.  1 trees retained.
 
Tree 0. min. branch lengths
      ┌──0 Base
      │              ┌──3 4
├──┤      ┌──99┴──5 2
      └──0┤       ┌──10 3
             └──99┴──16 1
Мы видим, что 99 выдавило из рёбер старые 3 и 1 в листочки, вес увеличился, но зато заданные ветки встали жёстко.

Попробуем то же на трёх ветвях:

Пример 8
xread
9 7
Base     355555000
1     123419001
2     366666010
3     955551010
4     355552001
5     244444100
6     233333100 
;
ccode +0.5;
ccode /99 6.8;

Best score: 342.  1 trees retained.
Tree 0. min. branch lengths
     ┌──0 Base
     │             ┌──2 4
├──┤    ┌──99┴──16 1
     │     │                ┌──5 6
     └──0┤     ┌──104┴──0 5
            └──0┤       ┌──9 3
                   └──99┴──5 2
И при коллапсе имеем три полноценные ветви:

Best score: 342.  1 trees retained.
Tree 0. min. branch lengths
      ┌──0 Base
      │               ┌──5 6
├──┤     ┌──104┴──0 5
      │     │       ┌──2 4
      └──0┼──99┴──16 1
             │       ┌──9 3
             └──99┴──5 2
Тут как бы всё логично.


Если бы до конца был понятен принцип вычисления длины или стоимости, то в принципе была бы решаема обратная филогении задача – синтез гаплотипов по известному дереву. Хотя и это не факт. Например, попытка синтеза гаплотипов из пары заданных ветвей привела к нерешаемой в целых числах с точки зрения числа мутаций задаче. Но программа справилась с ней легко, расставив длины ветвей чудным образом. Данное «судоку» (в смысле принципа определения длин ветвей)осталось нерешённым:

Пример 9
xread
6 5
Base     555555
1     444555
2     554445
3     533555
4     545345     
;
ccode +0.5;

Best score= Tree lengths=10.  1 trees retained.

 
Tree 0. min. branch lengths
      ┌──0 Base
      │            ┌──1 4
├──┤     ┌──2┴──0 2
     └──0┤     ┌──2 3
            └──0┴──1 1
В данном дереве непонятно даже почему его длина 10…
Когда число мутаций от базы 14…

Не вносит ясность и задание весов (то есть вес-то прибавляется, а цена ветвей так же непонятна):

Пример 10
Best score: 208.  1 trees retained.
Tree 0. min. branch lengths
     ┌──0 Base
     │                ┌──1 4
├──┤     ┌──101┴──0 2
     └──0┤       ┌──2 3
            └──99┴──1 1

Вопросов немало, но всё к одному. Вышеизложенное – что это было? Взвешенное-невзвешенное (branch lengths?)и где почитать про расчёт стоимости?

[Alexander]

Длина ветвей и деревьев ТНТ
Часть 2
Про длины и цены ветвей написано в работе J. Felsenstein «Interferring Phylogenies”, 2004. Только там, в основном, для Нетворка –Мурки (о сетях).
ТНТ же, по-видимому, программа строит бифуркационное дерево.
Справка из Вики:
«Укоренённое и неукоренённое филогетическое дерево может быть бифуркационным или небифуркационным. В бифуркационном дереве к каждому узлу подходят ровно три ветви (в случае укоренённого дерева — одна входящая ветвь и две исходящие). Таким образом бифуркационное дерево предполагает, что все эволюционные события состояли в происхождении от предкового объекта ровно двух потомков. К узлу небифуркационного дерева могут подходить четыре и более ветви».
То есть, запрограммировано из узла выходят только пара ветвей. Если случайно одно из рёбер окажется нулевой длины (без мутаций), то оно может быть сколлапсировано, и тогда из узла выйдет ветвей больше. Но такие случаи редки (в основном характерны для молодых выборок с большим количеством базовых).

Между тем, произвольная выборка не может быть 100%-но разложена на двузубцы с чёткой длиной ветвей, равной числу мутаций. Дробных же значений ТНТ не допускает. Поэтому длина ветвей определяется следующим образом. По максимуму возможного программа берёт и включает в расчёт локусы с числом мутаций, равным разности данного гаплотипа и базового. Но часть локусов (я бы назвал их «мутными» или «тёмными» аллелями) в расчёт не входят. Или входят с неким нейтральным значением.
Поясню на примере 67-маркерной выборки, состоящей из 4-х гаплотипов (включая базовый). Вид сразу ТНТ-шный (это разность значения аллеля и минимального в выборке):

xread
67 4
Modal   1122312111123123211133211121131442111111121123112121121132121112141
101     1122202141124133411032212111130322111121111134010122032232121113131
102     1122112132132031211164213121054332021111101112112121322233121121111
103     1122112132132031211164213121043632022111101122112121222233131111131
;
ccode +0.66;
 
Best score: 66.  1 trees retained
Tree 0. min. branch lengths
      ┌──0 Modal
├──┤     ┌──20 101
     └──0┤       ┌──4 103
            └──23┴──7 102

Число мутаций в 101 гаплотипе по отношению к базе 31. Однако мы видим на рис., что длина ветви равна всего 20. Это значит, что 11 мутаций, произошедших в «темных» аллелях неучтены. Если убирать по одному локусу и смотреть, какие длины и стоимости получаются, то можно найти такие локусы. Они в таблице:

Locus (number)   5   9   15   25   32   33   42   55   56   66
Base   3   1   2   1   4   4   2   1   1   4
Ivanov   2   4   3   2   3   2   1   2   2   3
Petrov   1   3   3   3   3   3   0   2   2   1
Sidorov   1   3   3   3   6   3   0   2   2   3
Branch length   no   yes=1   no   no   no   yes=1   no   no   no   no
Best Score   +1   +3   +1   +1   +2   +1   +1   +1   +1   +1
Видим, что в основном, к таковым относятся «перекошенные» маркеры. То есть отклонившиеся в одну сторону от базы, что бывает при базовом-первом, постороннем, но есть и иные случаи, по тем или иным признакам выбранные программой. Например, локус №32 – вроде бы «нормальный». Но в расчёт длины ветки не входит  И наоборот, №33 – «перекошенный», но в расчёт входит (правда, не с 2-мя мутациями, а только с 1-й).
Возможно какой-то слом происходит на середине гаплотипа (67 маркеров – 33-й как раз по середине).
Так же мутно эти локусы влияют на расчёт цены дерева, добавляя, в основном, по 1 ед. (но и больше тоже) вместо фактических мутаций (см. нижнюю строку таблицы). В остальных маркерах цена рассчитывается, как сумма длин ветвей.

Число «темных» локусов в данном случае 9/67=13%.

За счёт «темных» и длины ветвей и цена дерева всегда меньше (редко равна) истинной величины количества мутаций (к базе).

В более сложных случаях на больших выборках трудно определить «тёмных» из-за множественности деревьев. Как пример, могу привести один «тёмный» - маркер №25 (выпадающий из расчёта) маркер для выборки 90 гаплотипов+базовый:
4(базовый)3(первый)36442332533343451302233454475221011121001211221244231203221224544354435233454455544564554(девяностый).

Что за причина, вызвавшая невключение этого маркера в расчёт длин - неизвестно.

Следует также заметить, что из-за стремления к двузубцу, ТНТ принудительно вводит рёбра длиной 1 ед., которые не вывести в нулевые никакими средствами (и следовательно, не коллапсирующими), даже жёстким снипованием с большими весами.
При этом расчёт возраста предка методом выделенных пар наталкивается на сложность. Фактически, например, в выборке три-четыре или больше ветвей от корня. Но ТНТ предлагает только пару корневых ветвей. Остальные идут с единичными длинами явно принудительно.

Более того, эти единичные рёбра нарушают и расчёт длин ветвей.
Например, в 67-маркерной  выборке, из которой приведён вышестоящий пример (90 гаплотипов), есть удобный для слежения за длиной гаплотип №101 (также см. выше). И его длина на 67-маркерах равна 18 (как можно подсчитать, во всей выборке 31-18=13 «темных»). Цена всего дерева 1095.
При уменьшении числа маркеров, время от времени, скачкообразно уменьшается и длина пути, соединяющего №101 с базой. И так до момента, когда остаётся 25 маркеров в гаплотипах. При этом длина пути от №101 до базы равна 5+5 (два ребра)=10 ед.
Цена всего дерева падает до 389.
В этих 90 маркерах 3 известных ветви, которые довольно чётко разделяются (на 1 + 2 ветви через пресловутое единичное ребро).
Но при вычёркивании ещё всего одного маркера, когда остаётся 24 маркера, во первых ломается структура дерева и ветвей. Уже нет чёткого разделения. Вводится большое число единичных рёбер, а длина пути от №101 до базы увеличивается (!) снова до 18 – за счёт добавления в путь множества единичных рёбер.

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

[Alexander]

Вопрос по технике ТНТ.
Согласно рекомендациям для матрицы исходного файла берутся разности между значением аллели и её минимальным значением во всей выборке. При этом разность получается в диапазоне от 0 (совпадает с минимумом) до 9 (максимально возможный разброс значений). Если разность становится двузначной, например, 10, то файл не проходит, сбивается его структура. Как выходить из таких ситуаций?

[Овод]

Аллели более 9 должны кодироваться буквами. В руководстве написано, в какой именно последовательности.

[Alexander]

Аллели более 9 должны кодироваться буквами. В руководстве написано, в какой именно последовательности.

Спасибо, Овод. Это, ведь, я знал, когда описание переводил, и на тебе - забыл...
By default, up to 16 states are allowed by xread, using symbols 0-9 to indicate states 0-9, and A-F to indicate states 10-15.  This can be changed with the nstates command, which defines the number of states (and data type) to be subsequently read; in Windows versions, this can be changed with Format/DataType.  If the data are defined as DNA (with Format/DataType under menus, or using the command nstates dna;) the IUPAC nucleotide codes (including polymorphisms) are recognized and used throughout.  The states are internally stored as A=0, G=1, C=2, D=3, and gap=4.  If the data are instead defined as proteins (with the menus, or using nstates prot;), the IUPAC codes are recognized and used throughout.  For 32 states, non-protein data, the symbols 0-9 and A-V indicate states 0-31.