groupClade does not function correctly on single branches

[jasonbosch]

I am trying to colour branches on a tree based on their phylogenetic classification. I want to colour all branches from the same level, e.g. Class. I figured I should be able to do that with groupClade but I realised it does not function correctly. For single branches it only labels the last two nodes regardless of the node that is given to the command. I can't imagine this is the intended behaviour.

library(tidytree)
library(ggtree)
nwk <- "(((((((Terriglobusroseus)Terriglobus)Acidobacteriaceae)Acidobacteriales)Acidobacteriae)Acidobacteriota,(((((Corynebacteriumglutamicum)Corynebacterium)Mycobacteriaceae)Mycobacteriales,(((Nocardiopsisdassonvillei)Nocardiopsis)Streptosporangiaceae)Streptosporangiales)Actinomycetia)Actinobacteriota,(((((Echinicolavietnamensis)Echinicola)Cyclobacteriaceae)Cytophagales)Bacteroidia)Bacteroidota,(((((Meiothermus_Bsilvanus)Meiothermus_B)Thermaceae)Deinococcales)Deinococci)Deinococcota,(((((Hungateiclostridiumthermocellum)Hungateiclostridium)Acetivibrionaceae)Acetivibrionales,(((Clostridium_Pperfringens)Clostridium_P)Clostridiaceae)Clostridiales)Clostridia)Firmicutes_A,(((((Desulfosporosinusacidiphilus)Desulfosporosinus)Desulfitobacteriaceae)Desulfitobacteriales)Desulfitobacteriia)Firmicutes_B)Bacteria);"
tree <- read.tree(text=nwk)
tree_grouped <- groupClade(tree, c(35,30))
ggtree(tree_grouped, aes(colour = group)) +
  geom_text(aes(label=node))

As you can see in the example code, labelling is correct for node 35 onwards (two branches) but node 30 is not labelled and labelling starts later at node 33.

[jasonbosch]

I was looking at this again and found the problem is that the way the function works is by finding tips and then grouping the tips. That only works if there are several tips on the tree, otherwise it only takes the tip and the node just before it. That's obviously not ideal because you would want a clade to start from a specific node.

library(tidytree)
library(ggtree)
#One node beyond group
nwk <- "(test1(test2(test3,test4)),(test5(test6)));"
tree <- read.tree(text=nwk)
tree_grouped <- groupClade(tree, c(6,8))
ggtree(tree_grouped, aes(colour = group)) +
  geom_text(aes(label=node))

The above example works correctly because the clade with a single tip only has two nodes.

#Three nodes beyond group
nwk3 <- "(test1(test2(test3(test7(test10)),test4(test8(test11)))),(test5(test6(test9(test12)))));"
tree3 <- read.tree(text=nwk3)
tree_grouped3 <- groupClade(tree3, c(6,12))
ggtree(tree_grouped3, aes(colour = group)) +
  geom_text(aes(label=node)) 

This second example fails to correctly identify the clade because there are more than two nodes before the tip. So even though the clade starts at node 12, the labelling only begins at node 14.

To fix this, I made a new groupClade function which I think is more intuitive and works correctly in more situations.

#Custom function
groupClade2 <- function (.data, .node, group_name = "group", ...) 
{
  Grouped_tree <- as_tibble(.data)
  Grouped_tree[,group_name] <- "root"
  N_nodes <- length(.node)
  for (GROUP in 1:N_nodes) {
    NODES <- NULL
    FINISHED <- FALSE
    New_nodes <- .node[GROUP]
    while (FINISHED == FALSE) {
      if (sum(New_nodes%in%Grouped_tree$parent) == 0 & length(NODES) == 0) {
        NODES <- New_nodes
        FINISHED <- TRUE
      }
      if (sum(New_nodes%in%Grouped_tree$parent) >= 1) {
        New_nodes <- as.vector(Grouped_tree[Grouped_tree$parent%in%New_nodes,"node"])[[1]]
        NODES <- c(NODES,New_nodes)
      } else {
        FINISHED <- TRUE
      }
    }
    if (length(names(.node)) >= 1) {
      Grouped_tree[Grouped_tree$node%in%NODES,group_name] <- as.character(names(grouping)[GROUP])
    } else {
      Grouped_tree[Grouped_tree$node%in%NODES,group_name] <- as.character(GROUP)
    }
  }
  as.treedata(Grouped_tree)
}

If we now use the custom function to create the same clades...

#Three nodes beyond group
tree_grouped3b <- groupClade2(tree3, c(6,12,2))
ggtree(tree_grouped3b, aes(colour = group)) +
  geom_text(aes(label=node))

This has several advantages and functions more intuitively.

1. The groupClade2 function takes into account all the nodes in a clade, so even long branches will be correctly labelled.
2. I found the way the branches were coloured in the original to be unintuitive because the coloured branch extended to the left of the node. That doesn't make sense to me. In this version, only the branches to the right are coloured. This is not ideal when you display node labels because they will not be coloured but, in cases without node labels, I think it will make more sense.
3. The way branches are coloured differs if you are selecting tips. The previous function could not label tips as clades. But, sometimes, you might want to label a tip and this is possible with the groupClade2 function. However, in those cases, the preceding branch is coloured as otherwise there would be no visual indication of the tip clade. This is not ideal because it is an inconsistency but I think it is acceptable for edge cases. (PS: Because the clades are labelled sequentially, if a clade (or tip) is nested inside a larger clade, the nested clade must be provided last or it will be overwritten. That is the same behaviour as the original function.)

The groupClade2 function will also work with named lists.

#Also works for named lists
nwk <- "(((((((Terriglobusroseus)Terriglobus)Acidobacteriaceae)Acidobacteriales)Acidobacteriae)Acidobacteriota,(((((Corynebacteriumglutamicum)Corynebacterium)Mycobacteriaceae)Mycobacteriales,(((Nocardiopsisdassonvillei)Nocardiopsis)Streptosporangiaceae)Streptosporangiales)Actinomycetia)Actinobacteriota,(((((Echinicolavietnamensis)Echinicola)Cyclobacteriaceae)Cytophagales)Bacteroidia)Bacteroidota,(((((Meiothermus_Bsilvanus)Meiothermus_B)Thermaceae)Deinococcales)Deinococci)Deinococcota,(((((Hungateiclostridiumthermocellum)Hungateiclostridium)Acetivibrionaceae)Acetivibrionales,(((Clostridium_Pperfringens)Clostridium_P)Clostridiaceae)Clostridiales)Clostridia)Firmicutes_A,(((((Desulfosporosinusacidiphilus)Desulfosporosinus)Desulfitobacteriaceae)Desulfitobacteriales)Desulfitobacteriia)Firmicutes_B)Bacteria);"
grouping <- list("Clostridia" = 35,"Deinococci" = 30)
tree <- read.tree(text=nwk)
tree_grouped <- groupClade2(tree, grouping,group_name = "Taxon")
ggtree(tree_grouped, aes(colour = Taxon))

You can also see above how the new branch colouring works without node labels.

Hopefully this will be useful for anyone else trying to label clades. I'm not sure whether it will cope with all situations but it seems to work for these tests of the basic functionality.