#' Partition de.lorean object by cells
#'
#' @param dl de.lorean object
#' @param pieces How many pieces to partition into
#'
partition.de.lorean <- function(
dl,
pieces = 2
) {
partition <- (1:(dim(dl)[2]) %% pieces) + 1
cells <- sample(colnames(dl$expr))
get.piece <- function(p) {
partition.cells <- cells[partition == p]
cell.filter <- function(cells) cells %in% partition.cells
filter_cells(dl, cell.filter)
}
lapply(1:pieces, get.piece)
}
#' Test robustness of pseudotime estimation on subsets of de.lorean
#' object
#'
#' @param dl de.lorean object
#' @param pieces How many pieces to partition into
#'
test.robustness.de.lorean <- function(
dl,
pieces = 2
) {
# Partition de.lorean object into several pieces
partition <- partition.de.lorean(dl, pieces)
# Define function to fit each piece
run.model <- function(piece) {
piece <- prepare.for.stan(piece)
piece <- find.best.tau(piece)
piece <- fit.model(piece)
process.posterior(piece)
}
# Fit full de.lorean object
full.model <- run.model(dl)
# Get tau posterior
full.tau <- (full.model$samples.l$tau
%>% dplyr::select(-c)
%>% mutate(fit="part"))
# Fit each piece
piece.models <- lapply(partition, run.model)
# Get tau posterior from each piece and adjust mean to match full fit
# mean
get.tau.posterior <- function(piece.model) {
posterior <- (piece.model$samples.l$tau
%>% dplyr::select(-c)
%>% mutate(fit="full"))
cells <- unique(posterior$cell)
full.mean <- mean((full.tau %>% filter(cell %in% cells))$tau)
piece.mean <- mean(posterior$tau)
posterior$tau <- posterior$tau - piece.mean + full.mean
posterior
}
# Gather tau posteriors from each piece
pieces.tau <- do.call(rbind,
lapply(piece.models, get.tau.posterior))
# Combine fit pieces with fit full model
tau.posterior <- rbind(full.tau, pieces.tau)
# Sort by median tau
cells <- (tau.posterior
%>% group_by(cell)
%>% dplyr::summarise(median.tau=median(tau))
%>% arrange(median.tau)
)
# Create result list with plot
list(tau.posterior = tau.posterior,
gp = (ggplot(tau.posterior
%>% mutate(cell=factor(cell, levels=cells$cell)),
aes(x=cell, y=tau, fill=fit))
+ geom_boxplot()
))
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.