Merge branch 'R-devel' into R-release

Merge in stable 2022 dev build.

R/RcppExports.R

@@ -1099,6 +1099,14 @@ assess_label_enrichment <- function(G, M, thread_no = 0L) {
     .Call(`_ACTIONet_assess_label_enrichment`, G, M, thread_no)
 }
 
+autocorrelation_Moran_parametric_full <- function(G, scores, normalization_method = 4L, thread_no = 0L) {
+    .Call(`_ACTIONet_autocorrelation_Moran_parametric_full`, G, scores, normalization_method, thread_no)
+}
+
+autocorrelation_Moran_parametric <- function(G, scores, normalization_method = 4L, thread_no = 0L) {
+    .Call(`_ACTIONet_autocorrelation_Moran_parametric`, G, scores, normalization_method, thread_no)
+}
+
 roll_var <- function(X) {
     .Call(`_ACTIONet_roll_var`, X)
 }

R/utils_public.R

@@ -16,20 +16,21 @@ normalize.matrix <- function(S,
                              top_features_frac = 1.0,
                              scale_param = median,
                              transformation = "log",
-                             anchor_features = NULL) {
+                             anchor_features = NULL,
+                             post_rescale = FALSE) {
   if (!is.matrix(S) && !ACTIONetExperiment:::is.sparseMatrix(S)) {
     err <- sprintf("`S` must be `matrix` or `sparseMatrix`.\n")
     stop(err)
   }
 
   if (!is.null(anchor_features)) {
-    lib_sizes <- Matrix::colMeans(S[anchor_features, ])
+    lib_sizes <- Matrix::colSums(S[anchor_features, ])
   } else if (top_features_frac < 1.0) {
     universality <- Matrix::rowMeans(S != 0)
     selected_features <- which(universality > (1.0 - top_features_frac))
-    lib_sizes <- Matrix::colMeans(S[selected_features, ])
+    lib_sizes <- Matrix::colSums(S[selected_features, ])
   } else {
-    lib_sizes <- Matrix::colMeans(S)
+    lib_sizes <- Matrix::colSums(S)
   }
 
   if (is.matrix(S)) {
@@ -47,35 +48,51 @@ normalize.matrix <- function(S,
     S_scaled_norm <- S_scaled * kappa
   } else if (length(scale_param) == 1) {
     S_scaled_norm <- S_scaled * scale_param
-  } else {
+  } else if (length(scale_param) == ncol(S)) {
     if (is.matrix(S_scaled)) {
-      S_scaled_norm <- S_scaled %*% diag(lib_sizes)
+      S_scaled_norm <- S_scaled %*% diag(x = scale_param)
     } else {
       S_scaled_norm <- S_scaled %*% Diagonal(n = length(scale_param), x = scale_param)
     }
+  } else {
+    stop(sprintf("invalid scale_param in the normalize.matrix() function."))
   }
 
   if (transformation == "log") {
-    S_scaled_norm <- log1p(S_scaled_norm)
+    S_scaled_norm_trans <- log1p(S_scaled_norm)
   } else if (transformation == "tukey") {
     if (is.matrix(S_scaled_norm)) {
-      idx <- which(S_scaled_norm > 0)
-      vv <- S_scaled_norm[idx]
+      S_scaled_norm_trans <- S_scaled_norm
+      idx <- which(S_scaled_norm_trans > 0)
+      vv <- S_scaled_norm_trans[idx]
       vv_transformed <- sqrt(vv) + sqrt(1 + vv)
-      S_scaled_norm[idx] <- vv_transformed
+      S_scaled_norm_trans[idx] <- vv_transformed
     } else {
-      S_scaled_norm@x <- sqrt(S_scaled_norm@x) + sqrt(S_scaled_norm@x + 1)
+      S_scaled_norm_trans@x <- sqrt(S_scaled_norm_trans@x) + sqrt(S_scaled_norm_trans@x + 1)
     }
   } else if (transformation == "lsi") {
     if (is.matrix(S_scaled_norm)) {
       S_scaled_norm <- as(S_scaled_norm, "dMatrix")
     }
-    S_scaled_norm <- ACTIONet::LSI(S_sparse)
+    S_scaled_norm_trans <- ACTIONet::LSI(S_scaled_norm)
   }
 
-  return(S_scaled_norm)
+  if (post_rescale == TRUE) {
+    cs <- Matrix::colSums(S_scaled_norm_trans)
+    kappa <- median(cs) / cs
+
+    if (is.matrix(S_scaled_norm_trans)) {
+      S_scaled_norm_trans <- S_scaled_norm_trans %*% diag(x = kappa)
+    } else {
+      S_scaled_norm_trans <- S_scaled_norm_trans %*% Diagonal(n = length(kappa), x = kappa)
+    }
+  }
+
+  return(S_scaled_norm_trans)
 }
 
+
+
 .groupedRowSums <- function(S, group_vec) {
   if (ACTIONetExperiment:::is.sparseMatrix(S)) {
     mat <- compute_grouped_rowsums(S, sample_assignments = group_vec)

README.md

@@ -141,7 +141,7 @@ require(ACTIONet)
 ace = import.ace.from.10X.h5('pbmc_10k_v3.h5', prefilter = T, min_cells_per_feat = 0.01, min_feats_per_cell = 1000)
 ace = normalize.ace(ace)
 ace = reduce.ace(ace)
-ace = run.ACTIONet(ace)
+ace = runACTIONet(ace)
 
 # Annotate cell-types
 data("curatedMarkers_human")

inst/include/ACTIONet_RcppExports.h

@@ -2355,6 +2355,48 @@ namespace ACTIONet {
         return Rcpp::as<mat >(rcpp_result_gen);
     }
 
+    inline List autocorrelation_Moran_parametric_full(mat G, mat scores, int normalization_method = 4, int thread_no = 0) {
+        typedef SEXP(*Ptr_autocorrelation_Moran_parametric_full)(SEXP,SEXP,SEXP,SEXP);
+        static Ptr_autocorrelation_Moran_parametric_full p_autocorrelation_Moran_parametric_full = NULL;
+        if (p_autocorrelation_Moran_parametric_full == NULL) {
+            validateSignature("List(*autocorrelation_Moran_parametric_full)(mat,mat,int,int)");
+            p_autocorrelation_Moran_parametric_full = (Ptr_autocorrelation_Moran_parametric_full)R_GetCCallable("ACTIONet", "_ACTIONet_autocorrelation_Moran_parametric_full");
+        }
+        RObject rcpp_result_gen;
+        {
+            RNGScope RCPP_rngScope_gen;
+            rcpp_result_gen = p_autocorrelation_Moran_parametric_full(Shield<SEXP>(Rcpp::wrap(G)), Shield<SEXP>(Rcpp::wrap(scores)), Shield<SEXP>(Rcpp::wrap(normalization_method)), Shield<SEXP>(Rcpp::wrap(thread_no)));
+        }
+        if (rcpp_result_gen.inherits("interrupted-error"))
+            throw Rcpp::internal::InterruptedException();
+        if (Rcpp::internal::isLongjumpSentinel(rcpp_result_gen))
+            throw Rcpp::LongjumpException(rcpp_result_gen);
+        if (rcpp_result_gen.inherits("try-error"))
+            throw Rcpp::exception(Rcpp::as<std::string>(rcpp_result_gen).c_str());
+        return Rcpp::as<List >(rcpp_result_gen);
+    }
+
+    inline List autocorrelation_Moran_parametric(sp_mat G, mat scores, int normalization_method = 4, int thread_no = 0) {
+        typedef SEXP(*Ptr_autocorrelation_Moran_parametric)(SEXP,SEXP,SEXP,SEXP);
+        static Ptr_autocorrelation_Moran_parametric p_autocorrelation_Moran_parametric = NULL;
+        if (p_autocorrelation_Moran_parametric == NULL) {
+            validateSignature("List(*autocorrelation_Moran_parametric)(sp_mat,mat,int,int)");
+            p_autocorrelation_Moran_parametric = (Ptr_autocorrelation_Moran_parametric)R_GetCCallable("ACTIONet", "_ACTIONet_autocorrelation_Moran_parametric");
+        }
+        RObject rcpp_result_gen;
+        {
+            RNGScope RCPP_rngScope_gen;
+            rcpp_result_gen = p_autocorrelation_Moran_parametric(Shield<SEXP>(Rcpp::wrap(G)), Shield<SEXP>(Rcpp::wrap(scores)), Shield<SEXP>(Rcpp::wrap(normalization_method)), Shield<SEXP>(Rcpp::wrap(thread_no)));
+        }
+        if (rcpp_result_gen.inherits("interrupted-error"))
+            throw Rcpp::internal::InterruptedException();
+        if (Rcpp::internal::isLongjumpSentinel(rcpp_result_gen))
+            throw Rcpp::LongjumpException(rcpp_result_gen);
+        if (rcpp_result_gen.inherits("try-error"))
+            throw Rcpp::exception(Rcpp::as<std::string>(rcpp_result_gen).c_str());
+        return Rcpp::as<List >(rcpp_result_gen);
+    }
+
 }
 
 #endif // RCPP_ACTIONet_RCPPEXPORTS_H_GEN_

src/ACTIONet.cpp

@@ -2605,4 +2605,39 @@ mat assess_label_enrichment(sp_mat &G, mat &M, int thread_no = 0)
   mat logPvals = ACTIONet::assess_label_enrichment(G, M, thread_no);
 
   return (logPvals);
-}
+}
+
+// [[Rcpp::depends(RcppArmadillo)]]
+// [[Rcpp::export]]
+List autocorrelation_Moran_parametric_full(mat G, mat scores,
+                                           int normalization_method = 4,
+                                           int thread_no = 0)
+{
+  field<vec> out = ACTIONet::autocorrelation_Moran_parametric(G, scores, normalization_method, thread_no);
+
+  List res;
+  res["stat"] = out[0];
+  res["zscore"] = out[1];
+  res["mu"] = out[2];
+  res["sigma"] = out[3];
+
+  return (res);
+}
+
+// [[Rcpp::depends(RcppArmadillo)]]
+// [[Rcpp::export]]
+List autocorrelation_Moran_parametric(sp_mat G, mat scores,
+                                      int normalization_method = 4,
+                                      int thread_no = 0)
+{
+
+  field<vec> out = ACTIONet::autocorrelation_Moran_parametric(G, scores, normalization_method, thread_no);
+
+  List res;
+  res["stat"] = out[0];
+  res["zscore"] = out[1];
+  res["mu"] = out[2];
+  res["sigma"] = out[3];
+
+  return (res);
+}

