@@ -1,35 +1,33 @@

 # Generated by roxygen2: do not edit by hand

+export()

 export(DeconX)

 export(available_models)

 export(celda)

 export(celdaGridSearch)

 export(celdaHeatmap)

+export(celdaPerplexity)

 export(celdaProbabilityMap)

 export(celdaTsne)

 export(celda_C)

 export(celda_CG)

 export(celda_G)

 export(clusterProbability)

-export(clustering)

+export(clusters)

 export(compareCountMatrix)

-export(completeLogLik)

-export(countChecksum)

 export(differentialExpression)

 export(distinct_colors)

 export(factorizeMatrix)

 export(featureModuleLookup)

-export(finalLogLik)

-export(getPerplexity)

-export(initialSeed)

 export(logLikelihood)

 export(logLikelihood.celda_C)

 export(logLikelihood.celda_CG)

 export(logLikelihood.celda_G)

+export(logLikelihoodHistory)

 export(matrixNames)

-export(modelPriors)

 export(moduleHeatmap)

 export(normalizeCounts)

+export(params)

 export(perplexity)

 export(plotDimReduceCluster)

 export(plotDimReduceFeature)

@@ -77,17 +77,17 @@ moduleHeatmap <- function(counts, celda.mod, feature.module = 1, top.cells = 100

   cell_ix = match(colnames(filtered_norm.counts), celda.mod@names$column)

   z.to.plot = c()

   if(methods::.hasSlot(celda.mod, "z")){

-    cell <- distinct_colors(length(unique(celda.mod@clustering$z)))[sort(unique(celda.mod@clustering$z[cell_ix]))]

-    names(cell) <- sort(unique(celda.mod@clustering$z[cell_ix]))

+    cell <- distinct_colors(length(unique(celda.mod@clusters$z)))[sort(unique(celda.mod@clusters$z[cell_ix]))]

+    names(cell) <- sort(unique(celda.mod@clusters$z[cell_ix]))

     anno_cell_colors <- list(cell = cell)

-    z.to.plot = celda.mod@clustering$z[cell.indices]

+    z.to.plot = celda.mod@clusters$z[cell.indices]

   }else{

     anno_cell_colors <- NULL

   }

   plotHeatmap(

     filtered_norm.counts,

     z = z.to.plot,

-    y = celda.mod@clustering$y[gene_ix],

+    y = celda.mod@clusters$y[gene_ix],

     scale.row = scale.row,

     color.scheme = "divergent",

     show.names.feature = show_featurenames,

@@ -1,108 +1,23 @@

 setClass("celdaModel", 

-         representation(completeLogLik = "numeric", 

-                        finalLogLik = "numeric",

-                        seed = "numeric", 

-                        count.checksum = "character",

+         representation(params = "list", # K, L, model priors, seed, checksum

                         names = "list",

-                        clustering = "list",    # K, L, z, y, etc.

-                        modelPriors = "list"))  # alpha, beta, delta, etc.

-

-#' @title Get log-likelihood history

-#' @description Retrieves the complete log-likelihood from all iterations of Gibbs sampling used to generate a celda model.

-#' 

-#' @return Numeric. The log-likelihood at each step of Gibbs sampling used to generate the model.

-#' @examples

-#' completeLogLik(celda.CG.mod)

-#' @export

-setGeneric("completeLogLik",

-           function(celda.mod){ standardGeneric("completeLogLik") })

-setMethod("completeLogLik",

-           signature=c(celda.mod="celdaModel"),

-           function(celda.mod){  celda.mod@completeLogLik  })

-

-# @title Set log-likelihood history on a celda model

-# @description Set the complete log-likelihood from all iterations of Gibbs sampling used to generate a celda model.

-# @return A celdaModel object. The provided input object, with the updated completeLogLik property.

-# @examples

-# completeLogLik(celda.CG.mod) = c(0.00, 0.01, 0.02)  # Lenth must match num.iter

-#@export

-# setGeneric("completeLogLik<-",

-#            function(celda.mod, value){ standardGeneric("completeLogLik<-") })

-# setReplaceMethod("completeLogLik", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@completeLogLik = value

-#                    celda.mod

-#                  })

+                        completeLogLik = "numeric", 

+                        finalLogLik = "numeric",

+                        clusters = "list"))  # z and or y

-#' @title Get final log-likelihood 

-#' @description Retrieves the final log-likelihood from all iterations of Gibbs sampling used to generate a celda model.

-#' 

-#' @return Numeric. The log-likelihood at the final step of Gibbs sampling used to generate the model.

+#' @title Get parameter values provided for celda model creation

+#' @description Retrieves the K/L, model priors (e.g. alpha, beta), random seed, and count matrix checksum parameters provided during the creation of the provided celda model.

+#' @return List. Contains the model-specific parameters for the provided celda model object depending on its class.

 #' @examples

-#' finalLogLik(celda.CG.mod)

+#' params(celda.CG.mod)

 #' @export

-setGeneric("finalLogLik",

-           function(celda.mod){ standardGeneric("finalLogLik") })

-setMethod("finalLogLik",

-           signature=c(celda.mod="celdaModel"),

-           function(celda.mod){  celda.mod@finalLogLik  })

+setGeneric("params",

+           function(celda.mod){ standardGeneric("params") })

+setMethod("params",

+          signature=c(celda.mod="celdaModel"),

+          function(celda.mod){  celda.mod@params  })

-# @title Set log-likelihood history on a celda model

-# @description Set the final log-likelihood of Gibbs sampling used to generate a celda model.

-# @return A celdaModel object. The provided input object, with the updated completeLogLik property.

-# @examples

-# finalLogLik(celda.CG.mod) = 0.97 

-#@export

-# setGeneric("finalLogLik<-",

-#            function(celda.mod, value){ standardGeneric("finalLogLik<-") })

-# setReplaceMethod("finalLogLik", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@finalLogLik = value

-#                    celda.mod

-#                  })

-

-

-#' @title Get seed used to generate model

-#' @description Retrieves the random seed used to generate a celda model.

-#' @return Numeric. The random seed used to generate the provided celda model.

-#' @examples

-#' initialSeed(celda.CG.mod)

-#' @export

-setGeneric("initialSeed",

-           function(celda.mod){ standardGeneric("initialSeed") })

-setMethod("initialSeed",

-           signature=c(celda.mod="celdaModel"),

-           function(celda.mod){  celda.mod@seed  })

-#@export

-# setGeneric("initialSeed<-",

-#            function(celda.mod, value){ standardGeneric("initialSeed<-") })

-# setReplaceMethod("initialSeed", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@seed = value

-#                    celda.mod

-#                  })

-

-

-#' @title Get count matrix checksum for comparison

-#' @description Retrieves the MD5 checksum of the count matrix used to generate the provided celda mdoel.

-#' @return Character. The MD5 hash of the count matrix used to generate the provided celda model.

-#' @examples

-#' countChecksum(celda.CG.mod)

-#' @export

-setGeneric("countChecksum",

-           function(celda.mod){ standardGeneric("countChecksum") })

-setMethod("countChecksum",

-           signature=c(celda.mod="celdaModel"),

-           function(celda.mod){  celda.mod@count.checksum  })

-#@export

-# setGeneric("countChecksum<-",

-#            function(celda.mod, value){ standardGeneric("countChecksum<-") })

-# setReplaceMethod("countChecksum", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@count.checksum = value

-#                    celda.mod

-#                  })

 #' @title Get feature, cell and sample names from a celda model

 #' @description Retrieves the row, column, and sample names used to generate a celda model.

@@ -113,78 +28,57 @@ setMethod("countChecksum",

 setGeneric("matrixNames",

            function(celda.mod){ standardGeneric("matrixNames") })

 setMethod("matrixNames",

-           signature=c(celda.mod="celdaModel"),

-#@export

-# setGeneric("matrixNames<-",

-#            function(celda.mod, value){ standardGeneric("matrixNames<-") })

+          signature=c(celda.mod="celdaModel"),

           function(celda.mod){  celda.mod@names  })

-# setReplaceMethod("matrixNames", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@names = value

-#                    celda.mod

-#                  })

-#' @title Get clustering parameters and outcomes from a celda model.

-#' @description Returns the K/L parameters provided for modeling, as well as the corresponding z/y results.

-#' @return List. Contains K, z (for celda_C and celda_CG models), and/or L, y (for celda_G and celda_CG models.) 

+#' @title Get log-likelihood history

+#' @description Retrieves the complete log-likelihood from all iterations of Gibbs sampling used to generate a celda model.

+#' 

+#' @return Numeric. The log-likelihood at each step of Gibbs sampling used to generate the model.

 #' @examples

-#' clustering(celda.CG.mod)

+#' logLikelihoodHistory(celda.CG.mod)

 #' @export

-setGeneric("clustering",

-           function(celda.mod){ standardGeneric("clustering")})

-setMethod("clustering", signature=c(celda.mod="celdaModel"),

-          function(celda.mod){

-            return(celda.mod@clustering)

-          })

-#@export

-# setGeneric("clustering<-",

-#            function(celda.mod, value){ standardGeneric("clustering<-") })

+setGeneric("logLikelihoodHistory",

+           function(celda.mod){ standardGeneric("logLikelihoodHistory") })

+setMethod("logLikelihoodHistory",

+           signature=c(celda.mod="celdaModel"),

+           function(celda.mod){  celda.mod@completeLogLik  })

-#' @title Get model prior parameters from a celda model.

-#' @description Returns the model priors (e.g. alpha, beta) provided at model creation for a given celda model.

-#' @return List. Contains alpha, beta (for celda_C and celda_CG models), or delta, gamma (for celda_G and celda_CG models).

+#' @title Get final log-likelihood 

+#' @description Retrieves the final log-likelihood from all iterations of Gibbs sampling used to generate a celda model.

+#' @return Numeric. The log-likelihood at the final step of Gibbs sampling used to generate the model.

 #' @examples

-#' modelPriors(celda.CG.mod)

+#' logLikelihood(celda.CG.mod)

+#' @export

+setGeneric("logLikelihood",

+           function(celda.mod){ standardGeneric("logLikelihood") })

+setMethod("logLikelihood",

+           signature=c(celda.mod="celdaModel"),

+           function(celda.mod){  celda.mod@finalLogLik  })

+

+

+#' @title Get clustering outcomes from a celda model

+#' @description Returns the z / y results corresponding to the cell / gene cluster labels determined by the provided celda model.

+#' @return List. Contains z (for celda_C and celda_CG models) and/or y (for celda_G and celda_CG models)

+#' @examples

+#' clusters(celda.CG.mod)

 #' @export

-setGeneric("modelPriors",

-           function(celda.mod){ standardGeneric("modelPriors")})

-setMethod("modelPriors", signature=c(celda.mod="celdaModel"),

+setGeneric("clusters",

+           function(celda.mod){ standardGeneric("clusters")})

+setMethod("clusters", signature=c(celda.mod="celdaModel"),

           function(celda.mod){

-            return(celda.mod@modelPriors)

+            return(celda.mod@clusters)

           })

-#@export

-# setGeneric("modelPriors<-",

-#            function(celda.mod, value){ standardGeneric("modelPriors<-") })

 setClass("celda_C",

          representation(sample.label = "factor"),

          contains = "celdaModel")

-# setReplaceMethod("clustering", "celda_C",

-#                  function(celda.mod, value){

-#                    lapply(names(value),

-#                           function(key) {

-#                             if (!is.element(key, c("K", "z"))) {

-#                               stop(paste0("Invalid parameter provided: ", key))

-#                             }

-#                             celda.mod@clustering[[key]] = value[[key]]

-#                           })

-#                  })

-# setReplaceMethod("modelPriors", "celda_C",

-#                  function(celda.mod, value){

-#                    lapply(names(value),

-#                           function(key) {

-#                             if (!is.element(key, c("alpha", "beta"))) {

-#                               stop(paste0("Invalid parameter provided: ", key))

-#                             }

-#                             celda.mod@modelPriors[[key]] = value[[key]]

-#                           })

-#                  })

-

-

-#' @title Get sample labels

+

+

+#' @title Get sample labels from a celda model

 #' @description Returns the sample labels for the count matrix provided for generation of a given celda model.

 #' @return Character. Contains the sample labels provided at model creation time, or those automatically generated by celda.

 #' @examples

@@ -195,65 +89,13 @@ setGeneric("sampleLabel",

 setMethod("sampleLabel",

            signature=c(celda.mod="celdaModel"),

            function(celda.mod){  celda.mod@sample.label  })

-#@export

-# setGeneric("sampleLabel<-",

-#            function(celda.mod, value){ standardGeneric("sampleLabel<-") })

-# setReplaceMethod("sampleLabel", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@sample.label = value

-#                    celda.mod

-#                  })

 setClass("celda_G",

          contains = "celdaModel")

-# setReplaceMethod("clustering", "celda_G",

-#                  function(celda.mod, value){

-#                    lapply(names(value),

-#                           function(key) {

-#                             if (!is.element(key, c("L", "y"))) {

-#                               stop(paste0("Invalid parameter provided: ", key))

-#                             }

-#                             celda.mod@clustering[[key]] = value[[key]]

-#                           })

-#                  })

-# setReplaceMethod("modelPriors", "celda_G",

-#                  function(celda.mod, value){

-#                    lapply(names(value),

-#                           function(key) {

-#                             if (!is.element(key, c("beta", "delta", "gamma"))) {

-#                               stop(paste0("Invalid parameter provided: ", key))

-#                             }

-#                             celda.mod@modelPriors[[key]] = value[[key]]

-#                           })

-#                  })

-

 setClass("celda_CG",

          contains = c("celda_C", "celda_G"))

-# setReplaceMethod("clustering", "celda_CG",

-#                  function(celda.mod, value){

-#                    lapply(names(value),

-#                           function(key) {

-#                             if (!is.element(key, c("K", "L", "y", "z"))) {

-#                               stop(paste0("Invalid parameter provided: ", key))

-#                             }

-#                             celda.mod@clustering[[key]] = value[[key]]

-#                           })

-#                  })

-# setReplaceMethod("modelPriors", "celda_CG",

-#                  function(celda.mod, value){

-#                    lapply(names(value),

-#                           function(key) {

-#                             if (!is.element(key, c("alpha", "beta", "delta", "gamma"))) {

-#                               stop(paste0("Invalid parameter provided: ", key))

-#                             }

-#                             celda.mod@modelPriors[[key]] = value[[key]]

-#                           })

-#                  })

-

-

-

 setClass("celdaList",

          representation(run.params = "data.frame",

@@ -261,6 +103,7 @@ setClass("celdaList",

                         count.checksum = "character",

                         perplexity = "matrix"))

+

 #' @title Get run parameters provided to `celdaGridSearch()`

 #' @description Returns details on the clustering parameters, model priors, and seeds provided to `celdaGridSearch()` when the provided celdaList was created.

 #' @return Data Frame. Contains details on the various K/L parameters, chain parameters, and final log-likelihoods derived for each model in the provided celdaList.

@@ -272,14 +115,7 @@ setGeneric("runParams",

 setMethod("runParams",

            signature=c(celda.mod="celdaList"),

            function(celda.mod){  celda.mod@run.params  })

-#@export

-# setGeneric("runParams<-",

-#            function(celda.mod, value){ standardGeneric("runParams<-") })

-# setReplaceMethod("runParams", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@run.params = value

-#                    celda.mod

-#                  })

+

 #' @title Get final celda models from a celdaList

 #' @description Returns all models generated during a `celdaGridSearch()` run.

@@ -292,38 +128,19 @@ setGeneric("resList",

 setMethod("resList",

            signature=c(celda.mod="celdaList"),

            function(celda.mod){  celda.mod@res.list  })

-#@export

-#setGeneric("resList<-",

-#           function(celda.mod, value){ standardGeneric("resList<-") })

-# setReplaceMethod("resList", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@resList = value

-#                    celda.mod

-#                  })

-

-setMethod("countChecksum",

-           signature=c(celda.mod="celdaList"),

-           function(celda.mod){  celda.mod@count.checksum  })

 #' @title Get perplexity for every model in a celdaList

 #' @description Returns perplexity for each model in a celdaList as calculated by `perplexity().`

 #' @return List. Contains one celdaModel object for each of the parameters specified in the `runParams()` of the provided celda list.

 #' @examples

-#' celda.CG.grid.model.perplexities = getPerplexity(celda.CG.grid.search.res)

+#' celda.CG.grid.model.perplexities = celdaPerplexity(celda.CG.grid.search.res)

 #' @export

-setGeneric("getPerplexity",

-           function(celda.mod){ standardGeneric("getPerplexity") })

-setMethod("getPerplexity",

+setGeneric("celdaPerplexity",

+           function(celda.mod){ standardGeneric("celdaPerplexity") })

+setMethod("celdaPerplexity",

            signature=c(celda.mod="celdaList"),

            function(celda.mod){  celda.mod@perplexity  })

-# setGeneric("setPerplexity<-",

-#            function(celda.mod, value){ standardGeneric("setPerplexity<-") })

-# setReplaceMethod("setPerplexity", "celdaModel",

-#                  function(celda.mod, value){

-#                    celda.mod@perplexity = value

-#                    celda.mod

-#                  })

 ################################################################################

@@ -154,12 +154,12 @@ celda_C = function(counts, sample.label=NULL, K, alpha=1, beta=1,

     logMessages(date(), ".. Finished chain", i, "with seed", current.seed, logfile=logfile, append=TRUE, verbose=verbose)

   }  

   best.result = methods::new("celda_C",

-                             clustering=list(z=best.result$z, K=best.result$K),

-                             modelPriors=list(alpha=best.result$alpha, beta=best.result$beta),

+                             clusters=list(z=best.result$z),

+                             params=list(K=best.result$K, alpha=best.result$alpha, beta=best.result$beta,

+                                         sample.label=best.result$sample.label, 

+                                         count.checksum=best.result$count.checksum),

                              completeLogLik=best.result$completeLogLik,

                              finalLogLik=best.result$finalLogLik, seed=best.result$seed,

-                             sample.label=best.result$sample.label, 

-                             count.checksum=best.result$count.checksum,

                              names=best.result$names)

   best.result = reorder.celda_C(counts = counts, res = best.result)

@@ -306,14 +306,16 @@ simulateCells.celda_C = function(model, S=5, C.Range=c(50, 100), N.Range=c(500,1

   cell.counts = processCounts(cell.counts) 

   names = list(row=rownames(cell.counts), column=colnames(cell.counts), 

                sample=unique(cell.sample.label))

-  result = methods::new("celda_C", clustering=list(z=z, K=K),

-                        modelPriors=list(alpha=alpha, beta=beta), seed=seed, 

-                        sample.label=cell.sample.label, names=names,

-                        count.checksum=digest::digest(cell.counts, algo="md5"))

+  result = methods::new("celda_C", clusters=list(z=z),

+                        params=list(alpha=alpha, beta=beta, seed=seed,

+                                    count.checksum=digest::digest(cell.counts, algo="md5"),

+                                    K=K),

+                        sample.label=cell.sample.label, 

+                        names=names)

   class(result) = "celda_C" 

   result = reorder.celda_C(counts = cell.counts, res = result)

-  return(list(z=result@clustering$z, counts=processCounts(cell.counts), 

+  return(list(z=result@clusters$z, counts=processCounts(cell.counts), 

               sample.label=cell.sample.label, K=K, alpha=alpha, 

               beta=beta, C.Range=C.Range, N.Range=N.Range, S=S))

 }

@@ -337,8 +339,8 @@ setMethod("factorizeMatrix",

             counts = processCounts(counts) 

             compareCountMatrix(counts, celda.mod)

-            K = celda.mod@clustering$K

-            z = celda.mod@clustering$z

+            K = celda.mod@params$K

+            z = celda.mod@clusters$z

             alpha = celda.mod@modelPriors$alpha

             beta = celda.mod@modelPriors$beta

             sample.label = celda.mod@sample.label

@@ -487,11 +489,11 @@ cC.reDecomposeCounts = function(counts, s, z, previous.z, n.G.by.CP, K) {

 setMethod("clusterProbability", 

           signature(celda.mod = "celda_C"),

           function(counts, celda.mod, log=FALSE, ...) {

-            z = celda.mod@clustering$z

+            z = celda.mod@clusters$z

             sample.label = celda.mod@sample.label

             s = as.integer(sample.label)

-            K = celda.mod@clustering$K

+            K = celda.mod@params$K

             alpha = celda.mod@modelPriors$alpha

             beta = celda.mod@modelPriors$beta

@@ -558,10 +560,10 @@ setMethod("perplexity",

 reorder.celda_C = function(counts, res){

-  if(res@clustering$K > 2 & isTRUE(length(unique(res@clustering$z)) > 1)) {

-    res@clustering$z = as.integer(as.factor(res@clustering$z))

+  if(res@params$K > 2 & isTRUE(length(unique(res@clusters$z)) > 1)) {

+    res@clusters$z = as.integer(as.factor(res@clusters$z))

     fm <- factorizeMatrix(counts = counts, celda.mod = res)

-    unique.z = sort(unique(res@clustering$z))

+    unique.z = sort(unique(res@clusters$z))

     d <- cosineDist(fm$posterior$module[,unique.z])

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

     res <- recodeClusterZ(res, from = h$order, to = 1:length(h$order))

@@ -586,7 +588,7 @@ setMethod("celdaHeatmap",

           signature(celda.mod = "celda_C"),

           function(counts, celda.mod, feature.ix, ...) {

             norm = normalizeCounts(counts, normalize="proportion", transformation.fun=sqrt)

-            plotHeatmap(norm[feature.ix,], z=celda.mod@clustering$z, ...)

+            plotHeatmap(norm[feature.ix,], z=celda.mod@clusters$z, ...)

           })

@@ -617,9 +619,9 @@ setMethod("celdaTsne",

             compareCountMatrix(counts, celda.mod)

             ## Checking if max.cells and min.cluster.size will work

-            if((max.cells < ncol(counts)) & (max.cells / min.cluster.size < celda.mod@clustering$K)) {

+            if((max.cells < ncol(counts)) & (max.cells / min.cluster.size < celda.mod@params$K)) {

               stop(paste0("Cannot distribute ", max.cells, " cells among ", 

-                          celda.mod@clustering$K, " clusters while maintaining a minumum of ", 

+                          celda.mod@params$K, " clusters while maintaining a minumum of ", 

                           min.cluster.size, 

                           " cells per cluster. Try increasing 'max.cells' or decreasing 'min.cluster.size'."))

             }

@@ -628,7 +630,7 @@ setMethod("celdaTsne",

             total.cells.to.remove = ncol(counts) - max.cells

             z.include = rep(TRUE, ncol(counts))

             if(total.cells.to.remove > 0) {

-          	z.ta = tabulate(celda.mod@clustering$z, celda.mod@clustering$K)

+          	z.ta = tabulate(celda.mod@clusters$z, celda.mod@params$K)

           	## Number of cells that can be sampled from each cluster without 

           	## going below the minimum threshold

@@ -644,7 +646,7 @@ setMethod("celdaTsne",

           	## Perform sampling for each cluster

           	for(i in which(cluster.n.to.sample > 0)) {

-          	  z.include[sample(which(celda.mod@clustering$z == i), cluster.n.to.sample[i])] = FALSE

+          	  z.include[sample(which(celda.mod@clusters$z == i), cluster.n.to.sample[i])] = FALSE

           	}

             }   

             cell.ix = which(z.include)

@@ -681,7 +683,7 @@ setMethod("celdaProbabilityMap",

             counts = processCounts(counts)

             compareCountMatrix(counts, celda.mod)

-            z.include = which(tabulate(celda.mod@clustering$z, celda.mod@clustering$K) > 0)

+            z.include = which(tabulate(celda.mod@clusters$z, celda.mod@params$K) > 0)

             level = match.arg(level)

             factorized <- factorizeMatrix(celda.mod = celda.mod, 

@@ -215,12 +215,15 @@ celda_CG = function(counts, sample.label=NULL, K, L,

   ## Peform reordering on final Z and Y assigments:

   best.result = methods::new("celda_CG", 

-                             clustering=list(z=z.best, y=y.best, K=K, L=L),

-                             modelPriors=list(alpha=alpha, beta=beta, delta=delta, gamma=gamma),

-                             completeLogLik=ll, finalLogLik=ll.best,

-                  				   seed=current.seed, 

-                  				   sample.label=sample.label, names=names,

-                  				   count.checksum=count.checksum)

+                             clusters=list(z=z.best, y=y.best),

+                             params=list(K=K, L=L, alpha=alpha, beta=beta, 

+                                         delta=delta, gamma=gamma, 

+                                         seed=current.seed,

+                                         count.checksum=count.checksum),

+                             completeLogLik=ll, 

+                             finalLogLik=ll.best,

+                  				   sample.label=sample.label, 

+                  				   names=names)

   best.result = reorder.celda_CG(counts = counts, res = best.result)

   end.time = Sys.time()

@@ -325,14 +328,15 @@ simulateCells.celda_CG = function(model, S=5, C.Range=c(50,100), N.Range=c(500,1

   names = list(row=rownames(cell.counts), column=colnames(cell.counts), 

                sample=unique(cell.sample.label))

   result = methods::new("celda_CG", 

-                        clustering=list(z=z, y=y, K=K, L=L),

-                        modelPriors=list(alpha=alpha, beta=beta, delta=delta, gamma=gamma),

-            				    seed =seed,  sample.label=cell.sample.label, names=names,

-            				    count.checksum=digest::digest(cell.counts, algo="md5"))

+                        clusters=list(z=z, y=y),

+                        params=list(K=K, L=L, alpha=alpha, beta=beta, delta=delta, 

+                                    gamma=gamma, seed=seed,

+                                    count.checksum=digest::digest(cell.counts, algo="md5")),

+            				    sample.label=cell.sample.label, names=names)

   result = reorder.celda_CG(counts = cell.counts, res = result)

-  return(list(z=result@clustering$z, y=result@clustering$y, sample.label=cell.sample.label, 

+  return(list(z=result@clusters$z, y=result@clusters$y, sample.label=cell.sample.label, 

               counts=cell.counts, K=K, L=L, C.Range=C.Range, 

               N.Range=N.Range, S=S, alpha=alpha, beta=beta, gamma=gamma, 

               delta=delta, seed=seed))

@@ -357,10 +361,10 @@ setMethod("factorizeMatrix",

             counts = processCounts(counts)

             compareCountMatrix(counts, celda.mod)

-            K = celda.mod@clustering$K

-            L = celda.mod@clustering$L

-            z = celda.mod@clustering$z

-            y = celda.mod@clustering$y

+            K = celda.mod@params$K

+            L = celda.mod@params$L

+            z = celda.mod@clusters$z

+            y = celda.mod@clusters$y

             alpha = celda.mod@modelPriors$alpha

             beta = celda.mod@modelPriors$beta

             delta = celda.mod@modelPriors$delta

@@ -566,10 +570,10 @@ setMethod("clusterProbability",

           function(counts, celda.mod, log=FALSE, ...) {

             s = as.integer(celda.mod@sample.label)

-            z = celda.mod@clustering$z

-            K = celda.mod@clustering$K  

-            y = celda.mod@clustering$y

-            L = celda.mod@clustering$L

+            z = celda.mod@clusters$z

+            K = celda.mod@params$K  

+            y = celda.mod@clusters$y

+            L = celda.mod@params$L

             alpha = celda.mod@modelPriors$alpha

             delta = celda.mod@modelPriors$delta

             beta = celda.mod@modelPriors$beta

@@ -651,10 +655,10 @@ setMethod("perplexity",

 reorder.celda_CG = function(counts, res){

   # Reorder K

-  if(res@clustering$K > 2 & isTRUE(length(unique(res@clustering$z)) > 1)) {

-    res@clustering$z = as.integer(as.factor(res@clustering$z))

+  if(res@params$K > 2 & isTRUE(length(unique(res@clusters$z)) > 1)) {

+    res@clusters$z = as.integer(as.factor(res@clusters$z))

     fm <- factorizeMatrix(counts = counts, celda.mod = res, type="posterior")

-    unique.z = sort(unique(res@clustering$z))

+    unique.z = sort(unique(res@clusters$z))

     d <- cosineDist(fm$posterior$cell.population[,unique.z])

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

@@ -662,10 +666,10 @@ reorder.celda_CG = function(counts, res){

   }  

   # Reorder L

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

-    res@clustering$y = as.integer(as.factor(res@clustering$y))

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

+    res@clusters$y = as.integer(as.factor(res@clusters$y))

     fm <- factorizeMatrix(counts = counts, celda.mod = res, type="posterior")

-    unique.y = sort(unique(res@clustering$y))

+    unique.y = sort(unique(res@clusters$y))

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

     d <- cosineDist(cs)

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

@@ -695,7 +699,7 @@ setMethod("celdaHeatmap",

             top = celda::topRank(fm$proportions$module, n=nfeatures)

             ix = unlist(top$index)

             norm = normalizeCounts(counts, normalize="proportion", transformation.fun=sqrt)

-            plotHeatmap(norm[ix,], z=celda.mod@clustering$z, y=celda.mod@clustering$y[ix], ...)

+            plotHeatmap(norm[ix,], z=celda.mod@clusters$z, y=celda.mod@clusters$y[ix], ...)

           })

@@ -723,9 +727,9 @@ setMethod("celdaTsne",

                      seed=12345, ...) {

             ## Checking if max.cells and min.cluster.size will work

-            if((max.cells < ncol(counts)) & (max.cells / min.cluster.size < celda.mod@clustering$K)) {

+            if((max.cells < ncol(counts)) & (max.cells / min.cluster.size < celda.mod@params$K)) {

               stop(paste0("Cannot distribute ", max.cells, " cells among ",

-                          celda.mod@clustering$K, " clusters while maintaining a minumum of ", 

+                          celda.mod@params$K, " clusters while maintaining a minumum of ", 

                           min.cluster.size, " cells per cluster. Try increasing 'max.cells' or decreasing 'min.cluster.size'."))

             }

@@ -745,7 +749,7 @@ setMethod("celdaTsne",

             total.cells.to.remove = ncol(counts) - max.cells

             z.include = rep(TRUE, ncol(counts))

             if(total.cells.to.remove > 0) {

-            	z.ta = tabulate(celda.mod@clustering$z, celda.mod@clustering$K)

+            	z.ta = tabulate(celda.mod@clusters$z, celda.mod@params$K)

             	## Number of cells that can be sampled from each cluster without 

             	## going below the minimum threshold

@@ -761,7 +765,7 @@ setMethod("celdaTsne",

             	## Perform sampling for each cluster

             	for(i in which(cluster.n.to.sample > 0)) {

-            	  z.include[sample(which(celda.mod@clustering$z == i), cluster.n.to.sample[i])] = FALSE

+            	  z.include[sample(which(celda.mod@clusters$z == i), cluster.n.to.sample[i])] = FALSE

             	}

             }   

             cell.ix = which(z.include)

@@ -799,8 +803,8 @@ setMethod("celdaProbabilityMap",

             level = match.arg(level)

             factorized <- factorizeMatrix(celda.mod = celda.mod, counts = counts)

-            z.include = which(tabulate(celda.mod@clustering$z, celda.mod@clustering$K) > 0)

-            y.include = which(tabulate(celda.mod@clustering$y, celda.mod@clustering$L) > 0)    

+            z.include = which(tabulate(celda.mod@clusters$z, celda.mod@params$K) > 0)

+            y.include = which(tabulate(celda.mod@clusters$y, celda.mod@params$L) > 0)    

             if(level == "cell.population") {

               pop <- factorized$proportions$cell.population[y.include,z.include,drop=FALSE]

@@ -859,7 +863,7 @@ setMethod("featureModuleLookup",

             }

             for(x in 1:length(feature)){

               if(feature[x] %in% rownames(counts)){

-                list[x] <- celda.mod@clustering$y[which(rownames(counts) == feature[x])]

+                list[x] <- celda.mod@clusters$y[which(rownames(counts) == feature[x])]

               }else{

                 list[x] <- paste0("No feature was identified matching '", feature[x], "'.")

               }

@@ -148,10 +148,10 @@ celda_G = function(counts, L, beta=1, delta=1, gamma=1,

   best.result = methods::new("celda_G", 

-                             clustering=list(y=y.best, L=L),

-                             modelPriors=list(beta=beta, delta=delta, gamma=gamma),

+                             clusters=list(y=y.best),

+                             params=list(L=L, beta=beta, delta=delta, gamma=gamma,

+                                         count.checksum=count.checksum, seed=current.seed),

                              completeLogLik=ll, finalLogLik=ll.best, 

-                             count.checksum=count.checksum, seed=current.seed, 

                              names=names)

   best.result = reorder.celda_G(counts = counts, res = best.result) 

@@ -310,13 +310,14 @@ simulateCells.celda_G = function(model, C=100, N.Range=c(500,1000), G=100,

   ## Peform reordering on final Z and Y assigments:

   cell.counts = processCounts(cell.counts)

   names = list(row=rownames(cell.counts), column=colnames(cell.counts))

-  result = methods::new("celda_G", clustering=list(y=y, L=L), 

-                        modelPriors=list(beta=beta, delta=delta, gamma=gamma), 

-                        seed=seed,  names=names, 

-                        count.checksum=digest::digest(cell.counts, algo="md5"))

+  result = methods::new("celda_G", clusters=list(y=y), 

+                        params=list(L=L, beta=beta, delta=delta, gamma=gamma,

+                                    seed=seed,

+                                    count.checksum=digest::digest(cell.counts, algo="md5")),

+                        names=names)

   result = reorder.celda_G(counts = cell.counts, res = result)  

-  return(list(y=result@clustering$y, counts=processCounts(cell.counts), L=L, 

+  return(list(y=result@clusters$y, counts=processCounts(cell.counts), L=L, 

               beta=beta, delta=delta, gamma=gamma, seed=seed))

 }

@@ -339,8 +340,8 @@ setMethod("factorizeMatrix",

             counts = processCounts(counts)

             compareCountMatrix(counts, celda.mod)

-            L = celda.mod@clustering$L

-            y = celda.mod@clustering$y

+            L = celda.mod@params$L

+            y = celda.mod@clusters$y

             beta = celda.mod@modelPriors$beta

             delta = celda.mod@modelPriors$delta

             gamma = celda.mod@modelPriors$gamma

@@ -508,8 +509,8 @@ cG.reDecomposeCounts = function(counts, y, previous.y, n.TS.by.C, n.by.G, L) {

 setMethod("clusterProbability",

            signature(celda.mod = "celda_G"),

            function(counts, celda.mod, log=FALSE, ...) {

-              y = celda.mod@clustering$y

-              L = celda.mod@clustering$L

+              y = celda.mod@clusters$y

+              L = celda.mod@params$L

               delta = celda.mod@modelPriors$delta

               beta = celda.mod@modelPriors$beta

               gamma = celda.mod@modelPriors$gamma

@@ -574,10 +575,10 @@ setMethod("perplexity",

 reorder.celda_G = function(counts, res) {

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

-    res@clustering$y = as.integer(as.factor(res@clustering$y))

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

+    res@clusters$y = as.integer(as.factor(res@clusters$y))

     fm <- factorizeMatrix(counts = counts, celda.mod = res)

-    unique.y = sort(unique(res@clustering$y))

+    unique.y = sort(unique(res@clusters$y))

     cs = prop.table(t(fm$posterior$cell[unique.y,]), 2)

     d <- cosineDist(cs)

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

@@ -606,7 +607,7 @@ setMethod("celdaHeatmap",

             top = topRank(fm$proportions$module, n=nfeatures)

             ix = unlist(top$index)

             norm = normalizeCounts(counts, normalize="proportion", transformation.fun=sqrt)

-            plotHeatmap(norm[ix,], y=celda.mod@clustering$y[ix], ...)

+            plotHeatmap(norm[ix,], y=celda.mod@clusters$y[ix], ...)

           })

 #' @title tSNE for celda_G

@@ -687,7 +688,7 @@ setMethod("featureModuleLookup",

             }

             for(x in 1:length(feature)){

               if(feature[x] %in% rownames(counts)){

-                list[x] <- celda.mod@clustering$y[which(rownames(counts) == feature[x])]

+                list[x] <- celda.mod@clusters$y[which(rownames(counts) == feature[x])]

               }else{

                 list[x] <- paste0("No feature was identified matching '", 

                                   feature[x], "'.")

@@ -103,10 +103,10 @@ recodeClusterZ = function(celda.mod, from, to) {

   if (length(setdiff(from, to)) != 0) {

     stop("All values in 'from' must have a mapping in 'to'")

   }

-  if (is.null(celda.mod@clustering$z)) {

+  if (is.null(celda.mod@clusters$z)) {

     stop("Provided celda.mod argument does not have a z attribute")

   }

-  celda.mod@clustering$z = plyr::mapvalues(celda.mod@clustering$z, from, to)

+  celda.mod@clusters$z = plyr::mapvalues(celda.mod@clusters$z, from, to)

   return(celda.mod)

 }

@@ -125,10 +125,10 @@ recodeClusterY = function(celda.mod, from, to) {

   if (length(setdiff(from, to)) != 0) {

     stop("All values in 'from' must have a mapping in 'to'")

   }

-  if (is.null(celda.mod@clustering$y)) {

+  if (is.null(celda.mod@clusters$y)) {

     stop("Provided celda.mod argument does not have a y attribute")

   }

-  celda.mod@clustering$y = plyr::mapvalues(celda.mod@clustering$y, from, to)

+  celda.mod@clusters$y = plyr::mapvalues(celda.mod@clusters$y, from, to)

   return(celda.mod)

 }

@@ -146,13 +146,13 @@ recodeClusterY = function(celda.mod, from, to) {

 #' @export

 compareCountMatrix = function(counts, celda.mod, error.on.mismatch=TRUE) {

   if (methods::.hasSlot(celda.mod, "y")) {

-    if (nrow(counts) != length(celda.mod@clustering$y)) {

+    if (nrow(counts) != length(celda.mod@clusters$y)) {

       stop(paste0("The provided celda object was generated from a counts matrix with a different number of features than the one provided."))

     }  

   }

   if (methods::.hasSlot(celda.mod, "z")) {  

-    if (ncol(counts) != length(celda.mod@clustering$z)) {

+    if (ncol(counts) != length(celda.mod@clusters$z)) {

       stop(paste0("The provided celda object was generated from a counts matrix with a different number of cells than the one provided."))

     }

   }

@@ -27,21 +27,21 @@ differentialExpression <- function(counts, celda.mod, c1, c2 = NULL, only.pos =

   compareCountMatrix(counts, celda.mod)

   if (is.null(c2)){

-    c2 <- sort(setdiff(unique(celda.mod@clustering$z),c1))

+    c2 <- sort(setdiff(unique(celda.mod@clusters$z),c1))

   }

   if (length(c1) > 1){

     cells1 <-

-      celda.mod@names$column[which(celda.mod@clustering$z %in% c1)]

+      celda.mod@names$column[which(celda.mod@clusters$z %in% c1)]

   }else{

     cells1 <-

-      celda.mod@names$column[which(celda.mod@clustering$z == c1)]

+      celda.mod@names$column[which(celda.mod@clusters$z == c1)]

   }

   if (length(c2) > 1){

     cells2 <-

-      celda.mod@names$column[which(celda.mod@clustering$z %in% c2)]

+      celda.mod@names$column[which(celda.mod@clusters$z %in% c2)]

   }else{

     cells2 <-

-      celda.mod@names$column[which(celda.mod@clustering$z == c2)]

+      celda.mod@names$column[which(celda.mod@clusters$z == c2)]

   }

   mat <- counts[,c(cells1,cells2)]

   log_normalized_mat <- normalizeCounts(mat, normalize="cpm", transformation.fun=log1p)

@@ -64,7 +64,7 @@ recursive.splitZ = function(counts, s, K, alpha, beta, min.cell = 3, seed=12345)

       if(isTRUE(cluster.split.flag[i])) {

         ix = which(overall.z == i)

         res = suppressMessages(.celda_C(counts[,ix], K=2, stop.iter = 1, split.on.iter=-1, split.on.last=FALSE, nchains=1, seed=seed, verbose=FALSE, initialize="random"))

-        cluster.splits[cbind(ix, i)] = res@clustering$z

+        cluster.splits[cbind(ix, i)] = res@clusters$z