@@ -172,7 +172,7 @@ setMethod("celda_G",

         verbose = TRUE) {

         # Convert to sparse matrix

-        x <- methods::as(x, "dgCMatrix")

+        x <- methods::as(x, "CsparseMatrix")

         ls <- list()

         ls[[useAssay]] <- x

@@ -11,7 +11,7 @@

 #'  Rows represent features and columns represent cells. Alternatively,

 #'  any matrix-like object that can be coerced to a sparse matrix of class

 #'  "dgCMatrix" can be directly used as input. The matrix will automatically be

-#'  converted to a \linkS4class{SingleCellExperiment} object. 

+#'  converted to a \linkS4class{SingleCellExperiment} object.

 #' @param useAssay A string specifying the name of the

 #'  \link{assay} slot to use. Default "counts".

 #' @param altExpName The name for the \link{altExp} slot

@@ -173,7 +173,7 @@ setMethod("celda_G",

         # Convert to sparse matrix

         x <- methods::as(x, "dgCMatrix")

-      

+

         ls <- list()

         ls[[useAssay]] <- x

         sce <- SingleCellExperiment::SingleCellExperiment(assays = ls)

@@ -172,7 +172,7 @@ setMethod("celda_G",

         verbose = TRUE) {

         # Convert to sparse matrix

-        x <- as(x, "dgCMatrix")

+        x <- methods::as(x, "dgCMatrix")

         ls <- list()

         ls[[useAssay]] <- x

@@ -6,10 +6,12 @@

 #'  it exists. Otherwise, the \code{useAssay}

 #'  \link{assay} slot in \code{x} will be used if

 #'  \code{x} is a \linkS4class{SingleCellExperiment} object.

-#' @param x A numeric \link{matrix} of counts or a

-#'  \linkS4class{SingleCellExperiment}

+#' @param x A \linkS4class{SingleCellExperiment}

 #'  with the matrix located in the assay slot under \code{useAssay}.

-#'  Rows represent features and columns represent cells.

+#'  Rows represent features and columns represent cells. Alternatively,

+#'  any matrix-like object that can be coerced to a sparse matrix of class

+#'  "dgCMatrix" can be directly used as input. The matrix will automatically be

+#'  converted to a \linkS4class{SingleCellExperiment} object. 

 #' @param useAssay A string specifying the name of the

 #'  \link{assay} slot to use. Default "counts".

 #' @param altExpName The name for the \link{altExp} slot

@@ -149,7 +151,7 @@ setMethod("celda_G",

 #' @rdname celda_G

 #' @export

 setMethod("celda_G",

-    signature(x = "matrix"),

+    signature(x = "ANY"),

     function(x,

         useAssay = "counts",

         altExpName = "featureSubset",

@@ -169,6 +171,9 @@ setMethod("celda_G",

         logfile = NULL,

         verbose = TRUE) {

+        # Convert to sparse matrix

+        x <- as(x, "dgCMatrix")

+      

         ls <- list()

         ls[[useAssay]] <- x

         sce <- SingleCellExperiment::SingleCellExperiment(assays = ls)

@@ -753,7 +753,7 @@ setMethod("celda_G",

     counts <- SummarizedExperiment::assay(sce, i = useAssay)

     counts <- .processCounts(counts)

-    y <- SummarizedExperiment::rowData(sce)$celda_feature_module

+    y <- as.integer(SummarizedExperiment::rowData(sce)$celda_feature_module)

     L <- S4Vectors::metadata(sce)$celda_parameters$L

     beta <- S4Vectors::metadata(sce)$celda_parameters$beta

     delta <- S4Vectors::metadata(sce)$celda_parameters$delta

@@ -753,6 +753,13 @@ setMethod("celda_G",

     counts <- SummarizedExperiment::assay(sce, i = useAssay)

     counts <- .processCounts(counts)

+    y <- SummarizedExperiment::rowData(sce)$celda_feature_module

+    L <- S4Vectors::metadata(sce)$celda_parameters$L

+    beta <- S4Vectors::metadata(sce)$celda_parameters$beta

+    delta <- S4Vectors::metadata(sce)$celda_parameters$delta

+    gamma <- S4Vectors::metadata(sce)$celda_parameters$gamma

+    cNames <- colnames(sce)

+    rNames <- rownames(sce)

     if (is.null(maxCells) || maxCells > ncol(counts)) {

         maxCells <- ncol(counts)

@@ -761,7 +768,16 @@ setMethod("celda_G",

         cellIx <- sample(seq(ncol(counts)), maxCells)

     }

-    fm <- factorizeMatrix(x = sce, useAssay = useAssay, type = "counts")

+    fm <- .factorizeMatrixG(

+        counts = counts,

+        y = y,

+        L = L,

+        beta = beta,

+        delta = delta,

+        gamma = gamma,

+        cNames = cNames,

+        rNames = rNames,

+        type = "counts")

     modulesToUse <- seq(nrow(fm$counts$cell))

     if (!is.null(modules)) {

@@ -48,9 +48,8 @@

 #' @param countChecksum Character. An MD5 checksum for the `counts` matrix.

 #'  Default NULL.

 #' @param logfile Character. Messages will be redirected to a file named

-#'  `logfile`. If NULL, messages will be printed to stdout.  Default NULL.

+#'  \code{logfile}. If NULL, messages will be printed to stdout. Default NULL.

 #' @param verbose Logical. Whether to print log messages. Default TRUE.

-#' @param ... Ignored. Placeholder to prevent check warning.

 #' @return A \linkS4class{SingleCellExperiment} object. Function

 #'  parameter settings are stored in the \link{metadata}

 #'  \code{"celda_parameters"} slot. Column \code{celda_feature_module} in

@@ -62,7 +61,25 @@

 #' data(celdaGSim)

 #' sce <- celda_G(celdaGSim$counts, L = celdaGSim$L, nchains = 1)

 #' @export

-setGeneric("celda_G", function(x, ...) {

+setGeneric("celda_G",

+    function(x,

+        useAssay = "counts",

+        altExpName = "featureSubset",

+        L,

+        beta = 1,

+        delta = 1,

+        gamma = 1,

+        stopIter = 10,

+        maxIter = 200,

+        splitOnIter = 10,

+        splitOnLast = TRUE,

+        seed = 12345,

+        nchains = 3,

+        yInitialize = c("split", "random", "predefined"),

+        countChecksum = NULL,

+        yInit = NULL,

+        logfile = NULL,

+        verbose = TRUE) {

     standardGeneric("celda_G")})

@@ -744,7 +744,7 @@ setMethod("celda_G",

         cellIx <- sample(seq(ncol(counts)), maxCells)

     }

-    fm <- .factorizeMatrixCelda_G(sce, useAssay = useAssay, type = "counts")

+    fm <- factorizeMatrix(x = sce, useAssay = useAssay, type = "counts")

     modulesToUse <- seq(nrow(fm$counts$cell))

     if (!is.null(modules)) {

@@ -813,7 +813,7 @@ setMethod("celda_G",

     SummarizedExperiment::colData(sce)["colnames"] <-

         celdaGMod@names$column

     SummarizedExperiment::rowData(sce)["celda_feature_module"] <-

-        celdaClusters(celdaGMod)$y

+        as.factor(celdaClusters(celdaGMod)$y)

     return(sce)

 }

@@ -295,7 +295,7 @@ setMethod("celda_G",

   start.time <- Sys.time()

   ## Error checking and variable processing

-  #counts <- .processCounts(counts)

+  counts <- .processCounts(counts)

   if (is.null(countChecksum)) {

     countChecksum <- .createCountChecksum(counts)

   }

@@ -260,7 +260,7 @@ setMethod("celda_G",

     return(sce)

 }

-

+#' @importFrom Matrix colSums

 .celda_G <- function(counts,

                      L,

                      beta = 1,

@@ -304,12 +304,7 @@ setMethod("celda_G",

   allChains <- seq(nchains)

   # Pre-compute lgamma values

-  if(inherits(counts, "dgCMatrix")) {

-    cs <- Matrix::colSums(counts)

-  } else {

-    cs <- .colSums(counts, nrow(counts), ncol(counts))  

-  }

-  

+  cs <- colSums(counts)

   lgbeta <- lgamma(seq(0, max(cs)) + beta)

   lggamma <- lgamma(seq(0, nrow(counts) + L) + gamma)

   lgdelta <- c(NA, lgamma((seq(nrow(counts) + L) * delta)))

@@ -603,7 +598,7 @@ setMethod("celda_G",

   ix <- sample(seq(nG))

   for (i in ix) {

     probs[, i] <- cG_CalcGibbsProbY(index = i,

-      counts = as.numeric(counts[i,]),

+      counts = as.numeric(counts[i, ]),

       nTSbyC = nTSByC,

       nbyTS = nByTS,

       nGbyTS = nGByTS,

@@ -691,19 +686,15 @@ setMethod("celda_G",

 # cells.

 # @param y Numeric vector. Denotes feature module labels.

 # @param L Integer. Number of feature modules.

+#' @importFrom Matrix rowSums

 .cGDecomposeCounts <- function(counts, y, L) {

   if (any(y > L)) {

     stop("Assigned value of feature module greater than the total number",

         " of feature modules!")

   }

-  

+

   nTSByC <- .rowSumByGroup(counts, group = y, L = L)

-  

-  if (inherits(counts, "dgCMatrix")) {

-    nByG <- as.integer(Matrix::rowSums(counts))

-  } else {

-    nByG <- .rowSums(counts, nrow(counts), ncol(counts))

-  }

+  nByG <- rowSums(counts)

   nByTS <- .rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L)

   nGByTS <- tabulate(y, L) + 1 ## Add pseudogene to each state

   nM <- ncol(counts)

@@ -722,7 +713,6 @@ setMethod("celda_G",

 .cGReDecomposeCounts <- function(counts, y, previousY, nTSByC, nByG, L) {

   ## Recalculate counts based on new label

-  

   nTSByC <- .rowSumByGroupChange(counts, nTSByC, y, previousY, L)

   nByTS <- .rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L)

   nGByTS <- tabulate(y, L) + 1

@@ -295,7 +295,7 @@ setMethod("celda_G",

   start.time <- Sys.time()

   ## Error checking and variable processing

-  counts <- .processCounts(counts)

+  #counts <- .processCounts(counts)

   if (is.null(countChecksum)) {

     countChecksum <- .createCountChecksum(counts)

   }

@@ -304,10 +304,13 @@ setMethod("celda_G",

   allChains <- seq(nchains)

   # Pre-compute lgamma values

-  lgbeta <- lgamma(seq(0, max(.colSums(

-    counts,

-    nrow(counts), ncol(counts)

-  ))) + beta)

+  if(inherits(counts, "dgCMatrix")) {

+    cs <- Matrix::colSums(counts)

+  } else {

+    cs <- .colSums(counts, nrow(counts), ncol(counts))  

+  }

+  

+  lgbeta <- lgamma(seq(0, max(cs)) + beta)

   lggamma <- lgamma(seq(0, nrow(counts) + L) + gamma)

   lgdelta <- c(NA, lgamma((seq(nrow(counts) + L) * delta)))

@@ -694,25 +697,14 @@ setMethod("celda_G",

         " of feature modules!")

   }

-  if (inherits(counts, "matrix") & is.integer(counts)) {

-    nTSByC <- .rowSumByGroup(counts, group = y, L = L)

-    nByG <- as.integer(.rowSums(counts, nrow(counts), ncol(counts)))

-    nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1),

-                                       group = y, L = L))

-  } else if (inherits(counts, "matrix") & is.numeric(counts)) {

-    nTSByC <- .rowSumByGroupNumeric(counts, group = y, L = L)

-    nByG <- as.integer(.rowSums(counts, nrow(counts), ncol(counts)))

-    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1),

-                                              group = y, L = L)

-  } else if (inherits(counts, "dgCMatrix")) {

-    nTSByC <- rowSumByGroupSparse(counts, group = y)

+  nTSByC <- .rowSumByGroup(counts, group = y, L = L)

+  

+  if (inherits(counts, "dgCMatrix")) {

     nByG <- as.integer(Matrix::rowSums(counts))

-    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1),

-                                              group = y, L = L)

   } else {

-    stop("'counts' must be an integer, numeric, or dgCMatrix matrix.")

+    nByG <- .rowSums(counts, nrow(counts), ncol(counts))

   }

-  

+  nByTS <- .rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L)

   nGByTS <- tabulate(y, L) + 1 ## Add pseudogene to each state

   nM <- ncol(counts)

   nG <- nrow(counts)

@@ -730,18 +722,9 @@ setMethod("celda_G",

 .cGReDecomposeCounts <- function(counts, y, previousY, nTSByC, nByG, L) {

   ## Recalculate counts based on new label

-  if (inherits(counts, "matrix") & is.integer(counts)) {

-    nTSByC <- .rowSumByGroupChange(counts, nTSByC, y, previousY, L)

-    nByTS <- .rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L)

-  } else if (inherits(counts, "matrix") & is.numeric(counts)) {

-    nTSByC <- .rowSumByGroupChangeNumeric(counts, nTSByC, y, previousY, L)

-    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1), group = y, L = L)

-  } else if (inherits(counts, "dgCMatrix")) {

-    nTSByC <- rowSumByGroupChangeSparse(counts, nTSByC, y, previousY)

-    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1), group = y, L = L)

-  } else {

-    stop("'counts' must be an integer, numeric, or dgCMatrix matrix.")

-  }

+  

+  nTSByC <- .rowSumByGroupChange(counts, nTSByC, y, previousY, L)

+  nByTS <- .rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L)

   nGByTS <- tabulate(y, L) + 1

   return(list(

@@ -697,16 +697,22 @@ setMethod("celda_G",

   if (inherits(counts, "matrix") & is.integer(counts)) {

     nTSByC <- .rowSumByGroup(counts, group = y, L = L)

     nByG <- as.integer(.rowSums(counts, nrow(counts), ncol(counts)))

+    nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1),

+                                       group = y, L = L))

   } else if (inherits(counts, "matrix") & is.numeric(counts)) {

     nTSByC <- .rowSumByGroupNumeric(counts, group = y, L = L)

     nByG <- as.integer(.rowSums(counts, nrow(counts), ncol(counts)))

+    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1),

+                                              group = y, L = L)

   } else if (inherits(counts, "dgCMatrix")) {

     nTSByC <- rowSumByGroupSparse(counts, group = y)

     nByG <- as.integer(Matrix::rowSums(counts))

+    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1),

+                                              group = y, L = L)

   } else {

     stop("'counts' must be an integer, numeric, or dgCMatrix matrix.")

   }

-  nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1),group = y, L = L))  

+  

   nGByTS <- tabulate(y, L) + 1 ## Add pseudogene to each state

   nM <- ncol(counts)

   nG <- nrow(counts)

@@ -726,15 +732,16 @@ setMethod("celda_G",

   ## Recalculate counts based on new label

   if (inherits(counts, "matrix") & is.integer(counts)) {

     nTSByC <- .rowSumByGroupChange(counts, nTSByC, y, previousY, L)

+    nByTS <- .rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L)

   } else if (inherits(counts, "matrix") & is.numeric(counts)) {

     nTSByC <- .rowSumByGroupChangeNumeric(counts, nTSByC, y, previousY, L)

+    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1), group = y, L = L)

   } else if (inherits(counts, "dgCMatrix")) {

-    nTSByC <- rowSumByGroupChangeSparse(counts, nTSByC, y, previousY) 

+    nTSByC <- rowSumByGroupChangeSparse(counts, nTSByC, y, previousY)

+    nByTS <- .rowSumByGroupNumeric(matrix(nByG, ncol = 1), group = y, L = L)

   } else {

     stop("'counts' must be an integer, numeric, or dgCMatrix matrix.")

   }

-  

-  nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L))

   nGByTS <- tabulate(y, L) + 1

   return(list(

@@ -693,11 +693,20 @@ setMethod("celda_G",

     stop("Assigned value of feature module greater than the total number",

         " of feature modules!")

   }

-  nTSByC <- .rowSumByGroup(counts, group = y, L = L)

-  nByG <- as.integer(.rowSums(counts, nrow(counts), ncol(counts)))

-  nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1),

-    group = y, L = L

-  ))

