Merge branch 'basics'

* basics:
  Add BASiCS to lists
  Run checks
  Convert to SingleCellExperiment
  Fix lun2Estimate error
  Fix lun2Estimate error
  Fix default MCMC iterations
  Add estimation without spike-ins
  Set default MCMC iterations
  Update BASiCSParams print method
  Add BASiCS to listSims
  Add require check to newBASiCSParams
  Add BASiCSEstimate
  Updata BASiCSParams docs
  Add batches to BASiCSSimulate
  Add batch parameters
  Add BASiCSSimulate function
  Add nSpikes parameter
  Add BASiCSParams tests
  Add BASiCSParams class

DESCRIPTION

@@ -54,7 +54,8 @@ Suggests:
     scran,
     mfa,
     phenopath,
-    zinbwave
+    zinbwave,
+    BASiCS
 biocViews: SingleCell, RNASeq, Transcriptomics, GeneExpression, Sequencing,
     Software
 URL: https://github.com/Oshlack/splatter

NAMESPACE

 # Generated by roxygen2: do not edit by hand
 
+S3method(BASiCSEstimate,SCESet)
+S3method(BASiCSEstimate,matrix)
 S3method(lun2Estimate,SingleCellExperiment)
 S3method(lun2Estimate,matrix)
 S3method(lunEstimate,SingleCellExperiment)
@@ -17,6 +19,8 @@ S3method(splatEstimate,SingleCellExperiment)
 S3method(splatEstimate,matrix)
 S3method(zinbEstimate,SingleCellExperiment)
 S3method(zinbEstimate,matrix)
+export(BASiCSEstimate)
+export(BASiCSSimulate)
 export(addGeneLengths)
 export(compareSCEs)
 export(diffSCEs)
@@ -32,6 +36,7 @@ export(makeDiffPanel)
 export(makeOverallPanel)
 export(mfaEstimate)
 export(mfaSimulate)
+export(newBASiCSParams)
 export(newLun2Params)
 export(newLunParams)
 export(newMFAParams)
@@ -56,6 +61,7 @@ export(splatSimulateSingle)
 export(summariseDiff)
 export(zinbEstimate)
 export(zinbSimulate)
+exportClasses(BASiCSParams)
 exportClasses(Lun2Params)
 exportClasses(LunParams)
 exportClasses(MFAParams)
@@ -109,6 +115,7 @@ importFrom(methods,new)
 importFrom(methods,slot)
 importFrom(methods,slotNames)
 importFrom(methods,validObject)
+importFrom(scater,newSCESet)
 importFrom(stats,dnbinom)
 importFrom(stats,median)
 importFrom(stats,model.matrix)

R/AllClasses.R

@@ -308,7 +308,7 @@ setClass("LunParams",
 
 #' The Lun2Params class
 #'
-#' S4 class that holds parameters for the Lun simulation.
+#' S4 class that holds parameters for the Lun2 simulation.
 #'
 #' @section Parameters:
 #'
@@ -461,6 +461,90 @@ setClass("SCDDParams",
                                varInflation = c(1, 1),
                                condition = "condition"))
 