src/RcppExports.cpp

@@ -4192,6 +4192,80 @@ RcppExport SEXP _ACTIONet_assess_label_enrichment(SEXP GSEXP, SEXP MSEXP, SEXP t
     UNPROTECT(1);
     return rcpp_result_gen;
 }
+// autocorrelation_Moran_parametric_full
+List autocorrelation_Moran_parametric_full(mat G, mat scores, int normalization_method, int thread_no);
+static SEXP _ACTIONet_autocorrelation_Moran_parametric_full_try(SEXP GSEXP, SEXP scoresSEXP, SEXP normalization_methodSEXP, SEXP thread_noSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::traits::input_parameter< mat >::type G(GSEXP);
+    Rcpp::traits::input_parameter< mat >::type scores(scoresSEXP);
+    Rcpp::traits::input_parameter< int >::type normalization_method(normalization_methodSEXP);
+    Rcpp::traits::input_parameter< int >::type thread_no(thread_noSEXP);
+    rcpp_result_gen = Rcpp::wrap(autocorrelation_Moran_parametric_full(G, scores, normalization_method, thread_no));
+    return rcpp_result_gen;
+END_RCPP_RETURN_ERROR
+}
+RcppExport SEXP _ACTIONet_autocorrelation_Moran_parametric_full(SEXP GSEXP, SEXP scoresSEXP, SEXP normalization_methodSEXP, SEXP thread_noSEXP) {
+    SEXP rcpp_result_gen;
+    {
+        Rcpp::RNGScope rcpp_rngScope_gen;
+        rcpp_result_gen = PROTECT(_ACTIONet_autocorrelation_Moran_parametric_full_try(GSEXP, scoresSEXP, normalization_methodSEXP, thread_noSEXP));
+    }
+    Rboolean rcpp_isInterrupt_gen = Rf_inherits(rcpp_result_gen, "interrupted-error");
+    if (rcpp_isInterrupt_gen) {
+        UNPROTECT(1);
+        Rf_onintr();
+    }
+    bool rcpp_isLongjump_gen = Rcpp::internal::isLongjumpSentinel(rcpp_result_gen);
+    if (rcpp_isLongjump_gen) {
+        Rcpp::internal::resumeJump(rcpp_result_gen);
+    }
+    Rboolean rcpp_isError_gen = Rf_inherits(rcpp_result_gen, "try-error");
+    if (rcpp_isError_gen) {
+        SEXP rcpp_msgSEXP_gen = Rf_asChar(rcpp_result_gen);
+        UNPROTECT(1);
+        Rf_error(CHAR(rcpp_msgSEXP_gen));
+    }
+    UNPROTECT(1);
+    return rcpp_result_gen;
+}
+// autocorrelation_Moran_parametric
+List autocorrelation_Moran_parametric(sp_mat G, mat scores, int normalization_method, int thread_no);
+static SEXP _ACTIONet_autocorrelation_Moran_parametric_try(SEXP GSEXP, SEXP scoresSEXP, SEXP normalization_methodSEXP, SEXP thread_noSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::traits::input_parameter< sp_mat >::type G(GSEXP);
+    Rcpp::traits::input_parameter< mat >::type scores(scoresSEXP);
+    Rcpp::traits::input_parameter< int >::type normalization_method(normalization_methodSEXP);
+    Rcpp::traits::input_parameter< int >::type thread_no(thread_noSEXP);
+    rcpp_result_gen = Rcpp::wrap(autocorrelation_Moran_parametric(G, scores, normalization_method, thread_no));
+    return rcpp_result_gen;
+END_RCPP_RETURN_ERROR
+}
+RcppExport SEXP _ACTIONet_autocorrelation_Moran_parametric(SEXP GSEXP, SEXP scoresSEXP, SEXP normalization_methodSEXP, SEXP thread_noSEXP) {
+    SEXP rcpp_result_gen;
+    {
+        Rcpp::RNGScope rcpp_rngScope_gen;
+        rcpp_result_gen = PROTECT(_ACTIONet_autocorrelation_Moran_parametric_try(GSEXP, scoresSEXP, normalization_methodSEXP, thread_noSEXP));
+    }
+    Rboolean rcpp_isInterrupt_gen = Rf_inherits(rcpp_result_gen, "interrupted-error");
+    if (rcpp_isInterrupt_gen) {
+        UNPROTECT(1);
+        Rf_onintr();
+    }
+    bool rcpp_isLongjump_gen = Rcpp::internal::isLongjumpSentinel(rcpp_result_gen);
+    if (rcpp_isLongjump_gen) {
+        Rcpp::internal::resumeJump(rcpp_result_gen);
+    }
+    Rboolean rcpp_isError_gen = Rf_inherits(rcpp_result_gen, "try-error");
+    if (rcpp_isError_gen) {
+        SEXP rcpp_msgSEXP_gen = Rf_asChar(rcpp_result_gen);
+        UNPROTECT(1);
+        Rf_error(CHAR(rcpp_msgSEXP_gen));
+    }
+    UNPROTECT(1);
+    return rcpp_result_gen;
+}
 // roll_var
 vec roll_var(vec& X);
 RcppExport SEXP _ACTIONet_roll_var(SEXP XSEXP) {
@@ -4333,6 +4407,8 @@ static int _ACTIONet_RcppExport_validate(const char* sig) {
         signatures.insert("List(*aggregate_genesets)(sp_mat&,sp_mat&,sp_mat&,int,double,int)");
         signatures.insert("mat(*run_harmony)(mat&,mat&,vec,int,double,double,int,int,double,double,double,double,double,bool,int)");
         signatures.insert("mat(*assess_label_enrichment)(sp_mat&,mat&,int)");
+        signatures.insert("List(*autocorrelation_Moran_parametric_full)(mat,mat,int,int)");
+        signatures.insert("List(*autocorrelation_Moran_parametric)(sp_mat,mat,int,int)");
     }
     return signatures.find(sig) != signatures.end();
 }
