## CYTO_SPILLOVER_SPREAD_COMPUTE -----------------------------------------------
#' Compute Spillover Spreading Matrix
#'
#' \code{cyto_spillover_spread_compute} computes the spillover spreading matrix
#' as described by Nguyen et al. (2013) using gated and compensated compensation
#' controls.
#'
#' @param x object of class \code{\link[flowCore:flowSet-class]{flowSet}} or
#' \code{\link[flowWorkspace:GatingSet-class]{GatingSet}} containing
#' compensated, transformed and gated compensation controls.
#' @param parent name of the population to use for downstream calculations when
#' a \code{GatingSet} object is supplied, set to the last gated population by
#' default (e.g. "Single Cells").
#' @param axes_trans object of class
#' \code{\link[flowWorkspace:transformerList]{transformerList}} generated by a
#' \code{cyto_transformer} which contains the transformer definitions to be
#' used internally to transform the fluorescent channels of the samples as
#' required.
#' @param channel_match name of csv file to associate a fluorescent channel with
#' each of the compensation controls. The \code{channel_match} file should
#' contain a "name" column with the names of the compensation controls and a
#' "channel" column to associate a fluorescent channel with each compensation
#' control. Users need not generate this file by hand as it will be created
#' following the channel selection process.
#' @param compensated logical indicating whether the supplied compensation
#' controls have been compensated. Compensation controls must be compensated
#' prior to calculation of spillover spreading matrix. If \code{compensated}
#' is set to FALSE, the compensation controls will be compensated internally
#' using the supplied \code{spillover} matrix.
#' @param spillover name of the output spillover matrix csv file to be used
#' internally to compensate the compensation controls. If no spillover matrix
#' is supplied \code{spillover_spread_compute}, the spillover matrix attached
#' to compensation controls will be applied.
#' @param spillover_spread name of the csv file to which the spillover spreading
#' matrix will be saved, set to \code{"date-Spillover-Spread-Matrix.csv"} by
#' default.
#' @param axes_limits options include \code{"auto"}, \code{"data"} or
#' \code{"machine"} to use optimised, data or machine limits respectively. Set
#' to \code{"machine"} by default to use entire axes ranges. Fine control over
#' axes limits can be obtained by altering the \code{xlim} and \code{ylim}
#' arguments.
#' @param ... additional arguments passed to \code{\link{cyto_plot}}.
#'
#' @references Nguyen R, Perfetto S, Mahnke YD, Chattopadhyay P & Roederer M,
#' (2013). Quantifying spillover spreading for comparing instrument
#' performance and aiding in multicolor panel design. Cytometry A,
#' 83(3):306-15.
#'
#' @importFrom flowWorkspace pData
#' @importFrom flowCore Subset parameters exprs
#' @importFrom stats quantile na.omit
#' @importFrom grDevices graphics.off dev.new
#' @importFrom methods is
#'
#' @author Dillon Hammill, \email{Dillon.Hammill@anu.edu.au}
#'
#' @seealso \code{\link{cyto_spillover_compute}}
#' @seealso \code{\link{cyto_spillover_edit}}
#'
#' @export
cyto_spillover_spread_compute <- function(x,
parent = NULL,
axes_trans = NULL,
channel_match = NULL,
compensated = FALSE,
spillover = NULL,
spillover_spread = NULL,
axes_limits = "machine",
...){
# PREPARE DATA ----------------------------------------------------------
# COPY
cyto_copy <- cyto_copy(x)
# CHANNELS
channels <- cyto_fluor_channels(cyto_copy)
# TRANSFORMATIONS
axes_trans <- cyto_transformer_extract(cyto_copy)
# GATINGSET COMPENSATION
if(is(cyto_copy, "GatingSet")){
spill <- cyto_spillover_extract(cyto_copy)
if(length(spill) == 0){
compensated <- FALSE
}else{
compensated <- TRUE
}
}
# DATA SHOULD BE COMPENSATED & TRANSFORMED
if(compensated == FALSE){
# UNTRANSFORMED DATA
if(.all_na(axes_trans) | is.null(axes_trans)){
cyto_copy <- cyto_compensate(cyto_copy, spillover = spillover)
axes_trans <- cyto_transformer_biex(cyto_copy,
channels = channels,
type = "biex",
plot = FALSE)
cyto_copy <- suppressWarnings(cyto_transform(cyto_copy,
trans = axes_trans,
plot = FALSE))
# TRANSFORMED DATA
}else{
# INVERSE TRANSFORM
suppressWarnings(cyto_transform(cyto_extract(cyto_copy),
trans = axes_trans,
inverse = TRUE,
plot = FALSE))
# COMPENSATE
cyto_compensate(cyto_copy, spillover = spillover)
# TRANSFORM
suppressWarnings(cyto_transform(cyto_extract(cyto_copy),
trans = axes_trans,
plot = FALSE))
}
}else if(compensated == TRUE){
# TRANSFORMATIONS
if(.all_na(axes_trans) | is.null(axes_trans)){
axes_trans <- cyto_transformer_biex(cyto_copy,
channels = channels,
type = "biex",
plot = FALSE)
suppressWarnings(cyto_transform(cyto_copy,
trans = axes_trans,
plot = FALSE))
}
}
## SPILLOVER SPREAD COMPUTATION ----------------------------------------------
# NEGATIVE AND POSITIVE POPULATIONS (TRANSFORMED)
pops <- .cyto_spillover_pops(cyto_copy,
parent = parent,
channel_match = channel_match,
axes_trans = axes_trans,
axes_limits = axes_limits,
...)
neg_pops <- pops[["negative"]]
pos_pops <- pops[["positive"]]
# UPDATE DETAILS IN X
cyto_details(x) <- cyto_details(cyto_copy)
# EXPERIMENT DETAILS
pd <- cyto_details(pos_pops)
# SAMPLE NAMES
nms <- cyto_names(pos_pops)
# INVERSE TRANSFORMATIONS
neg_pops <- cyto_transform(neg_pops,
trans = axes_trans,
inverse = TRUE,
plot = FALSE)
pos_pops <- cyto_transform(pos_pops,
trans = axes_trans,
inverse = TRUE,
plot = FALSE)
# Calculate spillover spread for each compensation control
SSM <- lapply(seq_along(neg_pops), function(z){
# Extract NIL and pop
NIL <- neg_pops[[z]]
pop <- pos_pops[[z]]
# Channel
chan <- pd$channel[pd$name == cyto_names(pop)]
# Calculate 50th and 84th percentile for NIL
NIL_50 <- LAPPLY(channels, function(channel){
quantile(exprs(NIL)[,channel], 0.5)
})
names(NIL_50) <- channels
NIL_84 <- LAPPLY(channels, function(channel){
quantile(exprs(NIL)[,channel], 0.84)
})
names(NIL_84) <- channels
# Combine into a list - NIL_stats
NIL_stats <- list("50%" = NIL_50, "84%" = NIL_84)
# Calculate 50th and 84th percentile for pop
pop_50 <- LAPPLY(channels, function(channel){
quantile(exprs(pop)[,channel], 0.5)
})
names(pop_50) <- channels
pop_84 <- LAPPLY(channels, function(channel){
quantile(exprs(pop)[,channel], 0.84)
})
names(pop_84) <- channels
# Combine into a list - pop_stats
pop_stats <- list("50%" = pop_50, "84%" = pop_84)
# Calculate pop SD and VARIANCE
pop_SD <- pop_stats[[2]] - pop_stats[[1]]
pop_VAR <- pop_SD^2
# Calculate NIL SD and VARIANCE
NIL_SD <- NIL_stats[[2]] - NIL_stats[[1]]
NIL_VAR <- NIL_SD^2
# Difference in VAR of NIL and pop in other detectors
VAR_diff <- pop_VAR - NIL_VAR
# Replace negative numbers with 0 - stain has less spread than reference
if(any(VAR_diff < 0)){
VAR_diff[VAR_diff < 0] <- 0
}
# Difference in SD of NIL and pop in other detectors
SD_diff <- sqrt(VAR_diff)
# Change in fluorescence in primary detector
signal <- pop_stats[[1]][match(chan, channels)] -
NIL_stats[[1]][match(chan, channels)]
# Spillover spread values
return(SD_diff/sqrt(signal))
})
# Make matrix
SSM <- do.call("rbind", SSM)
# Add rownames (channel associated with each control)
rownames(SSM) <- pd$channel[match(cyto_names(pos_pops), pd$name)]
# Replace diagonal with NA
lapply(seq(1, nrow(SSM), 1), function(z) {
SSM[z, match(rownames(SSM)[z], colnames(SSM))] <<- NA
})
# Sort SSM by column names
ind <- na.omit(match(colnames(SSM), rownames(SSM)))
SSM <- SSM[ind, ]
# Turn off pop-up graphics device (gating)
graphics.off()
dev.new()
# Heatmap label size
ylab_width <- max(nchar(rownames(SSM))) * 0.02
xlab_width <- max(nchar(colnames(SSM))) * 0.04
# Plot heatmap
superheat::superheat(SSM,
title = "Spillover Spreading Matrix \n",
title.size = 6,
title.alignment = "center",
left.label.size = ylab_width,
left.label.col = "white",
left.label.text.alignment = "right",
bottom.label.text.angle = 90,
bottom.label.size = xlab_width,
bottom.label.col = "white",
bottom.label.text.alignment = "right",
legend.height = 0.2,
legend.vspace = 0,
row.title = "Fluorochrome \n",
row.title.size = 6,
heat.na.col = "grey",
heat.pal = c("white",
"steelblue",
"steelblue2",
"steelblue3",
"steelblue4",
"purple",
"mediumorchid",
"orchid",
"magenta",
"deeppink2",
"deeppink4",
"red"),
heat.pal.values = c(0,
1,
2,
3,
6,
10,
20,
30,
65,
100,
200,
300))
# Default spillover spread file name
if(is.null(spillover_spread)){
spillover_spread <- paste0(format(Sys.Date(), "%d%m%y"),
"-Spillover-Spread.csv")
}
# Write to csv file
if (!is(spillover_spread, "character")) {
stop("'spillover_spread' should be the name of a csv file.")
} else {
spillover_spread <- file_ext_append(spillover_spread, ".csv")
write.csv(SSM, spillover_spread)
}
return(SSM)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.