How to correct p-value for multiple comparisons regarding taxonomic rank

[lfca...@gmail.com]

Hello,

Thank you for package which is very useful, especially the possibility to compare treatments and visualize it with multiple trees. 

If I understand well, when adjusting the Wilcoxon test, it is done taking into account all taxon. 

In fact, I would like to adjust the p-value regarding each taxononmic rank. For example, if there is only 5 differents phyla in my tree, I would like to adjust all p-value corresponding to this 5 phyla only and so on for Class, Order...

Is there a way to do it ?

Thanks you,

Laurent

[Zachary Foster]

Hi Laurent,

Thanks, I am glad you found it useful!

See the example below. I think it does what you want.

Best,

Zach

library(metacoder)

# Generate exampe data
x = parse_tax_data(hmp_otus, class_cols = "lineage", class_sep = ";",
                   class_key = c(tax_rank = "taxon_rank", tax_name = "taxon_name"),
                   class_regex = "^(.+)__(.+)$")

# Convert counts to proportions
x$data$otu_table <- calc_obs_props(x, data = "tax_data", cols = hmp_samples$sample_id)

# Get per-taxon counts
x$data$tax_table <- calc_taxon_abund(x, data = "otu_table", cols = hmp_samples$sample_id)

# Calculate difference between groups
x$data$diff_table <- compare_groups(x, data = "tax_table",
                                    cols = hmp_samples$sample_id,
                                    groups = hmp_samples$body_site)

# Add rank information to table
x$data$diff_table$rank <- taxon_ranks(x)[x$data$diff_table$taxon_id]
x$data$diff_table

# Adjust p values per rank
for (a_rank in unique(x$data$diff_table$rank)) {
   x$data$diff_table$wilcox_p_value[x$data$diff_table$rank == a_rank] <- p.adjust(x$data$diff_table$wilcox_p_value[x$data$diff_table$rank == a_rank], method = "fdr")
}

[lfca...@gmail.com]

Perfect, it does exactly what I wanted.

Thank you,

Laurent