unrooted gene trees treated as rooted in gCF analysis

[ryneches]

I'm trying out iqtree2's gene concordance factor analysis for the first time, and I'm a bit stumped by this.

If I pass in the boostrap trees, the gCF analysis runs as as expected, like so.

Reading tree /tmp/77-849/concat/concat.treefile ...
un-rooted tree with 111 taxa and 219 branches
Reading tree(s) file /tmp/77-849/concat/concat.ufboot ...
1000 tree(s) loaded (0 rooted and 1000 unrooted)
Computing gene concordance factor...
0.442 sec
Tree with concordance factors written to /tmp/77-849/concat/concord/concord.cf.tree
Annotated tree (best viewed in FigTree) written to /tmp/77-849/concat/concord/concord.cf.tree.nex
Tree with branch IDs written to /tmp/77-849/concat/concord/concord.cf.branch
Concordance factors per branch printed to /tmp/77-849/concat/concord/concord.cf.stat

However, if I pass in my gene trees, they are interpreted as rooted, which causes the analysis to fail.

Reading tree /tmp/77-849/concat/concat.treefile ...
un-rooted tree with 111 taxa and 219 branches
Reading tree(s) file /tmp/77-849/gene_concord/genes.treefile ...
137 tree(s) loaded (137 rooted and 0 unrooted)
Trees have different number of taxa
Computing gene concordance factor...

ERROR: Taxon not found in full tree: __root__

Which is weird to me. These trees were also made with iqtree2, and were not rooted. They were, however, manipulated with dendropy to swap the gene names for the taxon names as tip labels. However, I was careful to read them using denropy's rooting='force-unrooted' option, and to write them with suppress_rooting=True. Since that didn't work, I also tried the using ete3, hoping that this was some fresh new hell a subtlety of NEWICK formatting. Sadly, no luck with that. Or with ape.

It looks like assignBranchSupportNew() receives a parameter telling it to treat the MTreeSet instance as rooted or unrooted, but it's unclear to me how this parameter is set. Why does it think one bag-o-trees of trees is rooted, but another isn't?

Also, "best viewed in FigTree"? I appreciate the suggestion, but FigTree? Really?

[roblanf]

@ryneches I'm confused about the __root__ taxon that's been added to your trees. It's not just that the input trees are rooted - it looks like IQ-TREE is falling over because there are different numbers of taxa, and that this __root__ taxon is the cause.

Can you post an (subset of?) your data, particularly the concat.treefile, the concat.ufboot (for a comparison that works) and the genes.treefile.

Also I suspect a workaround here would be to swap names on the alignments themselves (i.e. before estimating gene trees in IQ-TREE), so that you circumvent the need for an intermediate tree processing step.

Finally, I am shocked to see the wonderful FigTree dissed. That said, I am open to suggestions of other amazing interactive tree viewers. (Particularly if we can format amazingly long and complicated branch labels and then have really flexible ways of viewing them that don't require attempting to use the excellent, obviously, various implementations of tree-viewing libraries in R and Python).

[ryneches]

The is no taxon named __root__, or root or any obvious such variation in my gene trees. I think it is induced by iqtree2 when it inferred that the gene trees are rooted. This is how I am invoking the gCF analysis :

iqtree2 -t $PREFIX/concat/concat.treefile --gcf $PREFIX/concat/concat.ufboot --prefix $PREFIX/concat/concord/concord

iqtree2 -t $PREFIX/concat/concat.treefile --gcf $PREFIX/gene_concord/genes.treefile --prefix $PREFIX/gene_concord/concord

Here are my test files (gzipped to make GitHub happy) :

concat.treefile.gz
genes.treefile.gz
concat.ufboot.gz

With any luck, this will turn out to be some embarrassing oversight on my part, and we can all go about our day.

I would like to avoid renaming the alignments. There are a lot of paralogous genes in this analysis, and rebuilding the trees with taxa represented by every combination of paralogs would be absolutely impossible.

Anyway, I was just joking around about FigTree. We have students in the lab who are younger than the codebase! When folks ask for an approachable tool for interactive tree exploration, I think the recommendation these days is usually iTOL. All of these tools have their downsides, but it's hard to beat "browse to this website" for ease of installation.

[roblanf]

OK, I have an explanation and a solution for you. I also wonder if @bqminh and I should think about force-unrooting gene trees for gcf calculation.

1. Your gene trees are still rooted (I think...)

The newick utilities manual says that newick assumes all trees are rooted, unless there's a trifurcation at the root. Your gene trees (at least, the one I looked at in the wonderful FigTree) have bifurcations at the root, so are strictly rooted. I suspect (but don't know) that this is why IQ-TREE thinks they are rooted.

