search

ssolo / ALE

Amalgamated likelihood estimation (ALE) is a probabilistic approach to exhaustively explore all reconciled gene trees that can be amalgamated as a combination of clades observed in a sample of gene trees. We implement the ALE approach in the context of a reconciliation model (cf. http://arxiv.org/abs/1211.4606 ), which allows for the duplication, transfer and loss of genes. We use ALE to efficiently approximate the sum of the joint likelihood over amalgamations and to find the reconciled gene tree that maximizes the joint likelihood among all such trees.
44 stars 15 forks source link

CCPscore returns -inf #39

Closed Gullumluvl closed 10 months ago

Gullumluvl commented 1 year ago

Hi,

I tried scoring a single consensus tree with CCPscore (ALE version from commit 6444e6e, 2022/05).

It returns -inf. This looks like an overflow, so I tried downsampling my list of trees to only 100, to no avail.

All steps (see attached files):

ALEobserve PF02881_n1.pb12.treelist.100states

CCPscore PF02881_n1.pb12.treelist.100states.ale PF02881_n1.besttree

The consensus tree is rooted and binary.

Ideally, I would find it useful if ALE returned the tree maximizing the CCPscore. I don't know how different that would be from the Maximum Clade Credibility tree. I wrote an implementation that does this, but now I would like to verify the result.

Thanks in advance.

PF02881_n1.besttree.unnamed.txt PF02881_n1.pb12.treelist.100states.unnamed.txt

ssolo commented 1 year ago

Hi,

There is a function to compute the tree that maximises the CCP score in the approx_posterior class, e.g. just add

cout << "# mpp tree from sample: "<<endl;
cout << ale->mpp_tree().first << endl;

(It's mpp for maximum posterior probability, as this trees maximises the posterior probability under the ccp approximation.)

to the end of ALEobserve.cpp and recompile and rerun and for you input I get:

((21:0.134249,22:0.0979763)1:0.0809005,(18:0.187046,(20:0.168851,(19:0.0917971,((10:0.643948,11:0.629308)1:0.497442,(((23:0.716154,25:0.359882)1:0.092694,(24:0.626521,(((6:0.0994393,7:0.197541)1:0.717951,(5:0.192851,(1:0.105243,(2:0.670927,(3:0.210146,4:0.0985567)1:0.194033)0.909091:0.365248)1:0.271608)0.909091:0.552951)1:0.240396,((17:0.47561,(16:0.34304,15:0.181203)1:0.101747)0.441176:0.141741,(12:0.234521,(13:0.336066,14:0.160136)1:0.101809)0.75:0.153546)0.479705:0.645076)1:0.113491)1:0.145007)1:0.046875,((130:0.0709174,(129:0.0528045,128:0.0908273)1:0.0569322)0.484536:0.185011,(((33:0.0843068,(32:0.0941656,(29:0.0490208,(34:0.0789278,((35:0.0470979,36:0.0816743)1:0.0791744,(30:0.05411,31:0.0758012)1:0.0510662)0.967742:0.0334181)0.57643:0.0368001)0.667897:0.0461948)0.909091:0.0435826)0.177766:0.377014,((134:0.439687,(133:0.360585,((138:0.345009,(139:0.258548,(140:0.229494,141:0.692808)1:0.202465)0.962025:0.117007)0.878274:0.0784986,((135:0.156067,(136:0.150782,137:0.220067)1:0.113608)1:0.215636,((145:0.106314,146:0.131382)1:0.428723,(142:0.587657,(143:0.0831965,144:0.0876018)1:0.0900797)0.846939:0.173562)0.977372:0.0742247)0.746341:0.0896217)0.430397:0.0814642)0.469573:0.281422)0.651118:0.981152,((44:0.232095,(45:0.314295,(42:0.172805,43:0.122021)1:0.0898215)0.675676:0.0679917)0.634267:0.0838166,(41:0.0754275,((39:0.119086,28:0.0221217)1:0.017536,(46:0.0858774,(40:0.120623,(38:0.148127,(37:0.382289,((51:0.757704,(53:0.0913335,(52:0.10284,(54:0.14005,55:0.056016)1:0.0745844)1:0.138154)0.759494:0.371968)0.705224:0.483625,(50:0.0673209,(49:0.152457,(47:0.113204,48:0.0656937)1:0.0631917)0.833333:0.0707361)0.382932:0.455877)0.517259:0.489474)1:0.106565)1:0.17795)1:0.007433)1:0.010254)1:0.075252)1:0.037596)0.75:0.046517)1:0.013267)0.727273:0.095821,((((131:0.0738362,132:0.127895)1:0.192675,(105:0.176068,148:0.417382)1:0.0988657)0.818182:0.076485,(26:0.203991,(106:0.421087,(147:0.258056,(127:0.635867,(107:0.203185,(((111:0.177161,112:0.105254)1:0.213114,(110:0.19343,(108:0.0896598,109:0.122441)1:0.808864)1:0.128686)0.78125:0.130815,((123:0.132247,124:0.0991636)1:0.0827883,(113:0.160542,((116:0.0762346,117:0.114092)1:0.0339506,((120:0.0524114,(122:0.0880276,(121:0.0589544,(114:0.13327,115:0.122946)1:0.196826)1:0.065704)0.901408:0.0690292)1:0.0483173,(119:0.0695685,118:0.0933233)1:0.0323091)0.375538:0.0293234)0.267513:0.0655198)0.293951:0.0797087)0.253264:0.110597)0.283448:0.246812)0.885448:0.494539)1:0.162008)1:0.071227)1:0.076364)1:0.138842)0.604317:0.0769084,((((57:0.0714295,58:0.0826888)1:0.151876,(56:0.0476507,(59:0.0697381,(60:0.101412,61:0.0699282)1:0.0417347)0.555556:0.0917654)0.963115:0.0687735)0.484193:0.135816,(62:0.32716,(84:0.330372,85:0.285925)1:0.157648)0.926829:0.21304)0.390297:0.287671,(27:0.248776,(125:0.150408,(126:0.534081,((102:0.262461,((8:0.225181,9:0.149354)1:0.660704,(103:0.15552,104:0.0801251)1:0.507017)1:0.22989)0.8:0.110956,((87:0.292345,86:0.519975)1:0.180993,(88:0.216108,((100:0.114125,(101:0.0937446,((97:0.178413,(98:0.162424,99:0.0967859)1:0.12732)0.82:0.221274,((90:0.0457969,91:0.326005)1:0.0585891,(92:0.209711,((95:0.0876531,96:0.0695858)1:0.0655484,(93:0.103776,94:0.0161083)1:0.0289781)0.492958:0.0597342)1:0.0519134)0.40645:0.097478)0.627338:0.121953)0.950392:0.0937556)0.45614:0.280349,(89:0.347951,((82:0.131053,83:0.130295)1:0.200198,(((75:0.256043,76:0.402787)1:0.169977,((68:0.0793565,69:0.274232)1:0.121475,(67:0.0724335,(70:0.163944,71:0.163605)1:0.119434)0.857143:0.107724)0.717608:0.310029)0.952381:0.106897,((72:0.104003,73:0.0674318)1:0.141105,((63:0.279372,(64:0.388942,65:0.136815)1:0.130845)0.788462:0.134051,((66:0.273221,(77:0.178793,(78:0.239996,79:0.153353)1:0.135173)0.953488:0.129992)0.867641:0.116634,(74:0.3218,(80:0.253832,81:0.138991)1:0.187447)1:0.097131)1:0.081931)1:0.144103)1:0.049528)0.38534:0.13317)0.353297:0.332194)1:0.210518)1:0.0735296)0.75:0.101307)1:0.0965437)1:0.103929)1:0.208235)1:0.052885)1:0.013101)1:0.01092)1:0.042893)1:0.041891)0.728072:0.265211)0.168123:0.174577)0.254569:0.245207)1:0.05434)1:0.0420149)0.422507:0.0809005)0.00156977:0.0809005;

I have not really tested this though.. so do let me know if it is not working as it should!

g

On 30 Nov 2022, at 19:56, Gullumluvl @.***> wrote:

Hi, I tried scoring a single consensus tree with CCPscore (ALE version from commit 6444e6e, 2022/05). I returns -inf. This looks like an overflow, so I tried downsampling my list of trees to only 100, to no avail. All steps (see attached files): ALEobserve PF02881_n1.pb12.treelist.100states

CCPscore PF02881_n1.pb12.treelist.100states.ale PF02881_n1.besttree

The consensus tree is rooted and binary. Ideally, I would find it useful if ALE returned the tree maximizing the CCPscore. I don't know how different that would be from the Maximum Clade Credibility tree. I wrote an implementation that does this, but now I would like to verify the result. Thanks in advance. PF02881_n1.besttree.unnamed.txt PF02881_n1.pb12.treelist.100states.unnamed.txt — Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you are subscribed to this thread.Message ID: @.***>

Gullumluvl commented 1 year ago

Amazing, thanks a lot!

I made it work on the full sample, and I found the exact same topology as with my code \o/

However, my branch lengths are different (I take the median), and my supports (I take the bipartition frequency) as well, but I am less interested in those anyway.

One remark: I tried to add these to lines into the CCPscore.cpp file instead, but then it produces a core dump upon execution.

The ALE MPP tree solved my concerns, this is great, thanks.