mmDS {muscat} | R Documentation |
Performs cluster-wise DE analysis by fitting cell-level models.
mmDS( x, coef = NULL, covs = NULL, method = c("dream", "vst", "poisson", "nbinom", "hybrid"), n_cells = 10, n_samples = 2, min_count = 1, min_cells = 20, n_threads = 8, verbose = TRUE, dup_corr = FALSE, trended = FALSE, vst = c("sctransform", "DESeq2"), bayesian = FALSE, blind = TRUE, REML = TRUE, ddf = c("Satterthwaite", "Kenward-Roger", "lme4") ) .mm_dream( x, coef = NULL, covs = NULL, dup_corr = FALSE, trended = FALSE, ddf = c("Satterthwaite", "Kenward-Roger"), n_threads = 1, verbose = FALSE ) .mm_vst( x, vst = c("sctransform", "DESeq2"), coef = NULL, covs = NULL, bayesian = FALSE, blind = TRUE, REML = TRUE, ddf = c("Satterthwaite", "Kenward-Roger", "lme4"), n_threads = 1, verbose = FALSE )
x |
|
coef |
character specifying the coefficient to test.
If NULL (default), will test the last level of |
covs |
character vector of |
method |
a character string.
Either |
n_cells |
number of cells per cluster-sample required to consider a sample for testing. |
n_samples |
number of samples per group required to consider a cluster for testing. |
min_count |
numeric. For a gene to be tested in a given cluster,
at least |
min_cells |
number (or fraction, if < 1) of cells with a count >
|
n_threads |
number of threads to use. |
verbose |
logical specifying whether messages on progress and a progress bar should be displayed. |
dup_corr |
logical; whether to use
|
trended |
logical; whether to use expression-dependent variance priors
in |
vst |
method to use as variance-stabilizing transformations.
|
bayesian |
logical; whether to use bayesian mixed models. |
blind |
logical; whether to ignore experimental design for the vst. |
REML |
logical; whether to maximize REML instead of log-likelihood. |
ddf |
character string specifying the method for estimating
the effective degrees of freedom. For |
.mm_dream
and .mm_vst
expect cells from a single cluster,
and do not perform filtering or handle incorrect parameters well.
Meant to be called by mmDS
with method = c("dream", "vst")
and
vst = c("sctransform", "DESeq2")
to be applied across all clusters.
method = "dream"
variancePartition
's voom-lme4-implementation
of mixed models for RNA-seq data; function dream
.
method = "vst"
vst = "sctransform"
lmer
or blmer
mixed models on
vst
transformed counts.
vst = "DESeq2"
varianceStabilizingTransformation
followed by lme4
mixed models.
a data.frame
.mm_dream
: see details.
.mm_vst
: see details.
Pierre-Luc Germain & Helena L Crowell
Crowell, HL, Soneson, C, Germain, P-L, Calini, D, Collin, L, Raposo, C, Malhotra, D & Robinson, MD: On the discovery of population-specific state transitions from multi-sample multi-condition single-cell RNA sequencing data. bioRxiv 713412 (2018). doi: https://doi.org/10.1101/713412
data(sce) # subset "B cells" cluster sce <- sce[, sce$cluster_id == "B cells"] sce$cluster_id <- droplevels(sce$cluster_id) # downsample to 100 genes cs_by_s <- split(colnames(sce), sce$sample_id) gs <- sample(nrow(sce), 100) sce <- sce[gs, ] res <- mmDS(sce, method = "dream", n_threads = 1, verbose = FALSE) head(res$`B cells`)