Merge pull request #3 from PNNL-Comp-Mass-Spec/devel

Changes for v0.1.0

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: MotrpacRatTraining6moWAT
 Type: Package
 Title: Helper Functions for the MotrpacRatTraining6moWATData Package
-Version: 0.0.0.9000
+Version: 0.1.0
 Authors@R: 
     c(
     person(given = "Tyler", family = "Sagendorf", 
@@ -10,7 +10,7 @@ Authors@R:
            role = c("aut", "cre")),
     person(given = "Zhenxin", family = "Hou",
            comment = c(ORCID = "0000-0002-1301-0799",
-                       "Provided WGCNA-related code."),
+                       "Provided foundation for run_WGCNA function."),
            role = "ctb"),
     person(given = "Vladislav", family = "Petyuk",
            comment = c(ORCID = "0000-0003-4076-151X",
@@ -23,26 +23,21 @@ License: MIT + file LICENSE
 Encoding: UTF-8
 Depends: 
     R (>= 3.5.0),
-    MSnbase
+    Biobase
 Imports: 
     BiocFileCache,
     circlize,
     ComplexHeatmap (>= 2.12.0),
-    cowplot,
     data.table,
     dynamicTreeCut,
     edgeR,
     fgsea,
-    ggbeeswarm,
-    ggplot2,
     graphics,
     grDevices,
     latex2exp,
     limma,
     msigdbr (>= 7.5.1),
     ontologyIndex,
-    patchwork,
-    scales,
     statmod,
     stats,
     tibble,

NAMESPACE

@@ -3,37 +3,25 @@
 export(enrichmat)
 export(fgsea2)
 export(fora2)
-export(get_ranking)
-export(ggplotMDS)
 export(label_signif)
 export(limma_full)
 export(msigdbr2)
-export(plot_ORA)
-export(plot_analyte)
-export(plot_eigenfeature)
-export(plot_upset)
-export(plot_volcano)
-export(pval_hist)
 export(range_extend)
+export(rank_genes)
 export(run_WGCNA)
-export(theme_pub)
 export(update_GO_names)
 import(ComplexHeatmap)
-import(ggplot2)
 import(limma)
 import(statmod)
+importFrom(Biobase,exprs)
+importFrom(Biobase,fData)
+importFrom(Biobase,pData)
 importFrom(BiocFileCache,BiocFileCache)
 importFrom(BiocFileCache,bfcadd)
 importFrom(BiocFileCache,bfccount)
 importFrom(BiocFileCache,bfcnew)
 importFrom(BiocFileCache,bfcquery)
 importFrom(BiocFileCache,bfcrpath)
-importFrom(MSnbase,`fData<-`)
-importFrom(MSnbase,`pData<-`)
-importFrom(MSnbase,exprs)
-importFrom(MSnbase,fData)
-importFrom(MSnbase,featureNames)
-importFrom(MSnbase,pData)
 importFrom(WGCNA,TOMsimilarity)
 importFrom(WGCNA,adjacency)
 importFrom(WGCNA,bicor)
@@ -44,20 +32,16 @@ importFrom(WGCNA,pickSoftThreshold)
 importFrom(WGCNA,plotDendroAndColors)
 importFrom(WGCNA,standardColors)
 importFrom(circlize,colorRamp2)
-importFrom(cowplot,get_legend)
-importFrom(cowplot,plot_grid)
 importFrom(data.table,`.N`)
 importFrom(data.table,`.SD`)
 importFrom(data.table,`:=`)
-importFrom(data.table,as.data.table)
 importFrom(data.table,dcast)
 importFrom(data.table,melt)
 importFrom(data.table,rbindlist)
 importFrom(data.table,setDF)
 importFrom(data.table,setDT)
 importFrom(data.table,setkeyv)
 importFrom(data.table,setnames)
-importFrom(data.table,setorder)
 importFrom(data.table,setorderv)
 importFrom(dynamicTreeCut,cutreeDynamic)
 importFrom(edgeR,DGEList)
@@ -66,8 +50,6 @@ importFrom(edgeR,cpm)
 importFrom(edgeR,filterByExpr)
 importFrom(fgsea,fgseaMultilevel)
 importFrom(fgsea,fora)
-importFrom(ggbeeswarm,geom_quasirandom)
-importFrom(ggbeeswarm,position_beeswarm)
 importFrom(grDevices,bmp)
 importFrom(grDevices,dev.off)
 importFrom(grDevices,jpeg)
@@ -83,13 +65,9 @@ importFrom(grid,grid.circle)
 importFrom(grid,grid.rect)
 importFrom(grid,unit)
 importFrom(latex2exp,TeX)
-importFrom(limma,plotMDS)
 importFrom(msigdbr,msigdbr)
 importFrom(msigdbr,msigdbr_collections)
 importFrom(ontologyIndex,get_OBO)
-importFrom(patchwork,wrap_plots)
-importFrom(scales,breaks_extended)
-importFrom(scales,extended_breaks)
 importFrom(stats,as.dist)
 importFrom(stats,cor)
 importFrom(stats,hclust)

R/enrichmat.R

@@ -250,7 +250,7 @@ layer_fun <- function(j, i, x, y, w, h, f)
                   apply(rmat, MARGIN = margin, max, na.rm = TRUE),
                   FUN = "/")
   } else {
-    rmat <- rmat/max(rmat, na.rm = TRUE)
+    rmat <- rmat / max(rmat, na.rm = TRUE)
   }
 
   grid.circle(
@@ -271,11 +271,14 @@ heatmap_color_fun <- function(NES_mat,
   NES_mat <- NES_mat[!is.na(NES_mat)]
   extended_range <- range_extend(NES_mat, nearest = 0.1)
 
-  if (all(sign(NES_mat) %in% c(0, +1))) {
-    breaks <- c(0, 1, extended_range[2])
+  # if (all(sign(NES_mat) %in% c(0, +1))) {
+  if (all(NES_mat >= -1)) {
+    breaks <- c(-1, 1, extended_range[2])
     colors <- c("white", "white", colors[2])
-  } else if (all(sign(NES_mat) %in% c(0, -1))) {
-    breaks <- c(extended_range[1], -1, 0)
+  }
+  # } else if (all(sign(NES_mat) %in% c(0, -1))) {
+  else if (all(NES_mat <= +1)) {
+    breaks <- c(extended_range[1], -1, 1)
     colors <- c(colors[1], "white", "white")
   } else {
     breaks <- c(extended_range[1], -1, 1, extended_range[2])

R/fgsea2.R

@@ -4,14 +4,14 @@
 #'   output of \code{\link{msigdbr2}}.
 #'
 #' @param pathways output of \code{\link{msigdbr2}}.
-#' @param stats output of \code{\link{get_ranking}}.
-#' @param gene_column character string; the name of a column in \code{pathways}
+#' @param stats output of \code{\link{rank_genes}}.
+#' @param gene_column character; the name of a column in \code{pathways}
 #'   containing the genes in each pathway.
-#' @param adjust.method character string; the p-value correction method. Can be
+#' @param adjust.method character; the p-value correction method. Can be
 #'   abbreviated. See \code{\link[stats]{p.adjust.methods}}.
 #' @param adjust.globally logical; should p-values from all contrasts be
-#'   adjusted together using \code{adjust.method}? Set to \code{FALSE} if the contrasts
-#'   being tested are not closely related.
+#'   adjusted together using \code{adjust.method}? Set to \code{FALSE} if the
+#'   contrasts being tested are not closely related.
 #' @param seed numeric or \code{NULL}; passed to \code{set.seed}.
 #' @param ... additional arguments passed to
 #'   \code{\link[fgsea]{fgseaMultilevel}}.
@@ -60,15 +60,17 @@ fgsea2 <- function(pathways,
   res <- rbindlist(res)
   setDT(res)
 
-
   pathways <- pathways[, list(gs_subcat, gs_exact_source, gs_description)]
 
   res <- merge(res, pathways, sort = FALSE,
                by.x = "pathway", by.y = "gs_exact_source")
 
   # p-value adjustment
   by <- "gs_subcat"
-  if (!adjust.globally) { by <- c(by, "contrast") }
+  if (!adjust.globally) {
+    by <- c(by, "contrast")
+  }
+
   res[, padj := p.adjust(pval, method = adjust.method), by = by]
 
   return(res)

R/fora2.R

@@ -3,8 +3,10 @@
 #' @description Custom wrapper for \code{\link[fgsea]{fora}}.
 #'
 #' @param pathways output of \code{link{msigdbr2}}.
-#' @param genes named list of genes.
-#' @param universe character vector of background genes.
+#' @param genes named list of genes. These are usually clusters of related
+#'   genes.
+#' @param universe character vector of background genes. Usually, this is all
+#'   unique elements of `genes`.
 #' @param gene_column character; the name of the column in pathways containing
 #'   the elements of each gene set.
 #' @param minSize integer; minimum size of gene sets allowed for testing.
@@ -16,6 +18,18 @@
 #'   adjusted together using \code{adjust.method}? Set to \code{FALSE} if the
 #'   contrasts being tested are not closely related.
 #'
+#' @details If \code{adjust.method = "scale"}, the function will calculate the
+#'   maximum overlap for each combination of \code{names(genes)} and (if
+#'   \code{adjust.globally = TRUE}), the entries of the "gs_subcat" column of
+#'   \code{pathways}. For each gene set, it will then calculate the ratio of the
+#'   size of that set's overlap to one plus the maximum overlap. The addition of
+#'   1 in the denominator is to penalize small maximum overlaps. The
+#'   \code{log10(pval)} is multiplied by this overlap ratio and then
+#'   back-transformed to obtain p-values adjusted by how well they describe each
+#'   group defined by \code{genes}.
+#'
+#' @seealso \code{\link[fgsea]{fora}}
+#'
 #' @importFrom stats p.adjust p.adjust.methods
 #' @importFrom data.table setDT rbindlist setorderv `:=`
 #' @importFrom tibble deframe
@@ -30,7 +44,7 @@ fora2 <- function(pathways,
                   universe,
                   gene_column = "entrez_gene",
                   minSize = 1,
-                  maxSize = Inf,
+                  maxSize = length(universe) - 1,
                   adjust.method = "scale",
                   adjust.globally = TRUE)
 {
@@ -63,13 +77,15 @@ fora2 <- function(pathways,
 
   # p-value adjustment
   by <- "module"
-  if (!adjust.globally) { by <- c(by, "gs_subcat") }
+  if (!adjust.globally) {
+    by <- c(by, "gs_subcat")
+  }
 
   # Transform p-values to account for overlap ratio
   if (adjust.method == "scale") {
     res[, maxOverlap := max(overlap), by = by]
-    res[, overlapRatio := overlap / maxOverlap]
-    res[, padj := 10^(log10(pval) * overlapRatio)]
+    res[, overlapRatio := overlap / (maxOverlap + 1)]
+    res[, padj := 10 ^ (log10(pval) * overlapRatio)]
   } else {
     res[, padj := p.adjust(pval, method = adjust.method), by = by]
   }

R/globals.R

@@ -42,7 +42,8 @@ utils::globalVariables(
     "moduleNum",
     "ME",
     "plotlist",
-    "x"
+    "x",
+    "rank_metric"
   )
 )