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

Add estimateSimpleParams function

parent 89c3624c
No related branches found
No related tags found
No related merge requests found
...@@ -32,7 +32,7 @@ setMethod("show", "SimpleParams", function(object) { ...@@ -32,7 +32,7 @@ setMethod("show", "SimpleParams", function(object) {
pp <- list("Mean:" = c("(Rate)" = "mean.rate", pp <- list("Mean:" = c("(Rate)" = "mean.rate",
"(Shape)" = "mean.shape"), "(Shape)" = "mean.shape"),
"Counts:" = c("(Dispersion)" = "count.disp")) "Counts:" = c("[Dispersion]" = "count.disp"))
callNextMethod() callNextMethod()
showPP(object, pp) showPP(object, pp)
......
estimateSimpleParams <- function(data, params = newSimpleParams()) {
UseMethod("estimateSimpleParams")
}
#' @rdname simpleEstimate
#' @export
estimateSimpleParams.SCESet <- function(data, params = newSimpleParams()) {
counts <- scater::counts(x)
estimateSimpleParams(counts, params)
}
#' @rdname simpleEstimate
#' @importFrom stats median
#' @export
estimateSimpleParams.matrix <- function(data, params = newSimpleParams()) {
checkmate::assertClass(params, "SimpleParams")
# Normalise for library size and remove all zero genes
lib.sizes <- colSums(data)
lib.med <- median(lib.sizes)
norm.counts <- t(t(data) / lib.sizes * lib.med)
norm.counts <- norm.counts[rowSums(norm.counts > 0) > 1, ]
means <- rowMeans(norm.counts)
means.fit <- fitdistrplus::fitdist(means, "gamma", method = "mme")
params <- setParams(params, nGenes = nrow(data), nCells = nrow(data),
mean.shape = unname(means.fit$estimate["shape"]),
mean.rate = unname(means.fit$estimate["rate"]))
return(params)
}
\ No newline at end of file
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