+  

+  if (inherits(counts, "matrix") & is.integer(counts)) {

+    nTSByC <- .rowSumByGroup(counts, group = y, L = L)

+    nByG <- as.integer(.rowSums(counts, nrow(counts), ncol(counts)))

+  } else if (inherits(counts, "matrix") & is.numeric(counts)) {

+    nTSByC <- .rowSumByGroupNumeric(counts, group = y, L = L)

+    nByG <- as.integer(.rowSums(counts, nrow(counts), ncol(counts)))

+  } else if (inherits(counts, "dgCMatrix")) {

+    nTSByC <- rowSumByGroupSparse(counts, group = y)

+    nByG <- as.integer(Matrix::rowSums(counts))

+  } else {

+    stop("'counts' must be an integer, numeric, or dgCMatrix matrix.")

+  }

+  nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1),group = y, L = L))  

   nGByTS <- tabulate(y, L) + 1 ## Add pseudogene to each state

   nM <- ncol(counts)

   nG <- nrow(counts)

@@ -715,10 +724,17 @@ setMethod("celda_G",

 .cGReDecomposeCounts <- function(counts, y, previousY, nTSByC, nByG, L) {

   ## Recalculate counts based on new label

-  nTSByC <- .rowSumByGroupChange(counts, nTSByC, y, previousY, L)

-  nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1),

-    group = y, L = L

-  ))