@@ -4450,6 +4526,8 @@ RcppExport SEXP _ACTIONet_RcppExport_registerCCallable() {
     R_RegisterCCallable("ACTIONet", "_ACTIONet_aggregate_genesets", (DL_FUNC)_ACTIONet_aggregate_genesets_try);
     R_RegisterCCallable("ACTIONet", "_ACTIONet_run_harmony", (DL_FUNC)_ACTIONet_run_harmony_try);
     R_RegisterCCallable("ACTIONet", "_ACTIONet_assess_label_enrichment", (DL_FUNC)_ACTIONet_assess_label_enrichment_try);
+    R_RegisterCCallable("ACTIONet", "_ACTIONet_autocorrelation_Moran_parametric_full", (DL_FUNC)_ACTIONet_autocorrelation_Moran_parametric_full_try);
+    R_RegisterCCallable("ACTIONet", "_ACTIONet_autocorrelation_Moran_parametric", (DL_FUNC)_ACTIONet_autocorrelation_Moran_parametric_try);
     R_RegisterCCallable("ACTIONet", "_ACTIONet_RcppExport_validate", (DL_FUNC)_ACTIONet_RcppExport_validate);
     return R_NilValue;
 }
@@ -4566,6 +4644,8 @@ static const R_CallMethodDef CallEntries[] = {
     {"_ACTIONet_aggregate_genesets", (DL_FUNC) &_ACTIONet_aggregate_genesets, 6},
     {"_ACTIONet_run_harmony", (DL_FUNC) &_ACTIONet_run_harmony, 15},
     {"_ACTIONet_assess_label_enrichment", (DL_FUNC) &_ACTIONet_assess_label_enrichment, 3},
+    {"_ACTIONet_autocorrelation_Moran_parametric_full", (DL_FUNC) &_ACTIONet_autocorrelation_Moran_parametric_full, 4},
+    {"_ACTIONet_autocorrelation_Moran_parametric", (DL_FUNC) &_ACTIONet_autocorrelation_Moran_parametric, 4},
     {"_ACTIONet_roll_var", (DL_FUNC) &_ACTIONet_roll_var, 1},
     {"_ACTIONet_computeSparseRowVariances", (DL_FUNC) &_ACTIONet_computeSparseRowVariances, 4},
     {"_ACTIONet_RcppExport_registerCCallable", (DL_FUNC) &_ACTIONet_RcppExport_registerCCallable, 0},