Regarding support values derived from the TreeSearch resampling results

[HS6986]

Dear @ms609,

I'm thinking about support values derived from the resampling results of TreeSearch. To obtain the support values for an optimal tree in TreeSearch, I calculate the strict consensus of most parsimonious trees (MPTs) obtained in each jackknife iteration and apply the strict consensus trees to the optimal tree, using TreeSearch::JackLabels() or TreeSearch::PresCont(). This approach seems to be also used, for example, in TNT.

However, my concern on this approach is that strict consensus trees which are unresolved with respect to a certain split may underestimate (?) the support value for the split. For example, consider the case in which we have a total of 100 strict consensus trees derived from a jackknife resampling. Among them, 10 support a certain split and 90 are unresolved with respect to the split. In this case, JackLabels(), for example, computes the support value for the split as 10%, despite that the majority of the trees that are unresolved with respect to the split likely doesn't decisively deny the split.

As such, I've come up a new approach for this that results in probably more informative and higher support values. However, I'm not sure whether the approach or something similar to that has been proposed or implemented somewhere (I'm not familiar with the theories and circumstances in this area).

My approach is the following:

1. Take an optimal tree (which may be bifurcating or multifurcating) and a list of multiphylo of MPTs which were derived from a jackknife resampling in TreeSearch.
2. For each split, for each iteration, determine whether the iteration supports the split (all the MPTs support the split), ambiguous as to the split (some of the MPTs support the split, while others deny the split), or deny the split (all the MPTs deny the split).
3. For each split, calculate (the number of the iterations supporting the split) / (the number of the iterations supporting or denying the split), which becomes the support value for the split.

My questions are the following:

1. Is this approach justified? Or is it a good idea?
2. Has this approach, or something similar, been proposed or implemented somewhere?
3. Can better approach(es) be considered to address this issue?

If you think my approach would be valuable, I would greatly appreciate it if you could consider implementing it in TreeSearch.

Thank you very much for your time and support.

[ms609]

If by 'iteration' you mean 'sampled most parsimonious tree', then this approach is invalid for the same reason that it is invalid to compute a majority rule consensus from a set of MPTs, i.e. the frequency of a split's occurrence is not related to its support. But perhaps I misunderstood your proposal?

[ms609]

I've put together a draft implementation in #198. Let me know what you think.

[HS6986]

Dear @ms609,

Thank you very much for your draft implementation.

I've tried JackLabels() for my sets of MPTs with the following code:

allConsensus$node.label <- JackLabels(allConsensus, allJacksS, plot = FALSE, showFraction = TRUE)
write.tree(allConsensus, file = "allConsensusSato.tre")

But there seemed to be an issue with the output tree (the fraction of the decisive replicates for each split was not annotated and strange NAs were inserted at the end):

(Choanoflagellata,(Porifera,(Placozoa,((((Xianguangia,Daihua)0.96551724137931,Calathites)0.796610169491525,Dinomischus_venustus,Dinomischus_isolatus,(((((((Nematostella,(Anemonia,(Aiptasia,Metridium)0.876288659793814)0.98936170212766)0.73469387755102,Zoanthus)0.829268292682927,(Corynactis,(Montastraea,Porites,Acropora)0.581818181818182)0.99)0.597938144329897,Antipathes)0.927083333333333,Cerianthus)0.99,(((Anthomastus,(Nephthyigorgia,Leptogorgia)1)1,Cornularia)1,Virgularia)1)0.99,(((Haliclystus,Craterolophus)1,((Carybdea,Chiropsalmus,Carukia)1,((Nausithoe,Atorella)1,(Aurelia,Rhizostoma)1)1)1)1,(((Hydra,(Physalia,Craseoa,Nanomia)1)0.559139784946237,(Candelabrum,Ectopleura)0.775510204081633,Hydractinia)0.870967741935484,(Clytia,Obelia)0.858695652173913,(Aegina,Halitrephes)1)1)1)0.99,(((((Entoprocta,(Mollusca,(Annelida,((Chaetognatha,Gnathifera)1,Gastrotricha)1,Nemertea,Platyhelminthes)0.421052631578947)0.202702702702703)0.168421052631579,(Bryozoa,(Brachiopoda,Phoronida)0.745762711864407)0.623655913978495)1,((Onychophora,Arthropoda,Tardigrada)0.923076923076923,((Nematoda,(Priapulida,Kinorhyncha,Loricifera)0.790697674418605)0.816091954022989,Nematomorpha)0.867469879518072)1)1,(((Urochordata,Vertebrata)1,Cephalochordata)1,Xenacoelomorpha)0.464646464646465)0.463414634146341,(Echinodermata,(Pterobranchia,Enteropneusta)1)1)0.71)0.452631578947368,Siphusauctum)0.516129032258065,Daihuoides,Fasciculus,Xanioascus,Ctenorhabdotus_capulus,Ctenorhabdotus_campanelliformis,Thalassostaphylos,Thaumactena,Galeactena,Maotianoascus,Batofasciculus,Trigoides,Gemmactena,Beroe,Euplokamis,Pleurobrachia,Haeckelia,(((Mnemiopsis,(Cestum,Ganesha)0.245283018867925)0.245283018867925,Ocyropsis)0.75,Thalassocalyce)0.901639344262295,Ctenoplana)0.988235294117647)NA)NA;

[ms609]

I've updated the output format to be more intuitive, consistent and (via a new format argument) customizable.

[HS6986]

Dear @ms609,

Thank you for updating the code. It seems to be working fine now, though I haven't conducted thorough testing.

[ms609]

Great. Let me know when you publish something on this and I'll add a citation to the function documentation.