2
$\begingroup$

This is really a programming question but ggtree might be too obscure for this question to be posted on StackOverFlow as a general R related question.

I am having trouble extrapolating this this tutorial and the %<+% operator in ggtree.

I have a Newick file and a table that has a number for every tip and node in the Newick file. I want to produce a final graphic that shows a tree where the numbers from the table are shown next to respective tips and nodes.

What I want to do is close to the tutorial, I want to take a tree object and use a data-frame where one column is node-names and the second column holds numerical information per node. I want to attach the information from the data-frame to the tree object.

The data-frame looks something like this:

Node / tips Count
A 12
B 14

So, I want to ultimately display the number 12 next to node A and number 14 next to node B while keeping the names intact. Now I have the code working up to the following point:

(I am really only just working my way through learning R's grammar of graphics and R in general, so assume I am doing something wrong.)

library(ggplot2)
library(ggtree)

x <- read.tree("tree.txt")
info <- read.csv("counts.csv",sep = '\t',header=TRUE)

p <- ggtree(x) %<+% info + xlim(-.1, 5)

p2 <- p + geom_tiplab(offset = 0.2) + geom_tippoint(aes(caption = Increase)) + theme_tree2()+geom_tiplab(align=TRUE, linesize=.5)+ xlim(0,25)

ggsave("test_2.pdf", width = 60, height = 80, units = "cm", limitsize = FALSE)

print(p2)

?geom_tippoint::aes() # is there a way to list all the option inside geom_tippoint::aes()

The above code overlays the Newick_label column with the node names instead of attaching the value in the second column.

But specifically, I am at a loss as to how I can list the options available within geom_tippoint(aes()). I can use ?geom_tippoint to access the documentation of geom_tippoint but I think my task needs to manipulate options within the aes() function and the the geom_tippoint() documentation does not list all the aes() options available. And am stuck trying different keywords for geom_tippoint(aes()).

Can that be done in ggtree or is this specific application of ggtree obscure enough that there is no way to do this and produce a legible graphic? (eg issue: the graphic getting mushy overlays). I have hundred of (not uniform) genes trees so splitting data for each layer (node and tips) might be not the best way.

The files if anyone wants to try: counts and tree

Clarification on the labels for nodes and tips, I want to label both if possible. The data contains information for both.

Is annotating tips AND nodes out of the ordinary? Tips and nodes are technically different but I have quantitative data at both levels and am looking for the option to annotate at both levels.

$\endgroup$
3
  • 2
    $\begingroup$ The answer is geom_tiplab() if I get the chance this week I will post a response on a hypothetical data set, but can't do this today. Oh wait, maybe not, there is still the tip versus node issue - those are different solutions - might be easier just to copy part of the tree file to see what 'A' and 'B' are. At present the table is called node/tips: a tip is a subset of the nodes. $\endgroup$
    – M__
    Commented Dec 27, 2022 at 17:55
  • $\begingroup$ I will look into geom_tiplab() for now, thanks! $\endgroup$
    – Sudoh
    Commented Dec 27, 2022 at 18:16
  • 1
    $\begingroup$ There's nothing 'wonky' with any procedure, it's just I don't quite know what you are looking for. Main point, geom_tiplab() only works on tips and the tutorial (and the code) you provided is just for tips. It's different code for nodes and the labels are inserted at different places. Naming internal nodes is not straight forward and I'm not sure how you are doing that. It could be you assume what I would call the tip, taxa or terminal node are what you are calling a node and my taxa label or tip label are what you call a tip. If that is the case thats easy. $\endgroup$
    – M__
    Commented Dec 27, 2022 at 23:14

1 Answer 1

1
$\begingroup$

I ended up using the following code to solve this issue:

Change the tsv file and the tree file name accordingly. Be mindful of the format the tsv file should be in (indicated in the code).

Also change the tree and the font scales as needed.

# Load required libraries
library(ape)
library(ggtree)

# Read the TSV file with species and node names and numbers
# Make sure the TSV file has headers: "name", and "number"
data <- read.table("counts.tsv", header = TRUE, sep = "\t", fill = TRUE, colClasses = "character")

# Read the Newick file
newick_file <- "species_tree.newick"
tree <- read.tree(newick_file)

# Add labels and numbers to the tree
add_labels <- function(tree, data) {
  missing_numbers <- c()
  
  for (i in 1:nrow(data)) {
    label_name <- data$name[i]
label_number <- data$number[i]
    if (label_name %in% tree$tip.label) {
  tip_index <- match(label_name, tree$tip.label)
  tree$tip.label[tip_index] <- paste(label_name, label_number, sep = ":")
} else if (label_name %in% tree$node.label) {
      node_index <- match(label_name, tree$node.label)
  tree$node.label[node_index] <- paste(label_name, label_number, sep = ":")
    } else {
      missing_numbers <- c(missing_numbers, label_name)
    }
  }
  
  return(list(tree = tree, missing_numbers = missing_numbers))
}

# Modify the tree object to include labels and numbers
result <- add_labels(tree, data)
tree <- result$tree
missing_numbers <- result$missing_numbers

# Write missing labels to a log file
if (length(missing_numbers) > 0) {
  writeLines(c("Missing numbers for the following labels:", missing_numbers), "missing_labels.log")
}

# Convert factors to characters
tree$tip.label <- as.character(tree$tip.label)
tree$node.label <- as.character(tree$node.label)

# Visualize the tree
p <- ggtree(tree, layout = "rectangular") +
  geom_tiplab(size = 3) +
  geom_nodelab(size = 3, na.rm = TRUE, nudge_x = 0.008) +
  theme_tree2() +
  theme(legend.position = "none")

# Save the tree as a PDF
ggsave("species_tree_with_numbers.pdf", p, width = 20, height = 20, dpi = 300)
```
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.