2. Newick Utilities unroots them to IQ-TREE's satisfaction

If you don't already know about newick utilities, it's incredible:

• conda recipe: https://bioconda.github.io/recipes/newick_utils/README.html
• docs: https://gensoft.pasteur.fr/docs/newick-utils/1.6/nwutils_tutorial.pdf

Here's how it solves your issue.

I can reproduce your error like this (IQ-TREE v. 2.3.4)

iqtree2 -t concat.treefile --gcf genes.treefile --prefix genes

giving:

Reading tree concat.treefile ...
un-rooted tree with 111 taxa and 219 branches
Reading tree(s) file genes.treefile ...
137 tree(s) loaded (137 rooted and 0 unrooted)
Trees have different number of taxa
Computing gene concordance factor...
ERROR: Taxon not found in full tree: __root__

But if I first unroot with nw_utils, it's hunky dory (a phrase at least as old as the FigTree codebase)

nw_reroot -d genes.treefile > genes_unrooted.treefile
iqtree2 -t concat.treefile --gcf genes_unrooted.treefile --prefix genes_unrooted

Runs fine, and produces output like genes_unrooted.cf.stat e.g. first few lines:

ID      gCF     gCF_N   gDF1    gDF1_N  gDF2    gDF2_N  gDFP    gDFP_N  gN      Label   Length
112     30.77   4       0       0       0       0       69.23   9       13      100     0.0771731
113     0       0       0       0       0       0       100     16      16      100     0.0638753
114     0       0       0       0       7.14    2       92.86   26      28      95      2.1619e-06
115     0       0       0       0       0       0       100     2       2       100     0.0264323
116     0       0       0       0       0       0       100     20      20      99      2.5404e-06
117     0       0       0       0       0       0       100     24      24      69      2.1964e-06
118     3.92    2       0       0       0       0       96.08   49      51      96      0.00450702
119     0       0       0       0       0       0       100     51      51      59      2.9075e-06
120     1.37    1       0       0       0       0       98.63   72      73      82      0.00337056

A comment

Since in your case the gCFs are very low because most genes end up in gDF_P (i.e. one of the four groups required to define a branch is paraphyletic in your trees), the qCF might be more informative. For example, if most of the reason that they are ending up in gDF_P is noise, then qCF will avoid that because it uses quartets, meaning that the four groups are assumed to be monophyletic by default.

More info on that in this preprint:

https://ecoevorxiv.org/repository/view/6484/

And an accompanying almost-finished how-to-get qCFs and map them along with gCFs and other things (you need ASTRAL for qCFs) here:

https://github.com/iqtree/iqtree2/wiki/Estimating-gene,-site,-and-quartet-concordance-vectors

Hope some of that helps.

Also yes - +10 for iTOL. I do need to get good at using it, it looks amazing.

[bqminh]

Yes, I can confirm @roblanf is right. When there is a bifurcation at the root in the newick string, IQ-TREE will classify the tree as rooted. If you have trifurcation or more, IQ-TREE will treat it as unrooted tree. I believe that's also a convention in other software.

Re message about __root__: This is a "virtual taxon" that IQ-TREE internally added to a rooted tree, so that it internally becomes an unrooted tree and the whole machinery works without having to change other codes.

[ryneches]

@roblanf, you are a scholar and a gentleman. Thank you!

You are absolutely right about these gene trees, and qCF will indeed be my next stop. I have several thousand of these consensus trees to analyze, so my intent is is to use the gCF analysis to get a sense of how much of a mess the gene trees are.

I promise to stop making fun of FigTree. I named my own phylogenetics package after a meme from 2013, so stones and glass houses, etc etc.

@bqminh I thought that was where __root__ must have come from! Of course it makes sense that a trifurcation at the root would imply an unrooted tree, but I had never thought of a bifrication as the definition of a rooted tree. In my case, I think what happened is that I've lost the trifurcations in the process of trimming paralogs. There are also many packages that enforce strictly bifurcating trees.

May I humbly request that an option to override the rooting inference in some future release?

[roblanf]

Added an enhancement issue here: #251

[ryneches]

One minor issue with using nw_reroot -d, though is that some of the trees are balanced. When that happens, it can't decide which child of the root node to pick as the outgroup, and gives up. It's not an insurmountable problem, but it's not quite a trivial workaround.

[bqminh]

Well, there is a simpler solution: you can force all trees to be rooted via -rooted option. Let me know if that works.