Creating a heat tree using the results of a glm analysis

[LeonardosMageiros]

Hi,

I was wandering if it is possible to visualize the results of my generalized linear model analysis using the heat_tree function.

Specifically I have an otu table (species vs samples) and by performing glm analysis I get something like that:

Is there a way of plotting a heat_tree with all the taxa present in my otu table and color them based on my glm reults?

If yes can you please explain how?

Thank you in advance for your time and help.

Best Leonardos

[zachary-foster]

Yea, that's not too hard I think. You would use something like x <- parse_tax_data(my_table, class_cols = 7:12) to read the table. Then you would need to convert the per-OTU data to per-taxon data, probably by looping over the result of the obs(x, 'tax_data') to get the row indexes of your input table for each taxon, combining the values somehow, and saving the result in a new per-taxon table in the taxmap object. You could then plot the per-taxon values with heat_tree. Give that a try and let me know if you have questions.

[LeonardosMageiros]

Thank you for your help.

I believe I have the first part correct as the following works fine: `.libPaths("C:/CustomR") library(metacoder) library(readr) library(taxa) library(dplyr)

setwd("C:/Monash/projects/PATCH/2_taxonomic_profiling/3_analysis/") glm_data <- read_tsv("tax_tree/glm_trees/test.txt") glm_obj <- parse_tax_data(glm_data, class_cols = 7:13) `

I even see a very basic tree by running heat_tree(glm_obj, node_label = taxon_names)

But I am not sure I understand exactly what you mean when you say I have to convert the per-OTU data to per-Taxon. Running obs(glm_obj , 'tax_data') gives me something like the following:

But I am not sure what I can do with it. Can you please give me an example of my next step?

What I would ideally like to do is color the tree using the direction of logFC and adjust the size based on the FDR value. For your ease I upload also my test input.

test.txt

Thank you in advance for your help and your great software.

Best Leonardos

[zachary-foster]

The obs command is returning the row indexes in the "tax_data" table (your input data) for each taxon in the taxonomy. For example, the taxon for "Firmicutes" would have the rows indexes for all OTUs within the "Firmicutes" taxon. These taxa are nested, so all the indexes for "Clostridia" taxon would also be in the "Firmicutes" taxon.

Since heat trees show taxon data, we need to convert the OTU data (rows in your input table) to taxon data. obs returns on value per taxon, and therefore allows us to convert per-OTU data to per-taxon data by giving the list of OTUs corresponding each taxon.

For example, you can get the mean log fold change for each taxon and add it to the taxmap object using:

library(purrr)
glm_obj$data$mean_lgc <- map_dbl(obj(glm_obj, "tax_data"), function(index) mean(glm_obj$data$tax_data$logFC[index]))

Then mean_lgc should be available as a plotting parameter. You can do the same with the "FDR" column. I like to use the p-values/FDR to set any differences that are not significant to 0 and then use color to plot the significant logFC values.

glm_obj$data$mean_lgc <- map_dbl(obj(glm_obj, "tax_data"), function(index) {
  taxon_lfc <- glm_obj$data$tax_data$logFC[index]
  taxon_pval <- glm_obj$data$tax_data$FDR[index]
  if (any(taxon_pval < 0.05)) {
    return(mean(taxon_lfc[taxon_pval < 0.05], na.rm = TRUE))
  } else {
    return(0)
  }
})
heat_tree(glm_obj,
                node_label = taxon_names,
                node_size = n_obs,
                node_color = mean_lgc)

You can see a similar example here and on other issues.

No problem, thanks!

[LeonardosMageiros]

Thank you very much for your help. This work perfectly fine now.

Just one last question mainly for the aesthetics part. Using you code I can generate a figure that looks like that:

Is there a way to just color the tips of the tree (ie only the species level) and have all the above levels colored as grey?

Thank you very much for all your help Leonardos

[zachary-foster]

Nice, I am glad it worked out. Yea, you can do something like node_color = ifelse(is_leaf, mean_lgc, "grey") in the heat_tree command.

[LeonardosMageiros]

Thx a lot for all your help in that. This last part was easier that what I thought. Can I ask one more little detail please? Now my trees are looking like that:

Is there any way to have only the leaf node colored and not the edge behind it as well?

Once again thank you for your great software.

Best Leonardos

[zachary-foster]

No problem! Try setting the edge color manually with edge_color = "grey"