### GenericCompassWrapper.R #######################################################
#' Run COMPASS once
#'
#' This function runs COMPASS once, saving output to disk.
#' Intended for use by the generic.compass.wrapper function below.
#' @param gs The GatingSet of GatingSetList
#' @param cnode Node on which to run COMPASS
#' @param nodemarkermap List mapping nodes to marker names
#' @param outdir Directory in which to save output, e.g. heatmaps
#' @param individuals pData column containing individual identifiers (rows of heatmap)
#' @param seed (Optional) Number to set seed to [default NULL]
#' @param grouping (Optional) pData columns on which to group rows in heatmap, as a character vector [default NULL]
#' @param uniqueruns (Optional) pData column identifying unique runs. Use if you need multiple runs. [default NULL]
#' @param lineplotxvar (Optional) pData column which defines groups along x-axis in FS-score line plot, e.g. Time [default NULL]
#' @param iter (Optional) Number of COMPASS iterations to perform on each repitition (8 repetitions total) [default 40,000]
#' @param lineplotgroupby (Optional) This should be specified if lineplotxvar is given. pData column which defines which values to connect in the line plot (usually something like "PTID")
#' @return Nothing
#' @keywords COMPASS
#' @import grid
#' @export
run.compass.once <- function(gs,
cnode,
individuals,
nodemarkermap,
iter=40000,
lineplotxvar=NULL,
run=NULL,
outdir,
uniqueruns=NULL,
grouping=NULL,
lineplotgroupby=NULL,
countFilterThreshold=0,
seed=NULL) {
if (!is.null(seed)) {
set.seed(seed)
}
# Create a COMPASSContainer from the GatingSetList.
CC <- COMPASS::COMPASSContainerFromGatingSet(gs, node=cnode, individual_id=individuals,
mp=nodemarkermap, countFilterThreshold=countFilterThreshold)
# Run COMPASS for this unique run
fit <- COMPASS::COMPASS( CC,
treatment=trt == "Treatment",
control=trt == "Control",
iterations=iter
)
FS <- COMPASS::FunctionalityScore(fit)
PFS <- COMPASS::PolyfunctionalityScore(fit)
# Initialize subdirectory to save output for this run
fileSuffix <- if (is.null(run)) { cnode } else { paste(cnode, run, sep="_") }
subDir <- fileSuffix
subDirPath <- file.path(outdir, subDir)
dir.create(subDirPath, showWarnings = FALSE)
setwd(subDirPath)
# Save the COMPASS run as an RDS file for safekeeping
saveRDS(fit, paste(c("COMPASSResult_", fileSuffix, ".rds"), collapse=""))
# Initialize output text files for FS and PFS
columnsFormat <- paste("CellSubset", uniqueruns, sep="\t")
fsFile <- paste(paste("FS", fileSuffix, sep="_"), ".txt", sep="")
pfsFile <- paste("P", fsFile, sep="")
write(paste(columnsFormat, paste(names(FS), collapse="\t"), sep="\t"), file=fsFile, append=TRUE)
write(paste(columnsFormat, paste(names(PFS), collapse="\t"), sep="\t"), file=pfsFile, append=TRUE)
# Write results to file, one for each statistic type
write(paste(paste(cnode, run, sep="\t"), paste(FS, collapse="\t"), sep="\t"), file=fsFile, append=TRUE)
write(paste(paste(cnode, run, sep="\t"), paste(PFS, collapse="\t"), sep="\t"), file=pfsFile, append=TRUE)
plotTitleSuffix <- paste(c(",\n", run, ", ", cnode, " Cells"), collapse="")
cytokine_annotation_colors <- c("black", "black", "black", "black", "black", "black", "black")
# Plot a heatmap of the mean probability of response, to visualize differences
# in expression for each category
png(filename=paste(c("HeatmapMeanProbResponse", "_", cnode, run, ".png"), collapse=""),
width=800, height=650)
try(grid::grid.draw(print(plot(fit, grouping, show_rownames = TRUE,
main = paste("Heatmap of Mean Probability of Response", plotTitleSuffix, sep=""),
fontsize=14, fontsize_row=13, fontsize_col=11,
cytokine_annotation_colors=cytokine_annotation_colors))))
dev.off()
# Log scale of previous data, smaller changes show up better
png(filename=paste(c("HeatmapLogPostResponse", "_", cnode, run, ".png"), collapse=""),
width=800, height=650)
try(grid::grid.draw(print(plot(fit, grouping, show_rownames = TRUE,
measure=COMPASS::PosteriorLogDiff(fit), threshold=0,
main = paste("Heatmap of Log Posterior Differences in Response", plotTitleSuffix, sep=""),
fontsize=14, fontsize_row=13, fontsize_col=11,
cytokine_annotation_colors=cytokine_annotation_colors))))
dev.off()
# If applicable, create line plot of Functionality Score vs. lineplotxvar
try(if (!is.null(lineplotxvar)) {
# First format FS into a data.table
# The order of rows in metadata is not necessarily the same as those of FSdf. Perform a table merge
FSdf <- data.frame(names(FS), FS, row.names=NULL)
names(FSdf) <- c(individuals, "FunctionalityScore")
metasub <- pData(gsListForCOMPASSsub)
fsplotdf <- merge(x=metasub[metasub[,uniqueruns]==run,], y=FSdf, by.x=individuals, by.y=individuals)
# Draw the line plot
lineplot <- ggplot2::ggplot(data=fsplotdf, aes(x=factor(get(lineplotxvar)), y=FunctionalityScore, group=factor(get(lineplotgroupby)))) +
ggplot2::geom_point() + ggplot2::geom_line() +
ggplot2::labs(title=paste(c("Functionality Score vs. ", lineplotxvar, " Line Plot,\n", cnode, " ", run), collapse=""),
x=lineplotxvar) +
ggplot2::theme_set(theme_gray(base_size = 25))
ggplot2::ggsave(filename=paste(c("LinePlot_FSv", lineplotxvar, "_", cnode, run, ".png"), collapse=""),
plot=lineplot,
width=6.66, height=7.85)
})
# Call the garbage collector to free up memory
gc()
}
#' COMPASS Wrapper
#'
#' This function runs COMPASS once for each of the unique values defined by the uniqueruns argument (if provided)
#' @param path The path to the folder in which the GatingSetList or GatingSet is saved
#' @param cnode Node on which to run COMPASS
#' @param nodemarkermap List mapping nodes to marker names
#' @param outdir Directory in which to save output, e.g. heatmaps
#' @param individuals pData column containing individual identifiers (rows of heatmap)
#' @param seed (Optional) Number to set seed to [default NULL]
#' @param grouping (Optional) pData columns on which to group rows in heatmap, as a character vector [default NULL]
#' @param uniqueruns (Optional) pData column identifying unique runs. Use if you need multiple runs. [default NULL]
#' @param lineplotxvar (Optional) pData column which defines groups along x-axis in FS-score line plot, e.g. Time [default NULL]
#' @param iter (Optional) Number of COMPASS iterations to perform on each repitition (8 repetitions total) [default 40,000]
#' @param lineplotgroupby (Optional) This should be specified if lineplotxvar is given. pData column which defines which values to connect in the line plot (usually something like "PTID")
#' @return Nothing
#' @keywords COMPASS
#' @import grid
#' @export
#' @examples
#' \dontrun{
#' generic.compass.wrapper(path="/home/path/to/GatingSetList",
#' seed=1,
#' outdir="/path/to/OutDirectory",
#' cnode="8+",
#' nodemarkermap=list("8+/154+" = "CD154",
#' "8+/IFNg+" = "IFNg",
#' "8+/IL4+" = "IL4",
#' "8+/TNFa+" = "TNFa",
#' "8+/IL22+" = "IL22",
#' "8+/IL17+" = "IL17a",
#' "8+/IL2+" = "IL2"),
#' individuals="PATIENT ID",
#' uniqueruns="Peptide")
#' }
generic.compass.wrapper <- function(path=NULL,
seed=NULL,
outdir=NULL,
cnode=NULL,
nodemarkermap=NULL,
individuals=NULL,
grouping=NULL,
uniqueruns=NULL,
lineplotxvar=NULL,
iter=40000,
lineplotgroupby=NULL,
countFilterThreshold=0) {
# Assumptions: pData trt column contains "Treatment" and "Control" labels
# TODO: run in parallel
# TODO: test for single run from command line
# Check that required arguments are provided
if (is.null(path)) {
stop("Path parameter must be provided.")
}
if (is.null(outdir)) {
stop("Outdir parameter must be provided.")
}
if (is.null(nodemarkermap)) {
stop("Nodemarkermap parameter must be provided.")
}
if (is.null(individuals)) {
stop("Individuals parameter must be provided.")
}
if (is.null(cnode)) {
stop("cnode parameter must be provided.")
}
if (!is.null(lineplotxvar) & is.null(lineplotgroupby)) {
stop("Please specify lineplotgroupby parameter, which must be provided if lineplotxvar is provided.")
}
cat(paste(as.character(Sys.time()), "Loading GatingSetList\n"))
if (!is.null(seed)) {
set.seed(seed)
}
# Load the saved GatingSetList or GatingSet:
loadGSListOrGS <- function (path) {
out <- try(flowWorkspace::load_gslist(path))
if (class(out) == "try-error") {
cat("Caught an error during flowWorkspace::load_gslist, trying flowWorkspace::load_gs.\n")
out <- flowWorkspace::load_gs(path)
}
out
}
gsListForCOMPASS <- loadGSListOrGS(path)
meta <- flowWorkspace::pData(gsListForCOMPASS)
# Run COMPASS once or multiple times depending on whether uniqueruns is given
if (is.null(uniqueruns)) {
cat(paste(c(as.character(Sys.time()), " Running COMPASS for", cnode, "cells", "...\n"), collapse=" "))
try(run.compass.once(gs=gsListForCOMPASS,
cnode=cnode,
individuals=individuals,
nodemarkermap=nodemarkermap,
iter=iter,
lineplotxvar=lineplotxvar,
run=NULL,
outdir=outdir,
uniqueruns=NULL,
grouping=grouping,
countFilterThreshold=countFilterThreshold))
} else {
# Get a list of values identifying unique runs from the "uniqueruns" column, minus the control value
uniqueRunsList <- unique(meta[meta[,"trt"]!="Control",][,uniqueruns])
# Obtain the value in the "uniqueruns" column that corresponds to the control value
controlval <- unique(meta[meta[,"trt"]=="Control",][,uniqueruns])
for (run in uniqueRunsList) {
cat(paste(c(as.character(Sys.time()), " Running COMPASS for", cnode, "cells ", run, "...\n"), collapse=" "))
# Subset the GatingSetList for the desired unique run
rowSubsetBooleans <- meta[,uniqueruns] %in% c(run, controlval)
rownamesAsChar <- as.character(rownames(meta[rowSubsetBooleans,]))
gsListForCOMPASSsub <- gsListForCOMPASS[rownamesAsChar]
try(run.compass.once(gs=gsListForCOMPASSsub,
cnode=cnode,
individuals=individuals,
nodemarkermap=nodemarkermap,
iter=iter,
lineplotxvar=lineplotxvar,
run=run,
outdir=outdir,
uniqueruns=uniqueruns,
grouping=grouping,
lineplotgroupby=lineplotgroupby,
countFilterThreshold=countFilterThreshold))
}
}
cat(paste(as.character(Sys.time()), " All COMPASS runs Done\n"))
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.