Merge pull request #179 from jr-leary7/dev

bumped to 0.7.9

.github/workflows/bioc-check.yaml

@@ -32,10 +32,10 @@ jobs:
           bc_res <- BiocCheck::BiocCheck()
           bc_gc_res <- BiocCheck::BiocCheckGitClone()
           n_errors <- length(bc_res$error) + length(bc_gc_res$error)
-          Sys.setenv("BC_ERRORS_TOTAL" = n_errors)
+          if (n_errors > 0) {
+            message("BiocCheck failed ...")
+            system("exit 1")
+          } else {
+            message("BiocCheck passed !")
+          }
         shell: Rscript {0}
-      - name: Check BiocCheck results
-        run: |
-          echo "total errors: $BC_ERRORS_TOTAL"
-          if [[ $BC_ERRORS_TOTAL -gt 0 ]]; then exit 1; else echo "BiocCheck passed!"; fi
-        shell: bash

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: scLANE
 Type: Package
 Title: Model gene expression dynamics with spline-based NB GLMs, GEEs, & GLMMs
-Version: 0.7.8
+Version: 0.7.9
 Authors@R: c(person(given = "Jack", family = "Leary", email = "[email protected]", role = c("aut", "cre"), comment = c(ORCID = "0009-0004-8821-3269")), 
              person(given = "Rhonda", family = "Bacher", email = "[email protected]", role = c("ctb", "fnd"), comment = c(ORCID = "0000-0001-5787-476X")))
 Description: This package uses truncated power basis spline models to build flexible, interpretable models of single cell gene expression over pseudotime or latent time. 

NEWS.md

@@ -1,3 +1,8 @@
+# Changes in version 0.7.9
+
++ Added `geneProgramDrivers()` function to compute & test correlations of expression with gene module scores.
++ Updated documentation & unit tests.
+
 # Changes in version 0.7.8
 
 + Added progress bar to `testDynamic()`. 

R/geneProgramDrivers.R

@@ -12,6 +12,7 @@
 #' @param gene.program A vector of program scores as returned by \code{\link{geneProgramScoring}}. Defaults to NULL.
 #' @param cor.method (Optional) The correlation method to be used. Defaults to "spearman".
 #' @param fdr.cutoff (Optional) The FDR threshold for determining statistical significance. Defaults to 0.01.
+#' @param p.adj.method (Optional) The method used to adjust \emph{p}-values for multiple hypothesis testing. Defaults to "holm".
 #' @return Either a \code{Seurat} or \code{SingleCellExperiment} object if \code{expr.mat} is in either form, or a data.frame containing per-cell program scores if \code{expr.mat} is a matrix.
 #' @seealso \code{\link{geneProgramScoring}}
 #' @seealso \code{\link[stats]{cor.test}}
