Skip to content
Snippets Groups Projects
Commit 9bdeca50 authored by Luke Zappia's avatar Luke Zappia
Browse files

Merge branch 'scDD-est'

* scDD-est:
  Add counts check to scDDEstimate
  Adjust nEE/nEP estimation
  Update scDD estimation
parents d3ca4008 22ba7c63
No related branches found
No related tags found
No related merge requests found
......@@ -405,24 +405,24 @@ setClass("Lun2Params",
#' The SCDD simulation uses the following parameters:
#'
#' \describe{
#' \item{\code{[nGenes]}}{The number of genes to simulate (not used).}
#' \item{\code{nGenes}}{The number of genes to simulate (not used).}
#' \item{\code{nCells}}{The number of cells to simulate in each condition.}
#' \item{\code{[seed]}}{Seed to use for generating random numbers.}
#' \item{\code{SCdat}}{
#' \code{\link[SummarizedExperiment]{SummarizedExperiment}} containing real
#' data.}
#' \item{\code{[nDE]}}{Number of DE genes to simulate.}
#' \item{\code{[nDP]}}{Number of DP genes to simulate.}
#' \item{\code{[nDM]}}{Number of DM genes to simulate.}
#' \item{\code{[nDB]}}{Number of DB genes to simulate.}
#' \item{\code{[nEE]}}{Number of EE genes to simulate.}
#' \item{\code{[nEP]}}{Number of EP genes to simulate.}
#' \item{\code{nDE}}{Number of DE genes to simulate.}
#' \item{\code{nDP}}{Number of DP genes to simulate.}
#' \item{\code{nDM}}{Number of DM genes to simulate.}
#' \item{\code{nDB}}{Number of DB genes to simulate.}
#' \item{\code{nEE}}{Number of EE genes to simulate.}
#' \item{\code{nEP}}{Number of EP genes to simulate.}
#' \item{\code{[sd.range]}}{Interval for fold change standard deviations.}
#' \item{\code{[modeFC]}}{Values for DP, DM and DB mode fold changes.}
#' \item{\code{[varInflation]}}{Variance inflation factors for each
#' condition. If all equal to 1 will be set to \code{NULL} (default)}
#' condition. If all equal to 1 will be set to \code{NULL} (default).}
#' \item{\code{[condition]}}{String giving the column that represents
#' biological group of interest}
#' biological group of interest.}
#' }
#'
#' The parameters not shown in brackets can be estimated from real data using
......
......@@ -86,12 +86,12 @@ setMethod("setParam", "SCDDParams", function(object, name, value) {
setMethod("show", "SCDDParams", function(object) {
pp <- list("Genes:" = c("[nDE]" = "nDE",
"[nDP]" = "nDP",
"[nDM]" = "nDM",
"[nDP]" = "nDP",
"[nEE]" = "nEE",
"[nEP]" = "nEP"),
pp <- list("Genes:" = c("(nDE)" = "nDE",
"(nDP)" = "nDP",
"(nDM)" = "nDM",
"(nDP)" = "nDP",
"(nEE)" = "nEE",
"(nEP)" = "nEP"),
"Fold change:" = c("[SD Range]" = "sd.range",
"[Mode FC]" = "modeFC"),
"Variance:" = c("[Inflation]" = "varInflation"),
......
......@@ -7,11 +7,15 @@
#' @param conditions Vector giving the condition that each cell belongs to.
#' Conditions can be 1 or 2.
#' @param params SCDDParams object to store estimated values in.
#' @param BPPARAM A \code{\link[BiocParallel]{BiocParallelParam}} instance
#' giving the parallel back-end to be used. Default is
#' \code{\link[BiocParallel]{SerialParam}} which uses a single core.
#'
#' @details
#' This function is just a wrapper around \code{\link[scDD]{preprocess}} that
#' takes the output and converts it to a SCDDParams object. See
#' \code{\link[scDD]{preprocess}} for details.
#' This function applies \code{\link[scDD]{preprocess}} to the counts then uses
#' \code{\link[scDD]{scDD}} to estimate the numbers of each gene type to
#' simulate. The output is then converted to a SCDDParams object. See
#' \code{\link[scDD]{preprocess}} and \code{\link[scDD]{scDD}} for details.
#'
#' @return SCDDParams object containing the estimated parameters.
#'
......@@ -20,14 +24,17 @@
#' conditions <- sample(1:2, ncol(sc_example_counts), replace = TRUE)
#' params <- scDDEstimate(sc_example_counts, conditions)
#' params
#' @importFrom BiocParallel SerialParam
#' @export
scDDEstimate <- function(counts, conditions, params = newSCDDParams()) {
scDDEstimate <- function(counts, conditions, params = newSCDDParams(),
BPPARAM = SerialParam()) {
UseMethod("scDDEstimate")
}
#' @rdname scDDEstimate
#' @export
scDDEstimate.SCESet <- function(counts, conditions, params = newSCDDParams()) {
scDDEstimate.SCESet <- function(counts, conditions, params = newSCDDParams(),
BPPARAM = SerialParam()) {
counts <- scater::counts(counts)
scDDEstimate(counts, conditions, params)
}
......@@ -35,13 +42,17 @@ scDDEstimate.SCESet <- function(counts, conditions, params = newSCDDParams()) {
#' @rdname scDDEstimate
#' @importFrom methods as
#' @export
scDDEstimate.matrix <- function(counts, conditions, params = newSCDDParams()) {
scDDEstimate.matrix <- function(counts, conditions, params = newSCDDParams(),
BPPARAM = SerialParam()) {
if (!requireNamespace("scDD", quietly = TRUE)) {
stop("The scDD simulation requires the 'scDD' package.")
}
checkmate::assertClass(params, "SCDDParams")
checkmate::assertMatrix(counts, mode = "numeric", any.missing = FALSE,
min.rows = 1, min.cols = 1, row.names = "unique",
col.names = "unique")
checkmate::assertIntegerish(conditions, len = ncol(counts), lower = 1,
upper = 2)
......@@ -59,8 +70,29 @@ scDDEstimate.matrix <- function(counts, conditions, params = newSCDDParams()) {
SCdat <- SummarizedExperiment::SummarizedExperiment(assays = assays,
colData = colData)
params <- setParams(params, nCells = round(dim(SCdat)[2] / 2),
SCdat = SCdat)
SCdat <- scDD::scDD(SCdat, testZeroes = FALSE, param = BPPARAM)
res <- scDD::results(SCdat)
res <- res[!is.na(res$DDcategory), ]
dd.cats <- table(res$DDcategory)
not.dd <- res$DDcategory == "NS"
nDE <- ifelse("DE" %in% names(dd.cats), dd.cats["DE"], 0)
nDP <- ifelse("DP" %in% names(dd.cats), dd.cats["DP"], 0)
nDM <- ifelse("DM" %in% names(dd.cats), dd.cats["DM"], 0)
nDB <- ifelse("DB" %in% names(dd.cats), dd.cats["DB"], 0)
nEP <- sum(res$Clusters.c1[not.dd] > 1 & res$Clusters.c2[not.dd] > 1)
nEE <- nrow(counts) - nDE - nDP - nDM - nDB - nEP
params <- setParams(params,
nCells = round(dim(SCdat)[2] / 2),
SCdat = SCdat,
nDE = nDE,
nDP = nDP,
nDM = nDM,
nDB = nDB,
nEE = nEE,
nEP = nEP)
return(params)
}
\ No newline at end of file
}
......@@ -14,24 +14,24 @@ S4 class that holds parameters for the scDD simulation.
The SCDD simulation uses the following parameters:
\describe{
\item{\code{[nGenes]}}{The number of genes to simulate (not used).}
\item{\code{nGenes}}{The number of genes to simulate (not used).}
\item{\code{nCells}}{The number of cells to simulate in each condition.}
\item{\code{[seed]}}{Seed to use for generating random numbers.}
\item{\code{SCdat}}{
\code{\link[SummarizedExperiment]{SummarizedExperiment}} containing real
data.}
\item{\code{[nDE]}}{Number of DE genes to simulate.}
\item{\code{[nDP]}}{Number of DP genes to simulate.}
\item{\code{[nDM]}}{Number of DM genes to simulate.}
\item{\code{[nDB]}}{Number of DB genes to simulate.}
\item{\code{[nEE]}}{Number of EE genes to simulate.}
\item{\code{[nEP]}}{Number of EP genes to simulate.}
\item{\code{nDE}}{Number of DE genes to simulate.}
\item{\code{nDP}}{Number of DP genes to simulate.}
\item{\code{nDM}}{Number of DM genes to simulate.}
\item{\code{nDB}}{Number of DB genes to simulate.}
\item{\code{nEE}}{Number of EE genes to simulate.}
\item{\code{nEP}}{Number of EP genes to simulate.}
\item{\code{[sd.range]}}{Interval for fold change standard deviations.}
\item{\code{[modeFC]}}{Values for DP, DM and DB mode fold changes.}
\item{\code{[varInflation]}}{Variance inflation factors for each
condition. If all equal to 1 will be set to \code{NULL} (default)}
condition. If all equal to 1 will be set to \code{NULL} (default).}
\item{\code{[condition]}}{String giving the column that represents
biological group of interest}
biological group of interest.}
}
The parameters not shown in brackets can be estimated from real data using
......
......@@ -6,11 +6,14 @@
\alias{scDDEstimate.matrix}
\title{Estimate scDD simulation parameters}
\usage{
scDDEstimate(counts, conditions, params = newSCDDParams())
scDDEstimate(counts, conditions, params = newSCDDParams(),
BPPARAM = SerialParam())
\method{scDDEstimate}{SCESet}(counts, conditions, params = newSCDDParams())
\method{scDDEstimate}{SCESet}(counts, conditions, params = newSCDDParams(),
BPPARAM = SerialParam())
\method{scDDEstimate}{matrix}(counts, conditions, params = newSCDDParams())
\method{scDDEstimate}{matrix}(counts, conditions, params = newSCDDParams(),
BPPARAM = SerialParam())
}
\arguments{
\item{counts}{either a counts matrix or an SCESet object containing count
......@@ -20,6 +23,10 @@ data to estimate parameters from.}
Conditions can be 1 or 2.}
\item{params}{SCDDParams object to store estimated values in.}
\item{BPPARAM}{A \code{\link[BiocParallel]{BiocParallelParam}} instance
giving the parallel back-end to be used. Default is
\code{\link[BiocParallel]{SerialParam}} which uses a single core.}
}
\value{
SCDDParams object containing the estimated parameters.
......@@ -28,9 +35,10 @@ SCDDParams object containing the estimated parameters.
Estimate simulation parameters for the scDD simulation from a real dataset.
}
\details{
This function is just a wrapper around \code{\link[scDD]{preprocess}} that
takes the output and converts it to a SCDDParams object. See
\code{\link[scDD]{preprocess}} for details.
This function applies \code{\link[scDD]{preprocess}} to the counts then uses
\code{\link[scDD]{scDD}} to estimate the numbers of each gene type to
simulate. The output is then converted to a SCDDParams object. See
\code{\link[scDD]{preprocess}} and \code{\link[scDD]{scDD}} for details.
}
\examples{
data("sc_example_counts")
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment