feat: create fetch_celldex_ref() and label significant genes in make_volcano_plot()

.Rbuildignore

@@ -17,3 +17,8 @@
 ^Meta$
 ^SCOT\.Rproj$
 ^README.html$
+^bucket$
+^status$
+^tests/testthat/bucket$
+^tests/testthat/status$
+^\.vscode$

.gitignore

@@ -9,5 +9,7 @@ inst/doc
 /doc/
 /Meta/
 /README.html
-
+bucket
+status
 /.quarto/
+.vscode

DESCRIPTION

@@ -22,6 +22,7 @@ Imports:
     DropletUtils,
     ggplot2,
     ggpubr,
+    ggrepel,
     grDevices,
     harmony,
     magrittr,
@@ -47,12 +48,14 @@ Suggests:
     R.utils,
     rmarkdown,
     roxygen2,
+    SeuratData,
     SingleCellExperiment,
     testthat (>= 3.0.0),
     usethis
 Remotes:
     chris-mcginnis-ucsf/DoubletFinder,
     rnabioco/djvdj@*release,
+    SeuratData=satijalab/seurat-data,
     SeuratWrappers=satijalab/seurat-wrappers
 biocViews:
 Config/Needs/check: rcmdcheck, lintr, goodpractice

R/add_VDJ_metadata.R

@@ -7,9 +7,6 @@
 #' @param so A Seurat single cell RNA object
 #' @param vdj_file A path to a V(D)J annotation directory from CellRanger
 #'
-#' @import djvdj
-#' @import Seurat
-#'
 #' @export
 #'
 #' @return A Seurat object with V(D)J metadata

R/cluster_metrics.R

@@ -26,7 +26,6 @@
 #'
 #' @export
 #'
