Merge pull request #371 from taraeicher/cmdstanr_patch

Made dependency on cmdstanr conditional, fixed other minor issues

Author: taraeicher

.github/workflows/main.yml

@@ -70,7 +70,7 @@ jobs:
           python3 -m venv env
           source env/bin/activate
           echo "VIRTUAL ENV:" $VIRTUAL_ENV
-
+  
       - name: Install reticulate
         run: |
             library("reticulate")

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: netZooR
 Type: Package
 Title: A menagerie of Methods for the inference and analysis of gene regulatory networks
-Version: 1.6.2
+Version: 1.6.3
 Date: 2024-11-26
 Authors@R: c(person("Tara", "Eicher",
                         email = "teicher@hsph.harvard.edu", role = c("aut"), comment = c(ORCID = "0000-0003-1809-4458")),
@@ -71,7 +71,6 @@ Imports:
     methods,
     dplyr,
     graphics,
-    cmdstanr,
     GeneNet,
     loo,
     rARPACK,
@@ -95,8 +94,9 @@ Suggests:
     rmarkdown,
     pkgdown,
     dorothea,
+    cmdstanr
 VignetteEngine: knitr
 VignetteBuilder: knitr
-RoxygenNote: 7.3.2
+RoxygenNote: 7.3.3
 BugReports: https://github.com/netZoo/netZooR/issues
 URL: https://github.com/netZoo/netZooR, https://netzoo.github.io/

R/CONDOR.R

@@ -25,7 +25,7 @@
 #' matrices in memory. Useful on a machine with low RAM. However, runtimes
 #' are (much) longer.
 #' @param deltaQmin convergence parameter determining the minimum required increase
-#' in the modularity for each iteration. Default is min(10^{-4},1/(number of edges)),
+#' in the modularity for each iteration. Default is min(10^(-4),1/(number of edges)),
 #' with number of edges determined by \code{nrow(condor.object$edges)}. User can
 #' set this parameter by passing a numeric value to deltaQmin.
 #' @return \code{condor.object} with \code{\link{condorModularityMax}} output

R/DRAGON.R

@@ -363,14 +363,12 @@ estimate_p_values_dragon = function(r, n, p1, p2, lambdas, kappa="estimate",seed
 #' @param gradient : method for estimating parameters of p-value distribution, applies only if p-val == TRUE. default = "finite_difference"; other option = "exact"
 #' @param verbose : verbosity level (TRUE/FALSE)
 #' @return A list of model results. cov : the shrunken covariance matrix
-#' \itemize{
-#'  \item{\code{cov}}{  the shrunken covariance matrix}
-#'  \item{\code{prec}}{  the shrunken precision matrix}
-#'  \item{\code{ggm}}{ the shrunken Gaussian graphical model; matrix of partial correlations. Self-edges (diagonal elements) are set to zero.}
-#'  \item{\code{lambdas}}{  Vector of omics-specific tuning parameters (lambda1, lambda2) for \code{layer1} and \code{layer2}}
-#'  \item{\code{gammas}}{  Reparameterized tuning parameters; gamma = 1 - lambda^2}
-#'  \item{\code{risk_grid}}{  Risk grid, for assessing optimization. Grid boundaries are in terms of gamma.}
-#' }
+#' \item{\code{cov}}{  the shrunken covariance matrix}
+#' \item{\code{prec}}{  the shrunken precision matrix}
+#' \item{\code{ggm}}{ the shrunken Gaussian graphical model; matrix of partial correlations. Self-edges (diagonal elements) are set to zero.}
+#' \item{\code{lambdas}}{  Vector of omics-specific tuning parameters (lambda1, lambda2) for \code{layer1} and \code{layer2}}
+#' \item{\code{gammas}}{  Reparameterized tuning parameters; gamma = 1 - lambda^2}
+#' \item{\code{risk_grid}}{  Risk grid, for assessing optimization. Grid boundaries are in terms of gamma.}
 #' 
 #' @export
 dragon = function(layer1,layer2,pval = FALSE,gradient = "finite_difference", verbose = FALSE)

R/LIONESS.R

@@ -164,7 +164,7 @@ lionessPy <- function(expr_file, motif_file=NULL, ppi_file=NULL, computing="cpu"
 #' Options include "pearson". 
 #' @param ncores int specifying the number of cores to be used. Default is 1. 
 #' (Note: constructing panda networks can be memory-intensive, and the number of cores should take into consideration available memory.)
-#' @param union Aggregation mode between three input networks: Union (default), intersection, legacy (maps on motif network).
+#' @param mode 'legacy' refers to the processing mode in netZooPy<=0.5, 'union': takes the union of all TFs and genes across priors and fills the missing genes in the priors with zeros; 'intersection': intersects the input genes and TFs across priors and removes the missing TFs/genes. Default values is 'union'.
 #' @param ... additional arguments for panda analysis
 #' @keywords keywords
 #' @importFrom matrixStats rowSds

R/SEAHORSE.R

@@ -9,9 +9,8 @@
 #'               The measures of association 
 #'               for a numerical phenotype is Pearson correlation and
 #'               for a categorical phenotype is the p-value of an ANOVA test
-#'                
-#'
-#' Inputs:
+#'              
+#' @import fgsea
 #' @param expression : gene expression matrix (normalized, and filtered) 
 #'                     with rows as genes and columns as samples.
 #'                     Row and column names must be present.
@@ -24,7 +23,6 @@
 #'                               Types can be either "numeric" or "categorical" 
 #' @param pathways : a list of pathways (e.g. KEGG, GO, Reactome etc. 
 #'                   downloaded from http://www.gsea-msigdb.org/gsea/msigdb/human/collections.jsp)
-#'
 #' Outputs:
 #' @return results    : a list containing three objects
 #'         results$coexpression: a gene x gene Pearson correlation matrix.
@@ -41,7 +39,7 @@
 #' colnames(phenotype_data) = c("sex", "height")
 #' rownames(phenotype_data) = colnames(expression_data)
 #' phenotype_data$sex = c(rep("male", nrow(phenotype_data)/2), rep("female", nrow(phenotype_data)/2))
-#' phenotype_data$height = 65 + sample.int(10, nrow(phenotype_data), replace = T)
+#' phenotype_data$height = 65 + sample.int(10, nrow(phenotype_data), replace = TRUE)
 #' 
 #' phenotype_dictionary = c("categorical", "numeric")
 #' 
@@ -53,15 +51,46 @@
 #' # Run seahorse
 #' results <- seahorse(expression_data, phenotype_data, phenotype_dictionary, pathways)
 #'  
-#' 
-#'  
-
+#' @export
+seahorse <- function(expression, phenotype, phenotype_dictionary, pathways){
+  set.seed(0)
+  
+  results = list()
+  
+  # Compute coexpression of genes
+  results$coexpression = cor(t(expression), use="pairwise.complete.obs")
+  
+  # Compute association of gene expression with phenotypes and run GSEA
+  results$phenotype_association = list()
+  results$GSEA = list()
+  
+  for (i in 1:ncol(phenotype)){
+    pheno = phenotype[,i]
+    pheno_name = colnames(phenotype)[i]
+    
+    if (phenotype_dictionary[i] == "numeric"){
+      output_seahorse = gsea_numeric(expression, pheno, pathways)
+    }else {output_seahorse = gsea_categorical(expression, pheno, pathways)}
+    results$phenotype_association[[pheno_name]] = output_seahorse$cor
+    results$GSEA[[pheno_name]] = output_seahorse$GSEA
+  }
+  return(results)
+}
 
-# Function to run GSEA for a numeric phenotype
+#' Function to run GSEA for a numeric phenotype
+#' @param expression : gene expression matrix (normalized, and filtered) 
+#'                     with rows as genes and columns as samples.
+#'                     Row and column names must be present.
+#'                     Row names must be HGNC symbols.
+#'                     Column names must match the row names of the phenotype matrix.
+#' @param pheno : phenotype matrix
+#'                    with rows as samples and columns as phenotype variables.
+#' @param pathways : a list of pathways (e.g. KEGG, GO, Reactome etc. 
+#'                   downloaded from http://www.gsea-msigdb.org/gsea/msigdb/human/collections.jsp)
 #' @export
 #' @import fgsea
 
-gsea_numeric <- function(expression, pheno, pathways, results){
+gsea_numeric <- function(expression, pheno, pathways){
   output_seahorse = list()
   output_seahorse$cor = list()
   output_seahorse$GSEA = list()
@@ -78,10 +107,19 @@ gsea_numeric <- function(expression, pheno, pathways, results){
   return(output_seahorse)
 }
 
-# Function to run GSEA for a categorical phenotype
+#' Function to run GSEA for a categorical phenotype
 #' @import fgsea
+#' @param expression : gene expression matrix (normalized, and filtered) 
+#'                     with rows as genes and columns as samples.
+#'                     Row and column names must be present.
+#'                     Row names must be HGNC symbols.
+#'                     Column names must match the row names of the phenotype matrix.
+#' @param pheno : phenotype matrix
+#'                    with rows as samples and columns as phenotype variables.
+#' @param pathways : a list of pathways (e.g. KEGG, GO, Reactome etc. 
+#'                   downloaded from http://www.gsea-msigdb.org/gsea/msigdb/human/collections.jsp)
 #' @export
-gsea_categorical <- function(expression, pheno, pathways, results){
+gsea_categorical <- function(expression, pheno, pathways){
   output_seahorse = list()
   output_seahorse$cor = list()
   output_seahorse$GSEA = list()
@@ -96,32 +134,4 @@ gsea_categorical <- function(expression, pheno, pathways, results){
   output_seahorse$GSEA = fgseaRes
 
   return(output_seahorse)
-}
-
-# Main SEAHORSE function
-#' @import fgsea
-#' @export
-seahorse <- function(expression, phenotype, phenotype_dictionary, pathways){
-  set.seed(0)
-
-  results = list()
-  
-  # Compute coexpression of genes
-  results$coexpression = cor(t(expression), use="pairwise.complete.obs")
-  
-  # Compute association of gene expression with phenotypes and run GSEA
-  results$phenotype_association = list()
-  results$GSEA = list()
-  
-  for (i in 1:ncol(phenotype)){
-    pheno = phenotype[,i]
-    pheno_name = colnames(phenotype)[i]
-    
-    if (phenotype_dictionary[i] == "numeric"){
-      output_seahorse = gsea_numeric(expression, pheno, pathways, results)
-    }else {output_seahorse = gsea_categorical(expression, pheno, pathways, results)}
-    results$phenotype_association[[pheno_name]] = output_seahorse$cor
-    results$GSEA[[pheno_name]] = output_seahorse$GSEA
-  }
-  return(results)
 }

R/TIGER.R

@@ -53,9 +53,12 @@
 #' @export
 #'
 #' @examples
-#' data(TIGER_expr)
-#' data(TIGER_prior)
-#' tiger(TIGER_expr,TIGER_prior)
+#' if (requireNamespace("cmdstanr", quietly = TRUE) && tryCatch({
+#' nzchar(cmdstanr::cmdstan_path())}, error = function(e) FALSE)) {
+#'    data(TIGER_expr)
+#'    data(TIGER_prior)
+#'    tiger(TIGER_expr,TIGER_prior)
+#' }
 tiger = function(expr,prior,method="VB",TFexpressed = TRUE,
                  signed=TRUE,baseline=TRUE,psis_loo = FALSE,
                  seed=123,out_path=NULL,out_size = 300,

R/zzz.R

@@ -1,9 +1,3 @@
 
 .onLoad <- function(...) {
-  cmdstan_version <- cmdstanr::cmdstan_version(error_on_NA = FALSE)
-  if (is.null(cmdstan_version)) {
-    #stop("No CmdStan installation found. Run cmdstanr::install_cmdstan() to install.", call. = FALSE)
-    cmdstanr::check_cmdstan_toolchain(fix = TRUE)
-    cmdstanr::install_cmdstan()
-  }
 }