+  if (inherits(counts, "matrix") & is.integer(counts)) {

+    nTSByC <- .rowSumByGroupChange(counts, nTSByC, y, previousY, L)

+  } else if (inherits(counts, "matrix") & is.numeric(counts)) {

+    nTSByC <- .rowSumByGroupChangeNumeric(counts, nTSByC, y, previousY, L)

+  } else if (inherits(counts, "dgCMatrix")) {

+    nTSByC <- rowSumByGroupChangeSparse(counts, nTSByC, y, previousY) 

+  } else {

+    stop("'counts' must be an integer, numeric, or dgCMatrix matrix.")

+  }

+  

+  nByTS <- as.integer(.rowSumByGroup(matrix(nByG, ncol = 1), group = y, L = L))

   nGByTS <- tabulate(y, L) + 1

   return(list(

@@ -415,7 +415,7 @@ setMethod("celda_G",

       )

       if (L > 2 & iter != maxIter &

         ((((numIterWithoutImprovement == stopIter &

-          !all(tempLl > ll))) & isTRUE(splitOnLast)) |

+          !all(tempLl >= ll))) & isTRUE(splitOnLast)) |

           (splitOnIter > 0 & iter %% splitOnIter == 0 &

             isTRUE(doGeneSplit)))) {

         .logMessages(date(),

@@ -599,9 +599,8 @@ setMethod("celda_G",

   probs <- matrix(NA, ncol = nG, nrow = L)

   ix <- sample(seq(nG))

   for (i in ix) {

-    probs[, i] <- cG_CalcGibbsProbY(

-      index = i,

-      counts = counts,

+    probs[, i] <- cG_CalcGibbsProbY(index = i,

+      counts = as.numeric(counts[i,]),

       nTSbyC = nTSByC,

       nbyTS = nByTS,

       nGbyTS = nGByTS,

@@ -43,7 +43,7 @@

 #'  features will be split into sqrt(L) modules and then each module will be

 #'  subsequently split into another sqrt(L) modules. With 'predefined', values

 #'  in `yInit` will be used to initialize `y`. Default 'split'.

-#' @param yInit Integer vector. Sets initial starting values of y. 

+#' @param yInit Integer vector. Sets initial starting values of y.

 #'  `yInit` can only be used when `yInitialize = 'predefined'`. Default NULL.

 #' @param countChecksum Character. An MD5 checksum for the `counts` matrix.

 #'  Default NULL.

@@ -43,9 +43,8 @@

 #'  features will be split into sqrt(L) modules and then each module will be

 #'  subsequently split into another sqrt(L) modules. With 'predefined', values

 #'  in `yInit` will be used to initialize `y`. Default 'split'.

-#' @param yInit Integer vector. Sets initial starting values of y. If NULL,

-#'  starting values for each feature will be randomly sampled from `1:L`.

-#'  `yInit` can only be used when `initialize = 'random'`. Default NULL.

+#' @param yInit Integer vector. Sets initial starting values of y. 

+#'  `yInit` can only be used when `yInitialize = 'predefined'`. Default NULL.

 #' @param countChecksum Character. An MD5 checksum for the `counts` matrix.

 #'  Default NULL.

 #' @param logfile Character. Messages will be redirected to a file named

@@ -731,20 +731,6 @@ setMethod("celda_G",

 }

-.reorderCeldaG <- function(counts, res) {

-    if (params(res)$L > 2 & isTRUE(length(unique(celdaClusters(res)$y)) > 1)) {

-        res@clusters$y <- as.integer(as.factor(celdaClusters(res)$y))

-        fm <- factorizeMatrix(counts, res)

-        uniqueY <- sort(unique(celdaClusters(res)$y))

-        cs <- prop.table(t(fm$posterior$cell[uniqueY, ]), 2)

-        d <- .cosineDist(cs)

-        h <- stats::hclust(d, method = "complete")

-        res <- .recodeClusterY(res, from = h$order, to = seq(length(h$order)))

-    }

-    return(res)

-}

-

-

 .prepareCountsForDimReductionCeldaG <- function(sce,

     useAssay,

     maxCells,

@@ -1,18 +1,18 @@

 #' @title Feature clustering with Celda

 #' @description Clusters the rows of a count matrix containing single-cell data

 #'  into L modules. The

-#'  \code{useAssay} \link[SummarizedExperiment]{assay} slot in

-#'  \code{altExpName} \link[SingleCellExperiment]{altExp} slot will be used if

+#'  \code{useAssay} \link{assay} slot in

+#'  \code{altExpName} \link{altExp} slot will be used if

 #'  it exists. Otherwise, the \code{useAssay}

-#'  \link[SummarizedExperiment]{assay} slot in \code{x} will be used if

+#'  \link{assay} slot in \code{x} will be used if

 #'  \code{x} is a \linkS4class{SingleCellExperiment} object.

 #' @param x A numeric \link{matrix} of counts or a

 #'  \linkS4class{SingleCellExperiment}

 #'  with the matrix located in the assay slot under \code{useAssay}.

 #'  Rows represent features and columns represent cells.

 #' @param useAssay A string specifying the name of the

-#'  \link[SummarizedExperiment]{assay} slot to use. Default "counts".

-#' @param altExpName The name for the \link[SingleCellExperiment]{altExp} slot

+#'  \link{assay} slot to use. Default "counts".

+#' @param altExpName The name for the \link{altExp} slot

 #'  to use. Default "featureSubset".

 #' @param L Integer. Number of feature modules.

 #' @param beta Numeric. Concentration parameter for Phi. Adds a pseudocount to

@@ -53,9 +53,9 @@

 #' @param verbose Logical. Whether to print log messages. Default TRUE.

 #' @param ... Ignored. Placeholder to prevent check warning.

 #' @return A \linkS4class{SingleCellExperiment} object. Function

-#'  parameter settings are stored in the \link[S4Vectors]{metadata}

+#'  parameter settings are stored in the \link{metadata}

 #'  \code{"celda_parameters"} slot. Column \code{celda_feature_module} in

-#'  \link[SummarizedExperiment]{rowData} contains feature modules.

+#'  \link{rowData} contains feature modules.

 #' @seealso \link{celda_C} for cell clustering and \link{celda_CG} for

 #'  simultaneous clustering of features and cells. \link{celdaGridSearch} can

 #'  be used to run multiple values of L and multiple chains in parallel.

@@ -51,6 +51,7 @@

 #' @param logfile Character. Messages will be redirected to a file named

 #'  `logfile`. If NULL, messages will be printed to stdout.  Default NULL.

 #' @param verbose Logical. Whether to print log messages. Default TRUE.

+#' @param ... Ignored. Placeholder to prevent check warning.

 #' @return A \linkS4class{SingleCellExperiment} object. Function

 #'  parameter settings are stored in the \link[S4Vectors]{metadata}

 #'  \code{"celda_parameters"} slot. Column \code{celda_feature_module} in

@@ -763,7 +763,7 @@ setMethod("celda_G",

         cellIx <- sample(seq(ncol(counts)), maxCells)

     }

-    fm <- .factorizeMatrixCelda_G(x = sce, useAssay = useAssay, type = "counts")

+    fm <- .factorizeMatrixCelda_G(sce, useAssay = useAssay, type = "counts")

     modulesToUse <- seq(nrow(fm$counts$cell))

     if (!is.null(modules)) {

@@ -1,13 +1,19 @@

 #' @title Feature clustering with Celda

 #' @description Clusters the rows of a count matrix containing single-cell data

-#'  into L modules.

+#'  into L modules. The

+#'  \code{useAssay} \link[SummarizedExperiment]{assay} slot in

+#'  \code{altExpName} \link[SingleCellExperiment]{altExp} slot will be used if

+#'  it exists. Otherwise, the \code{useAssay}

+#'  \link[SummarizedExperiment]{assay} slot in \code{x} will be used if

+#'  \code{x} is a \linkS4class{SingleCellExperiment} object.

 #' @param x A numeric \link{matrix} of counts or a

 #'  \linkS4class{SingleCellExperiment}

 #'  with the matrix located in the assay slot under \code{useAssay}.

 #'  Rows represent features and columns represent cells.

-#' @param useAssay A string specifying which \link[SummarizedExperiment]{assay}

-#'  slot to use if \code{x} is a

-#'  \linkS4class{SingleCellExperiment} object. Default "counts".

+#' @param useAssay A string specifying the name of the

+#'  \link[SummarizedExperiment]{assay} slot to use. Default "counts".

+#' @param altExpName The name for the \link[SingleCellExperiment]{altExp} slot

+#'  to use. Default "featureSubset".

 #' @param L Integer. Number of feature modules.

 #' @param beta Numeric. Concentration parameter for Phi. Adds a pseudocount to

 #'  each feature module in each cell. Default 1.

@@ -66,6 +72,7 @@ setMethod("celda_G",

     signature(x = "SingleCellExperiment"),

     function(x,

         useAssay = "counts",

+        altExpName = "featureSubset",

         L,

         beta = 1,

         delta = 1,

@@ -83,12 +90,24 @@ setMethod("celda_G",

         verbose = TRUE) {

         xClass <- "SingleCellExperiment"

-        counts <- SummarizedExperiment::assay(x, i = useAssay)

-        sce <- .celdaGWithSeed(counts = counts,

+        if (!altExpName %in% SingleCellExperiment::altExpNames(x)) {

+            stop(altExpName, " not in 'altExpNames(x)'. Run ",

+                "selectFeatures(x) first!")

+        }

+

+        altExp <- SingleCellExperiment::altExp(x, altExpName)

+

+        if (!useAssay %in% SummarizedExperiment::assayNames(altExp)) {

+            stop(useAssay, " not in assayNames(altExp(x, altExpName))")

+        }

+

+        counts <- SummarizedExperiment::assay(altExp, i = useAssay)

+

+        altExp <- .celdaGWithSeed(counts = counts,

             xClass = xClass,

             useAssay = useAssay,

-            sce = x,

+            sce = altExp,

             L = L,

             beta = beta,

             delta = delta,

@@ -104,7 +123,8 @@ setMethod("celda_G",

             yInit = yInit,

             logfile = logfile,

             verbose = verbose)

-        return(sce)

+        SingleCellExperiment::altExp(x, altExpName) <- altExp

+        return(x)

     }

 )

@@ -114,6 +134,8 @@ setMethod("celda_G",

 setMethod("celda_G",

     signature(x = "matrix"),

     function(x,

+        useAssay = "counts",

+        altExpName = "featureSubset",

         L,

         beta = 1,

         delta = 1,

@@ -130,14 +152,16 @@ setMethod("celda_G",

         logfile = NULL,

         verbose = TRUE) {

+        ls <- list()

+        ls[[useAssay]] <- x

+        sce <- SingleCellExperiment::SingleCellExperiment(assays = ls)

+        SingleCellExperiment::altExp(sce, altExpName) <- sce

         xClass <- "matrix"

-        useAssay <- NULL

-        sce <- SingleCellExperiment::SingleCellExperiment(

-            assays = list(counts = x))

-        sce <- .celdaGWithSeed(counts = x,

+

+        altExp <- .celdaGWithSeed(counts = x,

             xClass = xClass,

             useAssay = useAssay,

-            sce = sce,

+            sce = SingleCellExperiment::altExp(sce, altExpName),

             L = L,

             beta = beta,

             delta = delta,

@@ -153,6 +177,7 @@ setMethod("celda_G",

             yInit = yInit,

             logfile = logfile,

             verbose = verbose)

+        SingleCellExperiment::altExp(sce, altExpName) <- altExp

         return(sce)

     }

 )

@@ -738,7 +763,7 @@ setMethod("celda_G",

         cellIx <- sample(seq(ncol(counts)), maxCells)

     }

-    fm <- factorizeMatrix(x = sce, useAssay = useAssay, type = "counts")

+    fm <- .factorizeMatrixCelda_G(x = sce, useAssay = useAssay, type = "counts")

     modulesToUse <- seq(nrow(fm$counts$cell))

     if (!is.null(modules)) {

@@ -723,7 +723,10 @@ setMethod("celda_G",

     useAssay,

     maxCells,

     minClusterSize,

-    modules) {

+    modules,

+    normalize,

+    scaleFactor,

+    transformationFun) {

     counts <- SummarizedExperiment::assay(sce, i = useAssay)

     counts <- .processCounts(counts)

@@ -750,9 +753,9 @@ setMethod("celda_G",

     }

     norm <- t(normalizeCounts(fm$counts$cell[modulesToUse, cellIx],

-        normalize = "proportion",

-        transformationFun = sqrt

-    ))

+        normalize = normalize,

+        scaleFactor = scaleFactor,

+        transformationFun = transformationFun))

     return(list(norm = norm, cellIx = cellIx))

 }

@@ -7,7 +7,7 @@

 #'  Rows represent features and columns represent cells.

 #' @param useAssay A string specifying which \link[SummarizedExperiment]{assay}

 #'  slot to use if \code{x} is a

-#'  \link[SingleCellExperiment]{SingleCellExperiment} object. Default "counts".

+#'  \linkS4class{SingleCellExperiment} object. Default "counts".

 #' @param L Integer. Number of feature modules.

 #' @param beta Numeric. Concentration parameter for Phi. Adds a pseudocount to

 #'  each feature module in each cell. Default 1.

@@ -45,8 +45,10 @@

 #' @param logfile Character. Messages will be redirected to a file named

 #'  `logfile`. If NULL, messages will be printed to stdout.  Default NULL.

 #' @param verbose Logical. Whether to print log messages. Default TRUE.

-#' @return An object of class `celda_G` with the feature module clusters stored

-#'  in `y`.

+#' @return A \linkS4class{SingleCellExperiment} object. Function

+#'  parameter settings are stored in the \link[S4Vectors]{metadata}

+#'  \code{"celda_parameters"} slot. Column \code{celda_feature_module} in

+#'  \link[SummarizedExperiment]{rowData} contains feature modules.

 #' @seealso \link{celda_C} for cell clustering and \link{celda_CG} for

 #'  simultaneous clustering of features and cells. \link{celdaGridSearch} can

 #'  be used to run multiple values of L and multiple chains in parallel.

@@ -724,6 +724,7 @@ setMethod("celda_G",

     modules) {

     counts <- SummarizedExperiment::assay(sce, i = useAssay)

+    counts <- .processCounts(counts)

     if (is.null(maxCells) || maxCells > ncol(counts)) {

         maxCells <- ncol(counts)

@@ -54,7 +54,8 @@

 #' data(celdaGSim)

 #' sce <- celda_G(celdaGSim$counts, L = celdaGSim$L, nchains = 1)

 #' @export

-setGeneric("celda_G", function(x, ...) {standardGeneric("celda_G")})

+setGeneric("celda_G", function(x, ...) {

+    standardGeneric("celda_G")})

 #' @rdname celda_G

@@ -799,4 +800,3 @@ setMethod("celda_G",

     return(sce)

 }

-

@@ -497,7 +497,7 @@ setMethod("celda_G",

   }

   bestResult <- methods::new("celda_G",

-    celdaClusters = list(y = yBest),

+    clusters = list(y = yBest),

     params = list(

       L = as.integer(L),

       beta = beta,

@@ -703,10 +703,10 @@ setMethod("celda_G",

 .reorderCeldaG <- function(counts, res) {

-    if (params(res)$L > 2 & isTRUE(length(unique(res@celdaClusters$y)) > 1)) {

-        res@celdaClusters$y <- as.integer(as.factor(res@celdaClusters$y))

+    if (params(res)$L > 2 & isTRUE(length(unique(celdaClusters(res)$y)) > 1)) {

+        res@clusters$y <- as.integer(as.factor(celdaClusters(res)$y))

         fm <- factorizeMatrix(counts, res)

-        uniqueY <- sort(unique(res@celdaClusters$y))

+        uniqueY <- sort(unique(celdaClusters(res)$y))

         cs <- prop.table(t(fm$posterior$cell[uniqueY, ]), 2)

         d <- .cosineDist(cs)

         h <- stats::hclust(d, method = "complete")

@@ -789,13 +789,13 @@ setMethod("celda_G",

         verbose = verbose,

         completeLogLik = celdaGMod@completeLogLik,

         finalLogLik = celdaGMod@finalLogLik,

-        featureModuleLevels = sort(unique(celdaGMod@celdaClusters$y)))

+        featureModuleLevels = sort(unique(celdaClusters(celdaGMod)$y)))

     SummarizedExperiment::rowData(sce)["rownames"] <- celdaGMod@names$row

     SummarizedExperiment::colData(sce)["colnames"] <-

         celdaGMod@names$column

     SummarizedExperiment::rowData(sce)["celda_feature_module"] <-

-        celdaGMod@celdaClusters$y

+        celdaClusters(celdaGMod)$y

     return(sce)

 }

@@ -613,90 +613,6 @@ setMethod("celda_G",

 }

-.factorizeMatrixCelda_G <- function(sce, useAssay, type) {

-    counts <- SummarizedExperiment::assay(sce, i = useAssay)

-    counts <- .processCounts(counts)

-    # compareCountMatrix(counts, celdaMod)

-

-    L <- S4Vectors::metadata(sce)$celda_parameters$L

-    y <- celdaModules(sce)

-    beta <- S4Vectors::metadata(sce)$celda_parameters$beta

-    delta <- S4Vectors::metadata(sce)$celda_parameters$delta

-    gamma <- S4Vectors::metadata(sce)$celda_parameters$gamma

-

-    p <- .cGDecomposeCounts(counts = counts, y = y, L = L)

-    nTSByC <- p$nTSByC

-    nByG <- p$nByG

-    nByTS <- p$nByTS

-    nGByTS <- p$nGByTS

-    nM <- p$nM

-    nG <- p$nG

-    rm(p)

-

-    nGByTS[nGByTS == 0] <- 1

-    nGByTS <- matrix(0, nrow = length(y), ncol = L)

-    nGByTS[cbind(seq(nG), y)] <- nByG

-

-    LNames <- paste0("L", seq(L))

-    colnames(nTSByC) <- colnames(sce)

-    rownames(nTSByC) <- LNames

-    colnames(nGByTS) <- LNames

-    rownames(nGByTS) <- rownames(sce)

-    names(nGByTS) <- LNames

-

-    countsList <- c()

-    propList <- c()

-    postList <- c()

-    res <- list()

-

-    if (any("counts" %in% type)) {

-        countsList <- list(

-            cell = nTSByC,

-            module = nGByTS,

-            geneDistribution = nGByTS

-        )

-        res <- c(res, list(counts = countsList))

-    }

-

-    if (any("proportion" %in% type)) {

-        ## Need to avoid normalizing cell/gene states with zero cells/genes

-        uniqueY <- sort(unique(y))

-        tempNGByTS <- nGByTS

-        tempNGByTS[, uniqueY] <- normalizeCounts(tempNGByTS[, uniqueY],

-            normalize = "proportion"

-        )

-        tempNGByTS <- nGByTS / sum(nGByTS)

-

-        propList <- list(

-            cell = normalizeCounts(nTSByC,

-                normalize = "proportion"

-            ),

-            module = tempNGByTS,

-            geneDistribution = tempNGByTS

-        )

-        res <- c(res, list(proportions = propList))

-    }

-

-    if (any("posterior" %in% type)) {

-        gs <- nGByTS

-        gs[cbind(seq(nG), y)] <- gs[cbind(seq(nG), y)] + delta

-        gs <- normalizeCounts(gs, normalize = "proportion")

-        tempNGByTS <- (nGByTS + gamma) / sum(nGByTS + gamma)

-