+#' The BASiCSParams class
+#'
+#' S4 class that holds parameters for the BASiCS simulation.
+#'
+#' @section Parameters:
+#'
+#' The BASiCS simulation uses the following parameters:
+#' \describe{
+#'     \item{\code{nGenes}}{The number of genes to simulate.}
+#'     \item{\code{nCells}}{The number of cells to simulate.}
+#'     \item{\code{[seed]}}{Seed to use for generating random numbers.}
+#'     \item{\emph{Batch parameters}}{
+#'         \describe{
+#'            \item{\code{nBatches}}{Number of batches to simulate.}
+#'            \item{\code{batchCells}}{Number of cells in each batch.}
+#'         }
+#'     }
+#'     \item{\emph{Gene parameters}}{
+#'         \describe{
+#'             \item{\code{gene.params}}{A \code{data.frame} containing gene
+#'             parameters with two coloumns: \code{Mean} (mean expression for
+#'             each biological gene) and \code{Delta} (cell-to-cell
+#'             heterogeneity for each biological gene).}
+#'         }
+#'     }
+#'     \item{\emph{Spike-in parameters}}{
+#'         \describe{
+#'             \item{\code{nSpikes}}{The number of spike-ins to simulate.}
+#'             \item{\code{spike.means}}{Input molecules for each spike-in.}
+#'         }
+#'     }
+#'     \item{\emph{Cell parameters}}{
+#'         \describe{
+#'             \item{\code{cell.params}}{A \code{data.frame} containing gene
+#'             parameters with two coloumns: \code{Phi} (mRNA content factor for
+#'             each cell, scaled to sum to the number of cells in each batch)
+#'             and \code{S} (capture efficient for each cell).}
+#'         }
+#'     }
+#'     \item{\emph{Variability parameters}}{
+#'         \describe{
+#'             \item{\code{theta}}{Technical variability parameter for each
+#'             batch.}
+#'         }
+#'     }
+#' }
+#'
+#' The parameters not shown in brackets can be estimated from real data using
+#' \code{\link{BASiCSEstimate}}. For details of the BASiCS simulation see
+#' \code{\link{BASiCSSimulate}}.
+#'
+#' @name BASiCSParams
+#' @rdname BASiCSParams
+#' @aliases BASiCSParams-class
+#' @exportClass BASiCSParams
+setClass("BASiCSParams",
+         contains = "Params",
+         slots = c(nBatches = "numeric",
+                   batchCells = "numeric",
+                   gene.params = "data.frame",
+                   nSpikes = "numeric",
+                   spike.means = "numeric",
+                   cell.params = "data.frame",
+                   theta = "numeric"),
+         prototype = prototype(nBatches = 1,
+                               batchCells = 100,
+                               gene.params =
+                                   data.frame(
+                                       Mean = c(8.36, 10.65, 4.88, 6.29, 21.72,
+                                                12.93, 30.19),
+                                       Delta = c(1.29, 0.88, 1.51, 1.49, 0.54,
+                                                 0.40, 0.85)
+                               ),
+                               nSpikes = 5,
+                               spike.means = c(12.93, 30.19, 1010.72, 7.90,
+                                               31.59),
+                               cell.params =
+                                   data.frame(
+                                       Phi = c(1.00, 1.06, 1.09, 1.05, 0.80),
+                                       S = c(0.38, 0.40, 0.38, 0.39, 0.34)
+                               ),
+                               theta = 0.39)
+)
+
 #' The MFAParams class
 #'
 #' S4 class that holds parameters for the mfa simulation.