@@ -37,9 +38,12 @@ geneProgramDrivers <- function(expr.mat = NULL,
                                genes = NULL,
                                gene.program = NULL,
                                cor.method = "spearman",
-                               fdr.cutoff = 0.01) {
+                               fdr.cutoff = 0.01,
+                               p.adj.method = "holm") {
   # check inputs
   if (is.null(expr.mat) || is.null(genes) || is.null(gene.program)) { stop("Arguments to geneProgramDrivers() are missing.") }
+  cor.method <- tolower(cor.method)
+  if (!cor.method %in% c("pearson", "spearman", "kendall")) { stop("Please specify a valid correlation metric.") }
   # set up counts matrix
   if (inherits(expr.mat, "SingleCellExperiment")) {
     counts_matrix <- SingleCellExperiment::logcounts(expr.mat)
@@ -50,22 +54,23 @@ geneProgramDrivers <- function(expr.mat = NULL,
   } else if (inherits(expr.mat, "dgCMatrix")) {
     counts_matrix <- Matrix::Matrix(expr.mat, sparse = FALSE)
   }
-  # iteratively compute correlations
+  # iteratively compute correlations & p-values
   cor_tests <- purrr::map(genes, \(g) {
     cor_res <- stats::cor.test(counts_matrix[g, ],
                                gene.program,
-                               method = "spearman",
+                               method = cor.method,
                                exact = FALSE)
     cor_df <- data.frame(gene = g,
                          corr = unname(cor_res$estimate),
                          pvalue = cor_res$p.value)
     return(cor_df)
   })
   cor_tests <- purrr::reduce(cor_tests, rbind)
+  rownames(cor_tests) <- genes
   cor_tests <- dplyr::arrange(cor_tests,
                               pvalue,
                               dplyr::desc(abs(corr))) %>%
-               dplyr::mutate(pvalue_adj = stats::p.adjust(pvalue, method = "holm")) %>%
+               dplyr::mutate(pvalue_adj = stats::p.adjust(pvalue, method = p.adj.method)) %>%
                dplyr::filter(pvalue_adj < fdr.cutoff)
   return(cor_tests)
 }

R/geneProgramScoring.R

@@ -10,6 +10,7 @@
 #' @param program.labels (Optional) A character vector specifying a label for each gene cluster. Defaults to NULL.
 #' @param n.cores (Optional) The number of cores used under the hood in \code{\link[UCell]{ScoreSignatures_UCell}}. Defaults to 2.
 #' @return Either a \code{Seurat} or \code{SingleCellExperiment} object if \code{expr.mat} is in either form, or a data.frame containing per-cell program scores if \code{expr.mat} is a matrix.
+#' @seealso \code{\link[UCell]{ScoreSignatures_UCell}}
 #' @seealso \code{\link{geneProgramDrivers}}
 #' @export
 #' @examples

R/getResultsDE.R

@@ -9,8 +9,8 @@
 #' @importFrom tidyselect everything
 #' @importFrom stats p.adjust p.adjust.methods
 #' @param test.dyn.res The nested list returned by \code{\link{testDynamic}}. Defaults to NULL.
-#' @param p.adj.method The method used to adjust \emph{p}-values for multiple hypothesis testing. Defaults to "holm".
-#' @param fdr.cutoff The FDR threshold for determining statistical significance. Defaults to 0.01.
+#' @param p.adj.method (Optional) The method used to adjust \emph{p}-values for multiple hypothesis testing. Defaults to "holm".
+#' @param fdr.cutoff (Optional) The FDR threshold for determining statistical significance. Defaults to 0.01.
 #' @return A data.frame containing differential expression results & test statistics for each gene.
 #' @export
 #' @seealso \code{\link{testDynamic}}

R/plotModels.R

@@ -278,7 +278,8 @@ plotModels <- function(test.dyn.res = NULL,
     p <- ggplot2::ggplot(counts_df, ggplot2::aes(x = PT, y = COUNT, group = ID)) +
          ggplot2::geom_point(mapping = ggplot2::aes(color = ID),
                              alpha = 0.5,
-                             size = 0.5)
+                             size = 2,
+                             stroke = 0)
     if (requireNamespace("ggh4x", quietly = TRUE)) {
       p <- p + ggh4x::facet_nested_wrap(~paste0("Lineage ", LINEAGE) + MODEL + ID,
                                         nrow = length(levels(counts_df$MODEL)),
@@ -312,7 +313,8 @@ plotModels <- function(test.dyn.res = NULL,
   } else {
     p <- ggplot2::ggplot(counts_df, ggplot2::aes(x = PT, y = COUNT, color = LINEAGE)) +
          ggplot2::geom_point(alpha = 0.5,
-                             size = 0.5,
+                             size = 2,
+                             stroke = 0,
                              show.legend = ifelse(ncol(pt) > 1, TRUE, FALSE))
     if (requireNamespace("ggh4x", quietly = TRUE)) {
       p <- p + ggh4x::facet_nested_wrap(~paste0("Lineage ", LINEAGE) + MODEL,

inst/rmarkdown/templates/Bacher_Group_HTML/skeleton/skeleton.Rmd

@@ -18,7 +18,6 @@ output:
 ```{r setup, include=FALSE}
 knitr::opts_chunk$set(echo = TRUE, 
                       comment = NA, 
-                      message = FALSE, 
                       warning = FALSE, 
                       fig.align = "center", 
                       dev = "png", 

tests/testthat/test_scLANE.R

@@ -144,7 +144,8 @@ withr::with_output_sink(tempfile(), {
                          expr.mat = sim_data,
                          plot.null = TRUE,
                          plot.glm = TRUE,
-                         plot.gam = TRUE)
+                         plot.gam = TRUE) +
+              theme_scLANE()
   plot_gee <- plotModels(test.dyn.res = gee_gene_stats,
                          gene = "ABR",
                          pt = pt_test,