-#' @import clusterSim
 cluster_metrics <- function(so,
                             cluster_list,
                             dims = 1:20,
@@ -40,7 +39,7 @@ cluster_metrics <- function(so,
     silhouette <- FALSE
   }
   if (silhouette == TRUE) {
-    distance <- dist(embed_mat)
+    distance <- stats::dist(embed_mat)
     cluster_scores <- cbind(cluster_scores, vector(length = length(cluster_list)))
     colnames(cluster_scores)[3] <- "Silhouette"
   }

R/convert_mtx_to_h5.R

@@ -17,9 +17,6 @@
 #' @param barcodes_file Path to the barcodes.tsv.gz file
 #' @param features_file Path to the features.tsv.gz file
 #'
-#' @import DropletUtils
-#' @import Seurat
-#'
 #' @export
 convert_mtx_to_h5 <- function(sample_name,
                               mtx_file,

R/fetch_celldex_ref.R

@@ -0,0 +1,61 @@
+#' fetch_celldex_ref: Helper function to import annotation references from
+#' celldex
+#'
+#' @description Retrieves one of the 7 available references in celldex. Should
+#' be cache-agnostic by declaring the cache location as the local working
+#' directory
+#'
+#' @param ref_name Character string to be used for retrieving references from
+#' the celldex package
+#' @param cache Character string to identify cache download location
+#'
+#' @export
+#'
+#' @return A celldex single cell object in SCE data structure
+fetch_celldex_ref <- function(ref_name, cache = ".") {
+  ref <- switch(ref_name,
+    "hpca" = ,
+    "HumanPrimaryCellAtlasData" = celldex::fetchReference(
+      "hpca",
+      version = "2024-02-26",
+      realize.assays = TRUE, cache = "./"
+    ),
+    "blueprint_encode" = ,
+    "BP_encode" = ,
+    "bpencode" = ,
+    "bp_encode" = ,
+    "BlueprintEncodeData" = celldex::fetchReference(
+      "blueprint_encode", "2024-02-26",
+      realize.assays = TRUE, cache = "./"
+    ),
+    "monaco" = ,
+    "MonacoImmuneData" = celldex::fetchReference(
+      "monaco_immune", "2024-02-26",
+      realize.assays = TRUE, cache = "./"
+    ),
+    "immu_cell_exp" = ,
+    "DatabaseImmuneCellExpressionData" = ,
+    "dice" = celldex::fetchReference(
+      "dice", "2024-02-26",
+      realize.assays = TRUE, cache = "./"
+    ),
+    "hematopoietic" = ,
+    "NovershternHematopoieticData" = ,
+    "novershtern_hematopoietic" = ,
+    "novershtern" = celldex::fetchReference(
+      "novershtern_hematopoietic", "2024-02-26",
+      realize.assays = TRUE, cache = "./"
+    ),
+    "immgen" = ,
+    "ImmGenData" = celldex::fetchReference(
+      "immgen", "2024-02-26",
+      realize.assays = TRUE, cache = "./"
+    ),
+    "mouseRNAseq" = ,
+    "MouseRNAseqData" = celldex::fetchReference(
+      "mouse_rnaseq", "2024-02-26",
+      realize.assays = TRUE, cache = "./"
+    )
+  )
+  return(ref)
+}

R/filter_doublets.R

@@ -13,9 +13,6 @@
 #' "DoubletFinder" (default), "scDblFinder", "consensus" removal of doublets,
 #' or the "union" of total doublets identified in both algorithms
 #'
-#' @import DoubletFinder
-#' @import scDblFinder
-#'
 #' @export
 #'
 #' @return A subsetted Seurat object with doublets removed
@@ -56,8 +53,7 @@ filter_doublets <- function(so_in, doublet_finder_method = "DoubletFinder") {
     so_in <- subset(so_in, cells = colnames(so_in)[which(so_in$doubletFinder_label == "Singlet")])
   } else if (doublet_finder_method == "scDblFinder") {
     sce <- Seurat::as.SingleCellExperiment(so_in)
-    set.seed(42)
-    sce_dbl <- scDblFinder(sce) %>% suppressWarnings()
+    sce_dbl <- scDblFinder::scDblFinder(sce) %>% suppressWarnings()
     so_in$scDblFinder_label <- sce_dbl$scDblFinder.class
     so_in <- subset(so_in, cells = colnames(so_in)[which(so_in$scDblFinder_label == "singlet")])
   } else if (doublet_finder_method == "union") {
@@ -95,7 +91,6 @@ filter_doublets <- function(so_in, doublet_finder_method = "DoubletFinder") {
     so_in$doubletFinder_label <- dfso[[utils::tail(names(dfso@meta.data), 1)]]
 
     sce <- Seurat::as.SingleCellExperiment(so_in)
-    set.seed(42)
     sce_dbl <- scDblFinder::scDblFinder(sce) %>% suppressWarnings()
     so_in$scDblFinder_label <- sce_dbl$scDblFinder.class
     so_in <- subset(so_in, cells = colnames(so_in)[intersect(
@@ -137,8 +132,7 @@ filter_doublets <- function(so_in, doublet_finder_method = "DoubletFinder") {
     so_in$doubletFinder_label <- dfso[[utils::tail(names(dfso@meta.data), 1)]]
 
     sce <- Seurat::as.SingleCellExperiment(so_in)
-    set.seed(42)
-    sce_dbl <- scDblFinder(sce) %>% suppressWarnings()
+    sce_dbl <- scDblFinder::scDblFinder(sce) %>% suppressWarnings()
     so_in$scDblFinder_label <- sce_dbl$scDblFinder.class
     so_in <- subset(so_in, cells = colnames(so_in)[unique(c(
       which(so_in$doubletFinder_label == "Singlet"),

R/make_bubble_plot.R

@@ -1,4 +1,4 @@
-#' make a flexible bubble plot with Seurat::DotPlot
+#' make_bubble_plot: Flexible dotplot
 #'
 #' @description Creates a bubble plot that shows average expression percent
 #' expression for genes in single cell categories
@@ -8,56 +8,48 @@
 #' @param palette A color palette from ggplot2
 #' @param assay The counts assay to use for determining expression
 #' @param ident The categorical identity to classify groups of cells
+#' @param scale Boolean to scale average expression across identities
 #'
 #' @export
 #'
 #' @return A ggplot2 figure
 #'
-make_bubble_plot <- function(so,
-                             features,
-                             palette = "RdBu",
-                             assay = "SCT",
-                             ident = "seurat_clusters") {
-  # data variables must be initialized to silence the R CMD check note:
-  #    'no visible binding for global variable'
-  pct.exp <- id <- features.plot <- group <- Gene <- Pct.Exp <- AvgExp <- NULL
+make_bubble_plot <- function(
+  so, features, palette = "RdBu", assay = "SCT",
+  scale = FALSE, ident = "seurat_clusters"
+) {
+  pct.exp <- id <- features.plot <- NULL
+  Group <- Gene <- PctExp <- AvgExp <- NULL
 
   Seurat::Idents(so) <- ident
 
   dotplot <- Seurat::DotPlot(so, features = features)
   pct_expression <- dplyr::select(dotplot$data, pct.exp, id, features.plot)
 
-  avg_expression <- Seurat::AverageExpression(so, assays = assay)[[assay]]
+  avg_expression <- Seurat::AverageExpression(so, assay = assay, features = features)[[assay]]
   avg_expression <- as.matrix(avg_expression)
+  colnames(avg_expression) <- gsub("-", "_", colnames(avg_expression))
+  if (isTRUE(scale)) {
+    avg_expression <- t(scale(t(avg_expression)))
+  }
   avg_expression_df <- reshape2::melt(avg_expression)
-  colnames(avg_expression_df) <- c("group", "Gene", "AvgExp")
+  colnames(avg_expression_df) <- c("Gene", "Group", "AvgExp")
   avg_expression_df$PctExp <- 0
-
-  # TODO refactor this with map
-  for (i in seq_len(avg_expression_df)) {
+  # avg_expression_df$Group <- as.integer(gsub("[^0-9]", "", avg_expression_df$Group))
+  for (i in seq_len(nrow(avg_expression_df))) {
     avg_expression_df[i, 4] <- pct_expression$pct.exp[which(
-      pct_expression$id == avg_expression_df[ident, i] &
+      pct_expression$id == avg_expression_df[i, "Group"] &
         pct_expression$features.plot == avg_expression_df$Gene[i]
     )]
   }
 
-  bubble_plot <- ggplot2::ggplot(
-    avg_expression_df,
-    ggplot2::aes(
-      x = group,
-      y = Gene,
-      size = Pct.Exp,
-      color = AvgExp
-    )
-  ) +
+  bubble_plot <- ggplot2::ggplot(avg_expression_df, ggplot2::aes(
+    x = Group, y = Gene, size = PctExp, color = AvgExp
+  )) +
     ggplot2::geom_point() +
     ggplot2::scale_color_distiller(palette = palette) +
     ggplot2::theme_bw() +
-    ggplot2::theme(axis.text.x = ggplot2::element_text(
-      angle = 45,
-      vjust = 1,
-      hjust = 1
-    ))
+    ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 45, vjust = 1, hjust = 1))
 
   return(bubble_plot)
 }