TreeNode printing and colorizing the branches

52 Views Asked by At

I am trying to generate an RShiny app which will create an interactive tree for phylogenetic and this should should print a hierarchical gene tree.

We have the following Newick File, (".nw").

structure(list(edge = structure(c(29L, 29L, 30L, 31L, 32L, 33L, 
34L, 35L, 36L, 37L, 38L, 39L, 40L, 40L, 39L, 38L, 41L, 41L, 37L, 
42L, 42L, 36L, 43L, 44L, 44L, 43L, 45L, 45L, 35L, 46L, 47L, 47L, 
46L, 48L, 48L, 34L, 49L, 49L, 33L, 50L, 50L, 32L, 51L, 51L, 31L, 
52L, 52L, 30L, 53L, 54L, 54L, 53L, 29L, 1L, 30L, 31L, 32L, 33L, 
34L, 35L, 36L, 37L, 38L, 39L, 40L, 2L, 3L, 4L, 41L, 5L, 6L, 42L, 
7L, 8L, 43L, 44L, 9L, 10L, 45L, 11L, 12L, 46L, 47L, 13L, 14L, 
48L, 15L, 16L, 49L, 17L, 18L, 50L, 19L, 20L, 51L, 21L, 22L, 52L, 
23L, 24L, 53L, 54L, 25L, 26L, 27L, 28L), .Dim = c(53L, 2L)), 
    edge.length = c(0.4551618719, 0.1476818952, 0.0990632662, 
    0.0371950264, 0.0257631258, 0.0426175095, 0.0234673195, 0.1041031134, 
    0.201856087, 0.0295971431, 0.0991494995, 0.1517511292, 0.4329144553, 
    0.2889804783, 0.7634777058, 0.0559182566, 0.5295162511, 0.7517785775, 
    0.0329285227, 0.5334232241, 0.4201915908, 0.1933011716, 0.0454164475, 
    0.4629568672, 0.5638560565, 0.3516864509, 0.1291711045, 0.1186964534, 
    0.0382650443, 0.2105653192, 0.1396779487, 0.1585672366, 0.0257753046, 
    0.3523889832, 0.281605366, 0.0456293758, 0.637940256, 0.6928517091, 
    0.4233376613, 0.2326908147, 0.2639081081, 0.075627338, 0.3064141728, 
    0.5365908087, 0.0787105838, 0.7727970128, 1.8925986754, 0.1991477293, 
    0.0440221141, 0.4123077483, 0.5161154588, 0.5192014606, 1.3871365814
    ), Nnode = 26L, node.label = c("", "100/100/100", "100/100/100", 
    "100/100/100", "100/100/100", "100/100/100", "98/100/100", 
    "100/100/100", "100/100/100", "87/95/100", "100/100/100", 
    "100/100/100", "100/99/100", "87/79/56", "100/100/100", "55/63/70", 
    "100/100/100", "100/100/100", "100/100/100", "63/72/10", 
    "97/100/100", "100/100/100", "100/100/100", "99/83/100", 
    "100/100/100", "88/92/100"), tip.label = c("Patl", "Ncra", 
    "Tmel", "Ylip", "Tdef", "Spom", "Nirr", "Scom", "Umay", "Spro", 
    "Lbic", "Ccin", "Pbla", "Rory", "Mver", "Rirr", "Ccor", "Crev", 
    "Amac", "Cang", "Spun", "Gpro", "Rall", "Mdap", "Pchi", "Nvec", 
    "Cper", "Falb")), class = "phylo", order = "cladewise")

The following code is in app.R file.

## Load necessary libraries
library(shiny)
library(ggtree)
library(treeio)

# Define UI
ui <- fluidPage(
  titlePanel("Phylogenetic Tree Visualization"),
  sidebarLayout(
    sidebarPanel(
      fileInput("file", "Choose Newick file"),
      checkboxInput("colorize", "Colorize branches", value = FALSE),
      sliderInput("thickness", "Branch thickness", min = 0.1, max = 2, value = 0.5)
    ),
    mainPanel(
      plotOutput("tree_plot")
    )
  )
)

# Define server logic
server <- function(input, output) {
  # Load tree from file
  tree <- reactive({
    inFile <- input$file
    if (is.null(inFile))
      return(NULL)
    
    # Read the Newick file
    newick <- readLines(inFile$datapath)
    
    # Define a function to parse Newick format and calculate support values
    parse_newick <- function(newick_str) {
      tree <- read.tree(text = newick_str)
      return(tree)
    }
    
    tree <- parse_newick(newick)
    return(tree)
  })
  
  # Render the tree plot
  output$tree_plot <- renderPlot({
    tree_data <- tree()
    
    # Plot the tree
    p <- ggtree(tree_data) +
      geom_treescale() +
      theme_tree2() +
      theme(legend.position = "none")
    
    # Colorize branches based on gene names
    if (input$colorize) {
      # Calculate branch lengths as a proxy for branches
      tree_data$branch <- sapply(seq_along(tree_data$edge), function(i) {
        edge <- tree_data$edge[[i]]
        if (!is.na(edge[2]) && edge[2] > 0) {
          return(sum(tree_data$edge.length[1:edge[2]]))
        } else {
          return(0)
        }
      })
      
      # Assign unique colors to branches based on genes
      num_genes <- length(unique(unlist(strsplit(gsub("[\\(\\):]", "", tree_data$tip.label), ","))))
      branch_colors <- rainbow(num_genes)
      names(branch_colors) <- unique(unlist(strsplit(gsub("[\\(\\):]", "", tree_data$tip.label), ",")))
      
      p <- p + geom_tree(aes(color = as.factor(branch)), size = input$thickness) +
        scale_color_manual(values = branch_colors)
    } else {
      p <- p + geom_tree(size = input$thickness)
    }
    
    print(p)
  })
}

# Run the application
shinyApp(ui = ui, server = server)

The code can be divided into input and output components where the input is a file and being read in the app, and the output will be the graph tree, where we can change the thickness of the tree connections and colorise it.

There are two major issues which I'm facing-

  1. The code doesn't print the gene names, i.e. the graph entries, in the tree.
  2. The code colorises all the branches with the same colour and there is no difference in the colours for different branches.

Can someone please suggest different ways to make changes to the code such that the tree nodes are printed and colorising the branches is done efficiently.

PS: If someone can give ideas about what new components to be included in the RShiny interactive app, that will be helpful as well.

1

There are 1 best solutions below

0
hmmmm On BEST ANSWER

Based no the requirements: This is the code I tried implementing

library(ape)
library(dplyr)

# Load tree from file
tree <- read.tree("MCs70_Holomycota_phylobayesConsensusTree.nw")

# Calculate node depths
node_depths <- node.depth(tree)

# Get node labels and depths
node_labels <- data.frame(node = tree$node.label, depth = node_depths[match(tree$node.label, tree$node.label)]) %>%
  filter(!is.na(node))

# Render the tree plot
plot(tree, edge.color = rainbow(length(node_depths)), edge.width = 2, label.offset = 0, cex = 1.2)

# Set node styles
nodelabels(pch = 21, bg = rainbow(length(tree$tip.label))[as.numeric(tree$tip.label)], col = "black", cex = 1.2)

# Print depth of each branch
edgelabels(round(node.depth.edgelength(tree), 2), bg = "white")

# Show node entries if selected
if (TRUE) {
  cat("Node Entries:\n")
  cat(paste(node_labels$node, collapse = "\n"))
}

maybe you can try implementing it in the app. if you are unable to do so, please let me know, I'll try to help you out there.