@@ -468,7 +552,6 @@ setClass("SCDDParams",
 #' @section Parameters:
 #'
 #' The mfa simulation uses the following parameters:
-#'
 #' \describe{
 #'     \item{\code{nGenes}}{The number of genes to simulate.}
 #'     \item{\code{nCells}}{The number of cells to simulate.}

R/BASiCS-estimate.R

+#' Estimate BASiCS simulation parameters
+#'
+#' Estimate simulation parameters for the BASiCS simulation from a real dataset.
+#'
+#' @param counts either a counts matrix or an SCESet object containing count
+#'        data to estimate parameters from.
+#' @param spike.info data.frame describing spike-ins with two columns: "Name"
+#'        giving the names of the spike-in features (must match
+#'        \code{rownames(counts)}) and "Input" giving the number of input
+#'        molecules.
+#' @param batch vector giving the batch that each cell belongs to.
+#' @param n total number of MCMC iterations. Must be \code{>= max(4, thin)} and
+#' a multiple of \code{thin}.
+#' @param thin thining period for the MCMC sampler. Must be \code{>= 2}.
+#' @param burn burn-in period for the MCMC sampler. Must be in the range
+#' \code{1 <= burn < n} and a multiple of \code{thin}.
+#' @param params BASiCSParams object to store estimated values in.
+#' @param verbose logical. Whether to print progress messages.
+#' @param progress logical. Whether to print additional BASiCS progress
+#' messages.
+#' @param ... Optional parameters passed to \code{\link[BASiCS]{BASiCS_MCMC}}.
+#'
+#' @details
+#' This function is just a wrapper around \code{\link[BASiCS]{BASiCS_MCMC}} that
+#' takes the output and converts it to a BASiCSParams object. Either a set of
+#' spike-ins or batch information (or both) must be supplied. If only batch
+#' information is provided there must be at least two batches. See
+#' \code{\link[BASiCS]{BASiCS_MCMC}} for details.
+#'
+#' @return BASiCSParams object containing the estimated parameters.
+#'
+#' @examples
+#' \dontrun{
+#' data("sc_example_counts")
+#' spike.info <- data.frame(Name = rownames(sc_example_counts)[1:10],
+#'                          Input = rnorm(10, 500, 200),
+#'                          stringsAsFactors = FALSE)
+#' params <- BASiCSEstimate(sc_example_counts[1:50, 1:20],
+#'                          spike.info)
+#' params
+#' }
+#' @export
+BASiCSEstimate <- function(counts, spike.info = NULL, batch = NULL,
+                           n = 20000, thin = 10, burn = 5000,
+                           params = newBASiCSParams(), verbose = TRUE,
+                           progress = TRUE, ...) {
+    UseMethod("BASiCSEstimate")
+}
+
+#' @rdname BASiCSEstimate
+#' @export
+BASiCSEstimate.SCESet <- function(counts, spike.info = NULL, batch = NULL,
+                                  n = 20000, thin = 10, burn = 5000,
+                                  params = newBASiCSParams(), verbose = TRUE,
+                                  progress = TRUE, ...) {
+    counts <- BiocGenerics::counts(counts)
+    BASiCSEstimate(counts, params)
+}
+
+#' @rdname BASiCSEstimate
+#' @export
+BASiCSEstimate.matrix <- function(counts, spike.info = NULL, batch = NULL,
+                                  n = 20000, thin = 10, burn = 5000,
+                                  params = newBASiCSParams(), verbose = TRUE,
+                                  progress = TRUE, ...) {
+
+    checkmate::assertClass(params, "BASiCSParams")
+    checkmate::assertMatrix(counts, mode = "numeric", any.missing = FALSE,
+                            min.rows = 1, min.cols = 1, row.names = "unique",
+                            col.names = "unique")
+    if (is.null(spike.info) && is.null(batch)) {
+        stop("At least one of spike.info and batch must be provided")
+    }
+    if (!is.null(spike.info)) {
+        checkmate::assertDataFrame(spike.info, any.missing = FALSE,
+                                   min.rows = 1, ncols = 2)
+        if (!all(colnames(spike.info) == c("Name", "Input"))) {
+            stop("spike.info must have columns named 'Name' and 'Input'")
+        }
+        checkmate::assertCharacter(spike.info$Name, min.chars = 1,
+                                   unique = TRUE)
+        checkmate::assertNumeric(spike.info$Input, lower = 0, finite = TRUE)
+    } #else {
+    #    spike.info <- data.frame(Name = c(), Input = c())
+    #}
+    if (!is.null(batch)) {
+        checkmate::assertIntegerish(batch, lower = 0, any.missing = FALSE,
+                                    len = ncol(counts))
+        if (is.null(spike.info) && length(unique(batch)) == 1) {
+            stop("If spike.info is not provided there must be at least two ",
+                 "batches")
+        }
+    } else {
+        batch <- rep(1, ncol(counts))
+    }
+    checkmate::assertInt(thin, lower = 2)
+    checkmate::assertInt(n, lower = max(4, thin))
+    if ((n %% thin) != 0) {
+        stop("'n' must be a multiple of 'thin'")
+    }
+    checkmate::assertInt(burn, lower = 1, upper = n - 1)
+    if ((burn %% thin) != 0) {
+        stop("'burn' must be a multiple of 'thin'")
+    }
+
+    is.spike <- rownames(counts) %in% spike.info$Name
+    BASiCS.data <- suppressMessages(
+                       BASiCS::newBASiCS_Data(counts, is.spike, spike.info,
+                                              batch)
+    )
+
+    if (verbose) {
+        mcmc <- BASiCS::BASiCS_MCMC(Data = BASiCS.data, N = n, Thin = thin,
+                                    Burn = burn, PrintProgress = progress, ...)
+    } else {
+        mcmc <- suppressMessages(
+                    BASiCS::BASiCS_MCMC(Data = BASiCS.data, N = n, Thin = thin,
+                                        Burn = burn, PrintProgress = progress,
+                                        ...)
+        )
+    }
+
+    mcmc.summ <- BASiCS::Summary(mcmc)
+
+    means <- BASiCS::displaySummaryBASiCS(mcmc.summ, Param = "mu")[, 1]
+    deltas <- BASiCS::displaySummaryBASiCS(mcmc.summ, Param = "delta")[, 1]
+    phis <- BASiCS::displaySummaryBASiCS(mcmc.summ, Param = "phi")[, 1]
+    ss <- BASiCS::displaySummaryBASiCS(mcmc.summ, Param = "s")[, 1]
+    thetas <- BASiCS::displaySummaryBASiCS(mcmc.summ, Param = "theta")[, 1]
+
+    params <- setParams(params,
+                        nGenes = sum(!is.spike),
+                        batchCells = as.vector(table(batch)),
+                        gene.params = data.frame(Mean = means, Delta = deltas),
+                        nSpikes = sum(is.spike),
+                        spike.means = ifelse(!is.null(spike.info),
+                                             spike.info$Input, numeric()),
+                        cell.params = data.frame(Phi = phis, S = ss),
+                        theta = thetas)
+
+    return(params)
+}

R/BASiCS-simulate.R

+#' BASiCS simulation
+#'
+#' Simulate counts using the BASiCS method.
+#'
+#' @param params BASiCSParams object containing simulation parameters.
+#' @param verbose logical. Whether to print progress messages
+#' @param ... any additional parameter settings to override what is provided in
+#'        \code{params}.
+#'
+#' @details
+#' This function is just a wrapper around \code{\link[BASiCS]{BASiCS_Sim}} that
+#' takes a \code{\link{BASiCSParams}}, runs the simulation then converts the
+#' output to an \code{\link[scater]{SCESet}} object. See
+#' \code{\link[BASiCS]{BASiCS_Sim}} for more details of how the simulation
+#' works.
+#'
+#' @return SingleCellExperiment containing simulated counts
+#'
+#' @references
+#' Vallejos CA, Marioni JC, Richardson S. BASiCS: Bayesian Analysis of
+#' Single-Cell Sequencing data. PLoS Comput. Biol. (2015).
+#'
+#' Paper: \url{10.1371/journal.pcbi.1004333}
+#'
+#' Code: \url{https://github.com/catavallejos/BASiCS}
+#'
+#' @examples
+#' sim <- BASiCSSimulate()
+#' @export
+#' @importFrom scater newSCESet
+BASiCSSimulate <- function(params = newBASiCSParams(), verbose = TRUE, ...) {
+
+    checkmate::assertClass(params, "BASiCSParams")
+    params <- setParams(params, ...)
+    params <- expandParams(params)
+    validObject(params)
+
+    # Set random seed
+    seed <- getParam(params, "seed")
+    set.seed(seed)
+
+    if (verbose) {message("Getting parameters...")}
+    nGenes <- getParam(params, "nGenes")
+    nCells <- getParam(params, "nCells")
+    nSpikes <- getParam(params, "nSpikes")
+    nBatches <- getParam(params, "nBatches")
+    batch.cells <- getParam(params, "batchCells")
+    gene.params <- getParam(params, "gene.params")
+
+    # Sample gene.params if necessary
+    if (nrow(gene.params) != nGenes) {
+        warning("Number of gene.params not equal to nGenes, ",
+                "gene.params will be sampled.")
+        selected <- sample(nrow(gene.params), nGenes, replace = TRUE)
+        gene.params <- gene.params[selected, ]
+    }
+
+    mu <- gene.params$Mean
+    delta <- gene.params$Delta
+
+    if (nSpikes > 0) {
+        spike.mu <- getParam(params, "spike.means")
+        if (length(spike.mu) != nSpikes) {
+            warning("Number of spike-in means not equal to nSpikes, ",
+                    "spike.means will be sampled.")
+            selected <- sample(length(spike.mu), nSpikes, replace = TRUE)
+            spike.mu <- spike.mu[selected]
+        }
+    } else {
+        # Create dummy spike-ins to get around BASiCS_Sim...
+        spike.mu <- c(10, 10)
+    }
+
+    cell.params <- getParam(params, "cell.params")
+    if (nrow(cell.params) != nCells) {
+        warning("Number of cell.params not equal to nCells, ",
+                "cell.params will be sampled.")
+        selected <- sample(nrow(cell.params), nCells, replace = TRUE)
+        cell.params <- cell.params[selected, ]
+    }
+
+    thetas <- getParam(params, "theta")
+
+    batches <- lapply(seq_len(nBatches), function(i, b) {rep(i, b[i])},
+                      b = batch.cells)
+    batches <- unlist(batches)
+
+    if (verbose) {message("Simulating counts with BASiCS...")}
+    counts.list <- list()
+    for (batch in seq_len(nBatches)) {
+        batch.cells <- batches == batch
+        phi <- cell.params[batch.cells, "Phi"]
+        phi <- (phi / sum(phi)) * sum(batch.cells)
+        s <- cell.params[batch.cells, "S"]
+        theta <- thetas[batch]
+        BASiCS.sim <- suppressMessages(
+                          BASiCS::BASiCS_Sim(mu, spike.mu, delta, phi, s, theta)
+                      )
+        batch.counts <- SummarizedExperiment::assay(BASiCS.sim)
+        counts.list[[batch]] <- batch.counts
+    }
+
+    counts <- do.call(cbind, counts.list)
+
+    if (verbose) {message("Creating final dataset...")}
+    cell.names <- paste0("Cell", seq_len(nCells))
+    gene.names <- paste0("Gene", seq_len(nGenes))
+    if (nSpikes > 0) {
+        gene.names <- c(gene.names, paste0("Spike", seq_len(nSpikes)))
+    } else {
+        # Remove dummy spikes
+        counts <- counts[1:(nrow(counts) - 2), ]
+        spike.mu <- numeric()
+    }
+
+    rownames(counts) <- gene.names
+    colnames(counts) <- cell.names
+
+    cells <- data.frame(Cell = cell.names,
+                        Phi = cell.params[, "Phi"],
+                        S = cell.params[, "S"],
+                        Batch = batches,
+                        BatchTheta = thetas[batches])
+    rownames(cells) <- cell.names
+
+    features <- data.frame(Gene = gene.names,
+                           Mean = c(mu, spike.mu),
+                           Delta = c(delta, rep(NA, nSpikes)),
+                           IsSpike = c(rep(FALSE, nGenes), rep(TRUE, nSpikes)))
+    rownames(features) <- gene.names
+
+    sim <- SingleCellExperiment(assays = list(counts = counts),
+                                rowData = features,
+                                colData = cells,
+                                metadata = list(params = params))
+
+    if (verbose) {message("Done!")}
+    return(sim)
+}

R/BASiCSParams-methods.R

+#' @rdname newParams
+#' @importFrom methods new
+#' @export
+newBASiCSParams <- function(...) {
+
+    if (!requireNamespace("BASiCS", quietly = TRUE)) {
+        stop("The BASiCS simulation requires the 'BASiCS' package.")
+    }
+
+    params <- new("BASiCSParams")
+    params <- setParams(params, ...)
+
+    return(params)
+}
+
+#' @importFrom checkmate checkInt checkDataFrame checkNumeric
+setValidity("BASiCSParams", function(object) {
+
+    object <- expandParams(object)
+    v <- getParams(object, slotNames(object))
+
+    nCells <- v$nCells
+    nGenes <- v$nGenes
+    nBatches <- v$nBatches
+    checks <- c(seed = checkInt(v$seed, lower = 0),
+                nGenes = checkInt(v$nGenes, lower = 1),
+                nCells = checkInt(v$nCells, lower = 1),
+                nBatches = checkInt(v$nBatches, lower = 1),
+                batchCells = checkIntegerish(v$batchCells, lower = 1,
+                                             len = nBatches),
+                gene.params = checkDataFrame(v$gene.params,
+                                             types = "numeric",
+                                             any.missing = FALSE,
+                                             min.rows = 1, ncols = 2),
+                nSpikes = checkNumber(v$nSpikes, lower = 0, finite = TRUE),
+                spike.means = checkNumeric(v$spike.means, lower = 0,
+                                           finite = TRUE),
+                cell.params = checkDataFrame(v$cell.params,
+                                             types = "numeric",
+                                             any.missing = FALSE,
+                                             min.rows = 1, ncols = 2),
+                theta = checkNumeric(v$theta, lower = 0, len = nBatches,
+                                     finite = TRUE)
+                )
+
+    if (!all(colnames(v$gene.params) == c("Mean", "Delta"))) {
+        checks <- c(checks, gene.params = "Incorrect column names")
+    }
+
+    if (!all(colnames(v$cell.params) == c("Phi", "S"))) {
+        checks <- c(checks, cell.params = "Incorrect column names")
+    }
+
+    # Check batchCells matches nCells, nBatches
+    if (nCells != sum(v$batchCells) || nBatches != length(v$batchCells)) {
+        checks <- c(checks,
+                    "nCells, nBatches and batchCells are not consistent")
+    }
+
+    if (all(checks == TRUE)) {
+        valid <- TRUE
+    } else {
+        valid <- checks[checks != TRUE]
+        valid <- paste(names(valid), valid, sep = ": ")
+    }
+
+    return(valid)
+})
+
+#' @rdname setParam
+setMethod("setParam", "BASiCSParams",function(object, name, value) {
+    checkmate::assertString(name)
+
+    if (name == "nCells" || name == "nBatches") {
+        stop(name, " cannot be set directly, set batchCells instead")
+    }
+
+    if (name == "batchCells") {
+        object <- setParamUnchecked(object, "nCells", sum(value))
+        object <- setParamUnchecked(object, "nBatches", length(value))
+    }
+
+    object <- callNextMethod()
+
+    return(object)
+})
+
+setMethod("show", "BASiCSParams", function(object) {
+
+    pp <- list("Batches:"     = c("(Batches)"     = "nBatches",
+                                  "(Batch Cells)" = "batchCells"),
+               "Spike-ins:"   = c("(Number)"      = "nSpikes",
+                                  "(Means)"       = "spike.means"),
+               "Variability:" = c("(Theta)"       = "theta"))
+
+    callNextMethod()
+
+    gene.params <- getParam(object, "gene.params")
+    cell.params <- getParam(object, "cell.params")
+    cat("Genes:", "\n")
+    cat("(Params)", "\n")
+    cat("data.frame with", dim(gene.params)[1], "features\n")
+    print(head(gene.params, n = 3))
+    cat("  ...  ...\n\n")
+
+    cat("Cells:", "\n")
+    cat("(Params)", "\n")
+    cat("data.frame with", dim(cell.params)[1], "features\n")
+    print(head(cell.params, n = 3))
+    cat("  ...  ...\n\n")
+
+    showPP(object, pp)
+})
+
+#' @rdname expandParams
+setMethod("expandParams", "BASiCSParams", function(object) {
+
+    n <- getParam(object, "nBatches")
+
+    vectors <- c("theta")
+
+    object <- callNextMethod(object, vectors, n)
+
+    return(object)
+})

R/listSims.R

@@ -48,6 +48,11 @@ listSims <- function(print = TRUE) {
                    "The scDD simulation samples a given dataset and can
                    simulate differentially expressed and differentially
                    distributed genes between two conditions."),
+                 c("BASiCS", "BASiCS", "10.1371/journal.pcbi.1004333",
+                   "catavallejos/BASiCS",
+                   "The BASiCS simulation is based on a bayesian model used to
+                   deconvolve biological and technical variation and
+                   includes spike-ins and batch effects."),
                  c("mfa", "mfa", "10.12688/wellcomeopenres.11087.1",
                    "kieranrcampbell/mfa",
                    "The mfa simulation produces a bifurcating pseudotime

_pkgdown.yml

@@ -8,23 +8,25 @@ reference:
   - title: Parameters
     desc: Parameters functions and classes
     contents:
+      - '`newParams`'
       - '`getParam`'
       - '`getParams`'
+      - '`setParam`'
+      - '`setParams`'
+      - '`BASiCSParams`'
       - '`Lun2Params`'
       - '`LunParams`'
       - '`MFAParams`'
-      - '`newParams`'
       - '`Params`'
       - '`PhenoParams`'
       - '`SCDDParams`'
       - '`SimpleParams`'
       - '`SplatParams`'
-      - '`setParam`'
-      - '`setParams`'
       - '`ZINBParams`'
   - title: Estimation
     desc: Functions for estimating parameters
     contents:
+      - '`BASiCSEstimate`'
       - '`lun2Estimate`'
       - '`lunEstimate`'
       - '`mfaEstimate`'
@@ -41,6 +43,7 @@ reference:
   - title: Simulation
     desc: Functions for simulating datasets
     contents:
+      - '`BASiCSSimulate`'
       - '`lun2Simulate`'
       - '`lunSimulate`'
       - '`mfaSimulate`'

man/BASiCSEstimate.Rd

+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/BASiCS-estimate.R
+\name{BASiCSEstimate}
+\alias{BASiCSEstimate}
+\alias{BASiCSEstimate.SCESet}
+\alias{BASiCSEstimate.matrix}
+\title{Estimate BASiCS simulation parameters}
+\usage{
+BASiCSEstimate(counts, spike.info = NULL, batch = NULL, n = 20000,
+  thin = 10, burn = 5000, params = newBASiCSParams(), verbose = TRUE,
+  progress = TRUE, ...)
+
+\method{BASiCSEstimate}{SCESet}(counts, spike.info = NULL, batch = NULL,
+  n = 20000, thin = 10, burn = 5000, params = newBASiCSParams(),
+  verbose = TRUE, progress = TRUE, ...)
+
+\method{BASiCSEstimate}{matrix}(counts, spike.info = NULL, batch = NULL,
+  n = 20000, thin = 10, burn = 5000, params = newBASiCSParams(),
+  verbose = TRUE, progress = TRUE, ...)
+}
+\arguments{
+\item{counts}{either a counts matrix or an SCESet object containing count
+data to estimate parameters from.}
+
+\item{spike.info}{data.frame describing spike-ins with two columns: "Name"
+giving the names of the spike-in features (must match
+\code{rownames(counts)}) and "Input" giving the number of input
+molecules.}
+
+\item{batch}{vector giving the batch that each cell belongs to.}
+
+\item{n}{total number of MCMC iterations. Must be \code{>= max(4, thin)} and
+a multiple of \code{thin}.}
+
+\item{thin}{thining period for the MCMC sampler. Must be \code{>= 2}.}
+
+\item{burn}{burn-in period for the MCMC sampler. Must be in the range
+\code{1 <= burn < n} and a multiple of \code{thin}.}
+
+\item{params}{BASiCSParams object to store estimated values in.}
+
+\item{verbose}{logical. Whether to print progress messages.}
+
+\item{progress}{logical. Whether to print additional BASiCS progress
+messages.}
+
+\item{...}{Optional parameters passed to \code{\link[BASiCS]{BASiCS_MCMC}}.}
+}
+\value{
+BASiCSParams object containing the estimated parameters.
+}
+\description{
+Estimate simulation parameters for the BASiCS simulation from a real dataset.
+}
+\details{
+This function is just a wrapper around \code{\link[BASiCS]{BASiCS_MCMC}} that
+takes the output and converts it to a BASiCSParams object. Either a set of
+spike-ins or batch information (or both) must be supplied. If only batch
+information is provided there must be at least two batches. See
+\code{\link[BASiCS]{BASiCS_MCMC}} for details.
+}
+\examples{
+\dontrun{
+data("sc_example_counts")
+spike.info <- data.frame(Name = rownames(sc_example_counts)[1:10],
+                         Input = rnorm(10, 500, 200),
+                         stringsAsFactors = FALSE)
+params <- BASiCSEstimate(sc_example_counts[1:50, 1:20],
+                         spike.info)
+params
+}
+}

man/BASiCSParams.Rd

+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/AllClasses.R
+\docType{class}
+\name{BASiCSParams}
+\alias{BASiCSParams}
+\alias{BASiCSParams-class}
+\title{The BASiCSParams class}
+\description{
+S4 class that holds parameters for the BASiCS simulation.
+}
+\section{Parameters}{
+
+
+The BASiCS simulation uses the following parameters:
+\describe{
+    \item{\code{nGenes}}{The number of genes to simulate.}
+    \item{\code{nCells}}{The number of cells to simulate.}
+    \item{\code{[seed]}}{Seed to use for generating random numbers.}
+    \item{\emph{Batch parameters}}{
+        \describe{
+           \item{\code{nBatches}}{Number of batches to simulate.}
+           \item{\code{batchCells}}{Number of cells in each batch.}
+        }
+    }
+    \item{\emph{Gene parameters}}{
+        \describe{
+            \item{\code{gene.params}}{A \code{data.frame} containing gene
+            parameters with two coloumns: \code{Mean} (mean expression for
+            each biological gene) and \code{Delta} (cell-to-cell
+            heterogeneity for each biological gene).}
+        }
+    }
+    \item{\emph{Spike-in parameters}}{
+        \describe{
+            \item{\code{nSpikes}}{The number of spike-ins to simulate.}
+            \item{\code{spike.means}}{Input molecules for each spike-in.}
+        }
+    }
+    \item{\emph{Cell parameters}}{
+        \describe{
+            \item{\code{cell.params}}{A \code{data.frame} containing gene
+            parameters with two coloumns: \code{Phi} (mRNA content factor for
+            each cell, scaled to sum to the number of cells in each batch)
+            and \code{S} (capture efficient for each cell).}
+        }
+    }
+    \item{\emph{Variability parameters}}{
+        \describe{
+            \item{\code{theta}}{Technical variability parameter for each
+            batch.}
+        }
+    }
+}
+
+The parameters not shown in brackets can be estimated from real data using
+\code{\link{BASiCSEstimate}}. For details of the BASiCS simulation see
+\code{\link{BASiCSSimulate}}.
+}
+

man/BASiCSSimulate.Rd

+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/BASiCS-simulate.R
+\name{BASiCSSimulate}
+\alias{BASiCSSimulate}
+\title{BASiCS simulation}
+\usage{
+BASiCSSimulate(params = newBASiCSParams(), verbose = TRUE, ...)
+}
+\arguments{
+\item{params}{BASiCSParams object containing simulation parameters.}
+
+\item{verbose}{logical. Whether to print progress messages}
+
+\item{...}{any additional parameter settings to override what is provided in
+\code{params}.}
+}
+\value{
+SingleCellExperiment containing simulated counts
+}
+\description{
+Simulate counts using the BASiCS method.
+}
+\details{
+This function is just a wrapper around \code{\link[BASiCS]{BASiCS_Sim}} that
+takes a \code{\link{BASiCSParams}}, runs the simulation then converts the
+output to an \code{\link[scater]{SCESet}} object. See
+\code{\link[BASiCS]{BASiCS_Sim}} for more details of how the simulation
+works.
+}
+\examples{
+sim <- BASiCSSimulate()
+}
+\references{
+Vallejos CA, Marioni JC, Richardson S. BASiCS: Bayesian Analysis of
+Single-Cell Sequencing data. PLoS Comput. Biol. (2015).
+
+Paper: \url{10.1371/journal.pcbi.1004333}
+
+Code: \url{https://github.com/catavallejos/BASiCS}
+}

man/Lun2Params.Rd

@@ -6,7 +6,7 @@
 \alias{Lun2Params-class}
 \title{The Lun2Params class}
 \description{
-S4 class that holds parameters for the Lun simulation.
+S4 class that holds parameters for the Lun2 simulation.
 }
 \section{Parameters}{
 

man/MFAParams.Rd

@@ -12,7 +12,6 @@ S4 class that holds parameters for the mfa simulation.
 
 
 The mfa simulation uses the following parameters:
-
 \describe{
     \item{\code{nGenes}}{The number of genes to simulate.}
     \item{\code{nCells}}{The number of cells to simulate.}

man/expandParams.Rd

 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/AllGenerics.R, R/LunParams-methods.R,
-%   R/SplatParams-methods.R
+% Please edit documentation in R/AllGenerics.R, R/BASiCSParams-methods.R,
+%   R/LunParams-methods.R, R/SplatParams-methods.R
 \docType{methods}
 \name{expandParams}
 \alias{expandParams}
+\alias{expandParams,BASiCSParams-method}
 \alias{expandParams,LunParams-method}
 \alias{expandParams,SplatParams-method}
 \title{Expand parameters}
 \usage{
 expandParams(object, ...)
 
+\S4method{expandParams}{BASiCSParams}(object)
+
 \S4method{expandParams}{LunParams}(object)
 
 \S4method{expandParams}{SplatParams}(object)

man/newParams.Rd

 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/AllGenerics.R, R/Lun2Params-methods.R,
-%   R/LunParams-methods.R, R/MFAParams-methods.R, R/PhenoParams-methods.R,
-%   R/SCDDParams-methods.R, R/SimpleParams-methods.R, R/SplatParams-methods.R,
-%   R/ZINBParams-methods.R
+% Please edit documentation in R/AllGenerics.R, R/BASiCSParams-methods.R,
+%   R/Lun2Params-methods.R, R/LunParams-methods.R, R/MFAParams-methods.R,
+%   R/PhenoParams-methods.R, R/SCDDParams-methods.R, R/SimpleParams-methods.R,
+%   R/SplatParams-methods.R, R/ZINBParams-methods.R
 \name{newParams}
 \alias{newParams}
+\alias{newBASiCSParams}
 \alias{newLun2Params}
 \alias{newLunParams}
 \alias{newMFAParams}
@@ -15,6 +16,8 @@
 \alias{newZINBParams}
 \title{New Params}
 \usage{
+newBASiCSParams(...)
+
 newLun2Params(...)
 
 newLunParams(...)

man/setParam.Rd

 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/AllGenerics.R, R/Lun2Params-methods.R,
-%   R/LunParams-methods.R, R/Params-methods.R, R/PhenoParams-methods.R,
-%   R/SCDDParams-methods.R, R/SplatParams-methods.R, R/ZINBParams-methods.R
+% Please edit documentation in R/AllGenerics.R, R/BASiCSParams-methods.R,
+%   R/Lun2Params-methods.R, R/LunParams-methods.R, R/Params-methods.R,
+%   R/PhenoParams-methods.R, R/SCDDParams-methods.R, R/SplatParams-methods.R,
+%   R/ZINBParams-methods.R
 \docType{methods}
 \name{setParam}
 \alias{setParam}
+\alias{setParam,BASiCSParams-method}
 \alias{setParam,Lun2Params-method}
 \alias{setParam,LunParams-method}
 \alias{setParam,Params-method}
@@ -16,6 +18,8 @@
 \usage{
 setParam(object, name, value)
 
+\S4method{setParam}{BASiCSParams}(object, name, value)
+
 \S4method{setParam}{Lun2Params}(object, name, value)
 
 \S4method{setParam}{LunParams}(object, name, value)

man/zinbEstimate.Rd

@@ -61,7 +61,9 @@ the fitted model and inserts it into a \code{\link{ZINBParams}} object. See
 \code{\link[zinbwave]{zinbFit}} for details of the estimation procedure.
 }
 \examples{
+\dontrun{
 data("sc_example_counts")
 params <- zinbEstimate(sc_example_counts)
 params
 }
+}

tests/testthat/test-BASiCSParams.R

+context("BASiCSParams")
+
+test_that("gene.params checks work", {
+    params <- newBASiCSParams()
+    expect_error(setParam(params, "gene.params", data.frame(A = 1, B = 1)),
+                 "gene.params: Incorrect column names")
+    expect_error(setParam(params, "gene.params",
+                          data.frame(Mean = 1, Disp = "a")),
+                 "gene.params: May only contain the following types: numeric")
+})
+
+test_that("cell.params checks work", {
+    params <- newBASiCSParams()
+    expect_error(setParam(params, "cell.params", data.frame(A = 1, B = 1)),
+                 "cell.params: Incorrect column names")
+    expect_error(setParam(params, "cell.params",
+                          data.frame(Phi = 1, S = "a")),
+                 "cell.params: May only contain the following types: numeric")
+})
+
+test_that("nBatches checks work", {
+    params <- newBASiCSParams()
+    expect_error(setParam(params, "nCells", 1),
+                 "nCells cannot be set directly, set batchCells instead")
+    expect_error(setParam(params, "nBatches", 1),
+                 "nBatches cannot be set directly, set batchCells instead")
+})