R code to reproduce all figures

[haruosuz]

I wonder if the code to reproduce Fig 6 (geographic location data) and Fig 7 (gene ontology) in the manuscript is available.

https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5340466/ Data Availability Statement A manual with documentation and examples is provided: Foster ZSL, Grunwald NJ. Metacoder user documentation [Internet]. 2016. doi:10.5281/zenodo.158228. This manual also provides the code to reproduce all figures included in this manuscript. https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5340466/figure/pcbi.1005404.g006/ Fig 6. Metacoder can be used with any type of data that can be organized hierarchically. This plot shows the results of the 2016 Democratic primary election organized by region, division, state, and county. https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5340466/figure/pcbi.1005404.g007/ Fig 7. Another alternate use example: Visualizing gene expression data in a GO hierarchy.

When I run the code example (https://grunwaldlab.github.io/metacoder_documentation/example.html), I got the following messages.

This is metacoder verison 0.3.3 (stable)

> obj$data$tax_occ <- calc_n_samples(obj, "tax_abund", groups = hmp_samples$body_site)
NOTE: Using the "groups" option without the "cols" option can yeild incorrect results if the column order is different from the group order.

> heat_tree_matrix(obj,
+                  dataset = "diff_table",

*snip*

Warning message:
Use of "dataset" is depreciated. Use "data" instead. 

[grunwald]

Dear Harao,

The code for the figures are available at out github website: https://grunwaldlab.github.io/metacoder_documentation/publication--03--greengenes.html https://grunwaldlab.github.io/metacoder_documentation/publication--03--greengenes.html

See tab under publication shows each figure somewhere there under digital PCR.

Thanks, Nik

On Jan 1, 2020, at 12:04 AM, Haruo Suzuki [email protected] wrote:

I wonder if the code to reproduce Fig 6 (geographic location data) and Fig 7 (gene ontology) in the manuscript is available.

https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5340466/ https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5340466/ Data Availability Statement A manual with documentation and examples is provided: Foster ZSL, Grunwald NJ. Metacoder user documentation [Internet]. 2016. doi:10.5281/zenodo.158228. This manual also provides the code to reproduce all figures included in this manuscript.

Fig 6. Metacoder can be used with any type of data that can be organized hierarchically. This plot shows the results of the 2016 Democratic primary election organized by region, division, state, and county.

Fig 7. Another alternate use example: Visualizing gene expression data in a GO hierarchy.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/grunwaldlab/metacoder/issues/279?email_source=notifications&email_token=AA4CKOQT3ZU3CPVWUGD4XZLQ3RE7XA5CNFSM4KB2CNE2YY3PNVWWK3TUL52HS4DFUVEXG43VMWVGG33NNVSW45C7NFSM4IDRY3RQ, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4CKOW4V4UZ3YGH3Z7UDWLQ3RE7XANCNFSM4KB2CNEQ.

[zachary-foster]

Hi @haruosuz, were you able to find the information you needed?

Thanks for pointing out the warnings generated by the examples. I have updated them.

[haruosuz]

Thank you. Yes, I found the following information.

https://grunwaldlab.github.io/metacoder_documentation/publication--08--voting.html Graphing voting geography

https://grunwaldlab.github.io/metacoder_documentation/publication--09--gene_expression.html Human gene expression

[laurahspencer]

I cannot reproduce the figures in the gene expression example. I successfully generate the final list of significant genes with their GO terms (final res object)

> head(res)
log2 fold change (MLE): dex trt vs untrt 
Wald test p-value: dex trt vs untrt 
DataFrame with 6 rows and 7 columns
                 baseMean log2FoldChange     lfcSE      stat      pvalue        padj       go_id
                <numeric>      <numeric> <numeric> <numeric>   <numeric>   <numeric> <character>
ENSG00000001167  369.3428      -0.500364  0.120880  -4.13935 3.48295e-05 3.41724e-04  GO:0000785
ENSG00000003096  377.9773      -0.940979  0.143801  -6.54363 6.00449e-11 1.51112e-09  GO:0004842
ENSG00000003137  104.5163      -0.896586  0.292593  -3.06428 2.18195e-03 1.29774e-02  GO:0001709
ENSG00000003402 2546.6142       1.190433  0.120621   9.86918 5.66250e-23 4.80618e-21  GO:0002020
ENSG00000003987   25.5043       1.004806  0.372657   2.69633 7.01081e-03 3.45151e-02  GO:0004438
ENSG00000004799  914.3790       2.644931  0.628565   4.20789 2.57766e-05 2.62790e-04  GO:0004672

but when I execute the code to run the term_class function and then generate cc_res, the result is a dataframe with 0 rows:

> (cc_res <- res[rep(1:nrow(res), sapply(cc_class, length)), ])
log2 fold change (MLE): dex trt vs untrt 
Wald test p-value: dex trt vs untrt 
DataFrame with 0 rows and 7 columns

Here is the term_class function, taken from the gene expression page on the publications tab.

term_class <- function(x, current = x, all_paths = TRUE, type = GOCCPARENTS, verbose = TRUE, 
                       valid_relationships = c("is_a")) {
  # Get immediate children of current taxon
  parents = tryCatch({
    possible_parents <- as.list(type[x[1]])[[1]]
    if (! is.null(valid_relationships)) {
      possible_parents <- possible_parents[names(possible_parents) %in% valid_relationships]
    }
    names(AnnotationDbi::Term(possible_parents))
  }, error = function(e) {
    c()
  })
  
  # only go down one path if desired
  if (! all_paths) {
    parents <- parents[1]
  }
  parents <- parents[parents != "all"]
  
  if (is.null(parents)) {
    return(c())
  } else if (length(parents) == 0) {
    return(paste0(collapse = "|", AnnotationDbi::Term(x)))
  } else {
    next_x <- lapply(parents, function(y) c(y, x))
    
    # Run this function on them to get their output
    child_output <- lapply(next_x, term_class, all_paths = all_paths, type = type)
    output <- unlist(child_output, recursive = FALSE)
    
    return(output)
  }
}