ExampleRunThrough/4_PostCompassPlots_Example.R

In seshadrilab/COMPASSHelpers: Collection of helper functions for analyzing FCM data

#!/usr/bin/env Rscript
library(flowHelpers)
library(flowWorkspace)
library(here)

compassResultsCD4Paths <- c(here::here("ExampleRunThrough/output/CompassOutput/CD4/4+_ESAT-6/COMPASSResult_4+_ESAT-6.rds"),
                       here::here("ExampleRunThrough/output/CompassOutput/CD4/4+_CFP-10/COMPASSResult_4+_CFP-10.rds"),
                       here::here("ExampleRunThrough/output/CompassOutput/CD4/4+_Ag85A/COMPASSResult_4+_Ag85A.rds"),
                       here::here("ExampleRunThrough/output/CompassOutput/CD4/4+_Ag85B/COMPASSResult_4+_Ag85B.rds"),
                       here::here("ExampleRunThrough/output/CompassOutput/CD4/4+_TB10.4/COMPASSResult_4+_TB10.4.rds"))
compassResultsCD4 <- lapply(compassResultsCD4Paths, readRDS)
names(compassResultsCD4) <- c("ESAT-6", "CFP-10", "Ag85A", "Ag85B", "TB10.4")
outdir <- here::here("ExampleRunThrough/output/PostCompassPlots")

# Plot the mean probablity of response for all interesting CD4+ COMPASSResults
# Re-creating Fig. 4: http://www.jimmunol.org/content/195/10/4595
# NOTE: The plotCOMPASSResultStack is not yet a part of the main COMPASS repository,
# so this function will not work on a regular COMPASS install.
# https://github.com/malisas/COMPASS/blob/stack-heatmaps/R/plotMeanGammaMulti.R
# PATIENT ID is not included in row_annotation because it is already the individual_id of the COMPASSResult (so it gets automatically included)
png(filename=file.path(outdir, "StackedHeatmapMeanProbResponse_CD4.png"),
    width=1000, height=1000)
plotCOMPASSResultStack(compassResultsCD4,
  row_annotation = c("Antigen", "QuantiferonInTubeResult"),
  variable = "Antigen",
  main = "Heatmap of Mean Probability of Response to Antigens, CD4+",
  fontsize = 14, fontsize_row = 13, fontsize_col = 11)
dev.off()

# load in the GatingSetlist generated by step 2
gatingSetListDir <- "/media/malisa/SharedStorage/uw/20170331_TB-ICS_ACS/GatingSetList4COMPASS"

# Now pick an interesting, representative patient and highlight their polyfunctional cells in FACS plots
# PATIENT ID 03-0311, Peptide ESAT-6, positive for IFNg, CD154, IL2, and TNFa, negative for IL22 and IL17a
# left out IL4 to be consistent w/ the paper
highlight.boolean.subset.flow.plot(path=gatingSetListDir,
                                   outdir=outdir,
                                   individualsCol="PATIENT ID",
                                   individual="03-0311",
                                   conditioncol="Peptide",
                                   exp="ESAT-6",
                                   ctrl="DMSO",
                                   parentsubset="4+",
                                   boolsubset="4+/IFNg+&4+/154+&4+/IL2+&4+/TNFa+&!4+/IL22+&!4+/IL17+",
                                   xaxis="TNFa",
                                   yaxis="IFNg",
                                   geomTextY=250)

# Same as above, TB10.4
highlight.boolean.subset.flow.plot(path=gatingSetListDir,
                                   outdir=outdir,
                                   individualsCol="PATIENT ID",
                                   individual="03-0311",
                                   conditioncol="Peptide",
                                   exp="TB10.4",
                                   ctrl="DMSO",
                                   parentsubset="4+",
                                   boolsubset="4+/IFNg+&4+/154+&4+/IL2+&4+/TNFa+&!4+/IL22+&!4+/IL17+",
                                   xaxis="TNFa",
                                   yaxis="IFNg",
                                   geomTextY=250)

# Look at the same subset across all individuals, in a boxplot comparing M.tb-uninfected to M.tb-infected
# For ESAT-6 and TB10.4:
antigens <- c("ESAT-6", "TB10.4")
for (i in 1:length(antigens)) {
  boxplot.boolean.subset.proportions(path=gatingSetListDir,
                                     outdir=outdir,
                                     conditioncol="Peptide",
                                     exp=antigens[[i]],
                                     ctrl="DMSO",
                                     conditioncol2="QuantiferonInTubeResult",
                                     parentsubset="4+",
                                     boolsubset="4+/IFNg+&4+/154+&4+/IL2+&4+/TNFa+&!4+/IL22+&!4+/IL17+",
                                     uniqueSamplesCol="PATIENT ID",
                                     ylimits=c(0, 0.007))
}
# Note: this is NOT comparable to Figure 3C because Figure 3C doesn't account for the difference between DMSO and the peptide of interest