Generates an overlay plot containing all objects in the GatingSet
Source:R/Utility_DensityOverlay.R
Utility_DensityOverlay.Rd
Generates an overlay plot containing all objects in the GatingSet
Usage
Utility_DensityOverlay(
gs,
subset,
inverse.transform = FALSE,
TheX = NULL,
TheFill,
returntype,
outpath,
filename,
therows = 3,
thecolumns = 3,
width = 7,
height = 9
)
Arguments
- gs
A GatingSet object
- subset
The desired Gating Hierarchy node (ex. "lymphocytes")
- inverse.transform
Default is FALSE, reverses transformation applied to GatingSet as its converted to a CytoSet
- TheX
A desired marker to plot, leave NULL for all markers
- TheFill
A desired marker to color individual specimens by (as named in pData)
- returntype
Whether to return a "pdf", "patchwork" or "plots"
- outpath
Desired storage location
- filename
The file name for the new .pdf
- therows
The desired number of rows per page
- thecolumns
The desired number of columns per page
- width
The desired page width
- height
The desired page height
Examples
library(BiocGenerics)
library(flowCore)
library(flowWorkspace)
library(openCyto)
library(data.table)
File_Location <- system.file("extdata", package = "Luciernaga")
FCS_Files <- list.files(path = File_Location, pattern = ".fcs",
full.names = TRUE)
Unmixed_FullStained <- FCS_Files[grep("Unmixed", FCS_Files)]
UnmixedFCSFiles <- Unmixed_FullStained[1:2]
UnmixedCytoSet <- load_cytoset_from_fcs(UnmixedFCSFiles[1:2],
truncate_max_range = FALSE,transformation = FALSE)
UnmixedGatingSet <- GatingSet(UnmixedCytoSet)
Markers <- colnames(UnmixedCytoSet)
KeptMarkers <- Markers[-grep("Time|FS|SC|SS|Original|-W$|-H$|AF", Markers)]
MyBiexponentialTransform <- flowjo_biexp_trans(channelRange = 256,
maxValue = 1000000,pos = 4.5, neg = 0, widthBasis = -1000)
TransformList <- transformerList(KeptMarkers, MyBiexponentialTransform)
UnmixedGatingSet <- flowWorkspace::transform(UnmixedGatingSet, TransformList)
FileLocation <- system.file("extdata", package = "Luciernaga")
UnmixedGates <- fread(file.path(path = FileLocation, pattern = 'GatesUnmixed.csv'))
UnmixedGating <- gatingTemplate(UnmixedGates)
#> Adding population:singletsFSC
#> Adding population:singletsSSC
#> Adding population:singletsSSCB
#> Adding population:nonDebris
#> Adding population:lymphocytes
#> Adding population:live
gt_gating(UnmixedGating, UnmixedGatingSet)
#> Gating for 'singletsFSC'
#> done!
#> done.
#> Gating for 'singletsSSC'
#> done!
#> done.
#> Gating for 'singletsSSCB'
#> done!
#> done.
#> Gating for 'nonDebris'
#> done!
#> done.
#> Gating for 'lymphocytes'
#> The prior specification has no effect when usePrior=no
#> Using the serial version of flowClust
#> The prior specification has no effect when usePrior=no
#> Using the serial version of flowClust
#> done!
#> done.
#> Gating for 'live'
#> done!
#> done.
#> finished.
removestrings <- c("DTR_", ".fcs")
StorageLocation <- file.path("C:", "Users", "JohnDoe", "Desktop")
Condition <- data.frame(Condition=c("Ctrl", "Ctrl"))
pd <- pData(UnmixedGatingSet)
new_pd <- cbind(pd, Condition)
pData(UnmixedGatingSet) <- new_pd
Plot <- Utility_DensityOverlay(gs=UnmixedGatingSet, subset="lymphocytes",
TheX="APC-Fire 810-A",TheFill="Condition", returntype="plots",
outpath="C:/Users/JohnDoe/Desktop/", filename="CD38_Expression")