Skip to contents

Overlay plot two .fcs files from two different gs.

Usage

Utility_ParallelNbyNPlots(
  x,
  y,
  sample.name,
  removestrings,
  Override = FALSE,
  marginsubset,
  gatesubset,
  ycolumn,
  bins,
  clearance,
  colorX,
  colorY,
  gatelines,
  reference = NULL,
  outpath,
  pdf
)

Arguments

x

The first gs object

y

The second gs object

sample.name

The keyword under which sample names are stored

removestrings

A list of values to remove from sample names

Override

Exclude raw columns with -A$

marginsubset

The gs subset that defines the plot margin

gatesubset

The gs subset of interest

ycolumn

The desired y-column for the comparisons

bins

Desired number of hex bins

clearance

A multiplication factor for margin wiggle room (0.2)

colorX

Color for the x gs

colorY

Color for the y gs

gatelines

Whether to plot .csv specified gate lines

reference

Reference for .csv specified gate lines

outpath

Location which to store the output

pdf

Whether to return as a pdf

Value

Either list ggplot objects or a pdf object

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")

Plot <- Utility_ParallelNbyNPlots(x=UnmixedGatingSet[1], y = UnmixedGatingSet[2],
 sample.name = "GROUPNAME", removestrings = ".fcs", Override = FALSE,
 marginsubset = "lymphocytes", gatesubset = "live", ycolumn = "Spark Blue 550-A",
  bins = 120, clearance = 0.2, colorX = "lightblue", colorY = "orange", gatelines = FALSE,
  reference = NULL, outpath = StorageLocation,pdf = FALSE)