-
Luke Zappia authoredLuke Zappia authored
splatSimulate.Rd 4.75 KiB
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/splat-simulate.R
\name{splatSimulate}
\alias{splatSimulate}
\alias{splatSimulateGroups}
\alias{splatSimulatePaths}
\alias{splatSimulateSingle}
\title{Splatter simulation}
\usage{
splatSimulate(params = newSplatParams(), method = c("single", "groups",
"paths"), verbose = TRUE, ...)
splatSimulateSingle(params = newSplatParams(), verbose = TRUE, ...)
splatSimulateGroups(params = newSplatParams(), verbose = TRUE, ...)
splatSimulatePaths(params = newSplatParams(), verbose = TRUE, ...)
}
\arguments{
\item{params}{SplatParams object containing parameters for the simulation.
See \code{\link{SplatParams}} for details.}
\item{method}{which simulation method to use. Options are "single" which
produces a single population, "groups" which produces distinct groups
(eg. cell types) or "paths" which selects cells from continuous
trajectories (eg. differentiation processes).}
\item{verbose}{logical. Whether to print progress messages.}
\item{...}{any additional parameter settings to override what is provided in
\code{params}.}
}
\value{
SCESet object containing the simulated counts and intermediate
values.
}
\description{
Simulate count data from a fictional single-cell RNA-seq experiment using
the Splatter method.
}
\details{
Parameters can be set in a variety of ways. If no parameters are provided
the default parameters are used. Any parameters in \code{params} can be
overridden by supplying additional arguments through a call to
\code{\link{setParams}}. This design allows the user flexibility in
how they supply parameters and allows small adjustments without creating a
new \code{SplatParams} object. See examples for a demonstration of how this
can be used.
The simulation involves the following steps:
\enumerate{
\item Set up simulation object
\item Simulate library sizes
\item Simulate gene means
\item Simulate groups/paths
\item Simulate BCV adjusted cell means
\item Simulate true counts
\item Simulate dropout
\item Create final SCESet object
}
The final output is an \code{\link[scater]{SCESet}} object that contains the
simulated counts but also the values for various intermediate steps. These
are stored in the \code{\link[Biobase]{phenoData}} (for cell specific
information), \code{\link[Biobase]{featureData}} (for gene specific
information) or \code{\link[Biobase]{assayData}} (for gene by cell matrices)
slots. This additional information includes:
\describe{
\item{\code{phenoData}}{
\describe{
\item{Cell}{Unique cell identifier.}
\item{Group}{The group or path the cell belongs to.}
\item{ExpLibSize}{The expected library size for that cell.}
\item{Step (paths only)}{how far along the path each cell is.}
}
}
\item{\code{featureData}}{
\describe{
\item{Gene}{Unique gene identifier.}
\item{BaseGeneMean}{The base expression level for that gene.}
\item{OutlierFactor}{Expression outlier factor for that gene. Values
of 1 indicate the gene is not an expression outlier.}
\item{GeneMean}{Expression level after applying outlier factors.}
\item{DEFac[Group]}{The differential expression factor for each gene
in a particular group. Values of 1 indicate the gene is not
differentially expressed.}
\item{GeneMean[Group]}{Expression level of a gene in a particular
group after applying differential expression factors.}
}
}
\item{\code{assayData}}{
\describe{
\item{BaseCellMeans}{The expression of genes in each cell adjusted for
expected library size.}
\item{BCV}{The Biological Coefficient of Variation for each gene in
each cell.}
\item{CellMeans}{The expression level of genes in each cell adjusted
for BCV.}
\item{TrueCounts}{The simulated counts before dropout.}
\item{Dropout}{Logical matrix showing which values have been dropped
in which cells.}
}
}
}
Values that have been added by Splatter are named using \code{CamelCase} in
order to differentiate them from the values added by Scater which uses
\code{underscore_naming}.
}
\examples{
# Simulation with default parameters
sim <- splatSimulate()
# Simulation with different number of genes
sim <- splatSimulate(nGenes = 1000)
# Simulation with custom parameters
params <- newSplatParams(nGenes = 100, mean.rate = 0.5)
sim <- splatSimulate(params)
# Simulation with adjusted custom parameters
sim <- splatSimulate(params, mean.rate = 0.6, out.prob = 0.2)
# Simulate groups
sim <- splatSimulate(method = "groups")
# Simulate paths
sim <- splatSimulate(method = "paths")
}
\seealso{
\code{\link{splatSimLibSizes}}, \code{\link{splatSimGeneMeans}},
\code{\link{splatSimDE}}, \code{\link{splatSimCellMeans}},
\code{\link{splatSimBCVMeans}}, \code{\link{splatSimTrueCounts}},
\code{\link{splatSimDropout}}
}