xyplot,flowFrame,missing-method {flowViz} | R Documentation |
These functions create Trellis scatter plots (a.k.a. dot plots in the Flow Cytometry community) from flow cytometry data.
## S4 method for signature 'flowFrame,missing' xyplot(x, data, time, xlab, ylab = "", layout, prepanel = prepanel.xyplot.flowframe.time, panel = panel.xyplot.flowframe.time, type = "discrete", ...) prepanel.xyplot.flowframe.time(x, y, frame, time, xlim, ylim, ...) panel.xyplot.flowframe.time(x, y, frame, time, type = "discrete", nrpoints = 0, binSize = 100, ...) ## S4 method for signature 'formula,flowFrame' xyplot(x, data, filter = NULL, overlay = NULL, stats = FALSE, strip.text = NULL, ...) prepanel.xyplot.flowframe(frame, channel.x.name, channel.y.name, x, y, xlim, ylim, ...) panel.xyplot.flowframe(frame, filter = NULL, smooth = TRUE, margin = TRUE, outline = FALSE, channel.x.name, channel.y.name, pch = gp$flow.symbol$pch, alpha = gp$flow.symbol$alpha, cex = gp$flow.symbol$cex, col = gp$flow.symbol$col, gp, xbins = 0, binTrans = sqrt, stats = FALSE, pos = 0.5, digits = 2, abs = FALSE, overlay = NULL, checkName = TRUE, sample.ratio = 1, overlay.symbol = NULL, ...) ## S4 method for signature 'formula,flowSet' xyplot(x, data, ...) prepanel.xyplot.flowset(x, frames, channel.x.name, channel.y.name, xlim, ylim, ...) panel.xyplot.flowset(x, frames, filter = NULL, channel.x, channel.y, overlay = NULL, stats = FALSE, ...) ## S4 method for signature 'formula,view' xyplot(x, data, ...) ## S4 method for signature 'view,missing' xyplot(x, data, ...) ## S4 method for signature 'formula,gateView' xyplot(x, data, filter = NULL, par.settings, ...)
x |
A formula describing the structure of the plot and the variables to
be used in the display. In the |
data,y,frame |
a |
time |
A character string giving the name of the data column recording time. If not provided, we try to guess from the available parameters. |
xlab,ylab |
Labels for data axes, with suitable defaults taken from the formula. |
layout |
These arguments are passed unchanged to the corresponding methods in lattice, and are listed here only because they provide different defaults. See documentation for the original methods for details. |
prepanel |
The prepanel function. See
|
panel |
The panel function. See |
type |
type of rendering; see
|
xlim,ylim |
limits for data axes. If not given, they are taken from the ranges stored in flowFrame |
nrpoints |
The number of points plotted on the smoothed plot in sparse
regions. This is only listed here because we use a different default. See
|
binSize |
The size of a bin (i.e., the number of events within a bin) used for the smoothed average timeline plots. |
filter |
A |
overlay |
The extra cell events plotted on top of the current cell
population. It is a |
stats,pos,digits,abs |
Arguments to control statistics that is
associated with |
strip.text |
A |
channel.x.name,channel.y.name |
Character strings giving corresponding names used to match filter parameters if applicable. |
smooth |
Logical. If |
margin |
Logical indicating whether to truncate the density estimation
on the margins of the measurement range and plot margin events as lines if
|
outline |
Logical, specifying whether to add the boundaries of a gate
to the plot when |
pch,cex,col,alpha |
Graphical parameters used when |
gp |
A list of graphical parameters that are passed down to the low
level panel functions. This is for internal use only. The public user
interface to set graphical parameters is either |
xbins |
The argument passed to |
binTrans |
The argument passed to
|
checkName |
|
sample.ratio |
|
overlay.symbol |
list of the lattice graphic parameters to format the overlay points. |
frames |
An environment containing frame-specific data. |
channel.x,channel.y |
Expressions defining the x and y variables in terms of columns in the data. Can involve functions or multiple columns from the data, however this usage is discouraged. |
par.settings |
A list of lists of graphical parameters. See |
... |
marker.only More arguments, usually passed on to the underlying lattice methods. |
The implementation of xyplot
in flowViz
is very close to the
original lattice
version. Concepts like conditioning and the use of
panels apply directly to the flow cytometry data. The single fundamental
difference is that conditioning variables are not evaluated in the context
of the raw data, but rather in the phenoData
slot environment (only
for the flowSet
methods. Thus, we can directly condition on pheotypic
variables like sample groups, patients or treatments.
In the formula interface, the primary and secondary variables (separated by
the tilde) have to be valid parameter names. Please note that frequently
used variants like FSC-H
and SSC-H
are not syntactically
correct R symbols, and need to be wrapped in ` `
. E.g.,
`FSC-H`
. For flowSets
, the use of a conditioning variable is
optional. We implicitely condition on flowFrames
and the default is
to arrange panels by sample names.
signature(x = "flowFrame", data = "missing")
: Creates
diagnostic time series plots of flow parameter values against time. These
plots are useful to detect quality issues in the raw data. If not provided
explicitely via the tine
argument, the time parameter will be
automatically detected. The additional arguments xlab
, ylab
,
nrpoints
, and layout
are only listed because flowViz
provides different defaults. Internally, they are directly passed on to the
underlying lattice
functions
. Argument type
can be a
combination of any of the types allowed in lattice xyplots
, or
discrete
, in which case a smoothed average of the parameter against
time is plotted. binSize
controls the binning that is used for the
smoothing procedure.
signature(x = "formula", data = "flowFrame")
: Creates
scatter plots (a.k.a. dot plots) of a pair of FCM channels. Depending on the
setting of the smooth
argument, the data will be rendered as a
partially smoothed density estimate (smooth=TRUE
, the default) or as
a regular scatter plot with separate points for individual events. The
formula interface allows for fairly general plotting, however there are
certain limitations on the use of expressions as part of the formulae.
Unless you are sure about what you are doing, you should transform the raw
data in a separate step using one of the tools in the
flowCore
package rather than inline
using the formula interface. The method allows to superimpose gating results
though the filter
argument. If smooth=TRUE
, we try to add
spherical 2D representations of the gates if applicable. For
smooth=FALSE
, gates are indicated by a grouping mechanism using
different point shapes or colors (unless outline
is also TRUE
,
in which case the gate outlines are superimposed in addition to the
grouping). Argument margins
controls how events on the margins of the
measurement range are treated. The default (TRUE
) is to discard them
from any density estimation and later add them as separate glyphs. Argument par.settings
can be used
to supply lists of graphical parameters. See flowViz.par.set
for details on controlling graphical
parameters in these plots.
signature(x = "formula", data = "flowSet")
: Scatter
plots from a flowSet
object. We allow for conditioning on variables
in the phenoData
slot of the flowSet
. All additional arguments
that apply to the flowFrame
method are also valid for
flowSets
.
F. Hahne, D. Sarkar
Not all standard lattice arguments will have the intended effect, but many should. For a fuller description of possible arguments and their effects, consult documentation on lattice.
data(GvHD) GvHD <- GvHD[pData(GvHD)$Patient %in% 5:6] ## a bivariate scatterplot ## by default ('smooth=TRUE') panel.smoothScatter is used xyplot(`FSC-H` ~ `SSC-H`, GvHD[["s5a05"]], nbin = 100, main="A single flowFrame") ## A non-smooth version of the same data xyplot(`FSC-H` ~ `SSC-H`, GvHD[["s5a05"]], nbin = 100, main="A single flowFrame", smooth=FALSE) ## A non-smooth version of the same data with customerized color scheme require(IDPmisc) colramp <- colorRampPalette(IDPcolorRamp(21)) xyplot(`FSC-H` ~ `SSC-H`, GvHD[["s5a05"]], nbin = 100, main="A single flowFrame", smooth=FALSE, colramp=colramp, pch=20, cex=0.1) ## A hexbin version of non-smooth scatter plot xyplot(`FSC-H` ~ `SSC-H`, GvHD[["s5a05"]], xbin = 128 ,main="A single flowFrame", smooth=FALSE) ## Visual artifacts created by the pileup of margin events xyplot(`FSC-H` ~ `SSC-H`, GvHD[["s5a05"]], nbin = 100, main="A single flowFrame", margin=FALSE) ## simple bivariate scatter plot (a.k.a. dot plot) ## for the whole flowSet, conditioning on Patient and ## Visit xyplot(`SSC-H` ~ `FSC-H` | Patient:Visit, data = GvHD) ## Same bivariate scatter plot with replacing default color require(IDPmisc) cols <- colorRampPalette(IDPcolorRamp(21)) xyplot(`SSC-H` ~ `FSC-H` | Patient:Visit, data = GvHD, colramp=cols) ## several examples with time on the X axis ## first for a flowFrame xyplot(GvHD[[1]]) ## and for flowSets xyplot(`FSC-H` ~ Time | Visit, GvHD, smooth = FALSE, type = "l", subset = (Patient == 5), xbin = 32) xyplot(`FSC-H` ~ Time | Patient+Visit, GvHD, smooth = FALSE, type = "a", strip = FALSE, strip.left = TRUE, aspect = "xy", xbin = 32) ## combine plots for two channels ssc.time <- xyplot(`SSC-H` ~ Time | factor(Patient):factor(Visit), GvHD, smooth = FALSE, type = "a", strip = FALSE, strip.left = strip.custom(horizontal = TRUE), par.strip.text = list(lines = 3), between = list(y = rep(c(0, 0.5), c(6, 1))), scales = list(x = list(axs = "i"), y = list(draw = FALSE)), layout = c(1, 14), xbin = 32) fsc.time <- xyplot(`FSC-H` ~ Time | factor(Patient):factor(Visit), GvHD, smooth = FALSE, type = "a", strip = FALSE, strip.left = strip.custom(horizontal = TRUE), par.strip.text = list(lines = 3), between = list(y = rep(c(0, 0.5), c(6, 1))), scales = list(x = list(axs = "i"), y = list(draw = FALSE)), layout = c(1, 14), xbin = 32) plot(fsc.time, split = c(1, 1, 2, 1)) plot(ssc.time, split = c(2, 1, 2, 1), newpage = FALSE) ## saving plots as variables allows more manipulation plot(update(fsc.time[8:14], layout = c(1, 7)), split = c(1, 1, 1, 2)) plot(update(ssc.time[8:14], layout = c(1, 7)), split = c(1, 2, 1, 2), newpage = FALSE) ## displaying filters n2gate <- norm2Filter("SSC-H", "FSC-H") xyplot(`SSC-H` ~ `FSC-H` | Patient:Visit, data = GvHD, filter=n2gate, subset=Patient==5) xyplot(`SSC-H` ~ `FSC-H` | Patient:Visit, data=transform("SSC-H"=asinh,"FSC-H"=asinh) %on% GvHD, smooth=FALSE, filter=n2gate, subset=Patient == 5, xbin = 32) ## displaying filters with stats n2gate.results <- filter(GvHD, n2gate) xyplot(`SSC-H` ~ `FSC-H` | Visit, data=GvHD, subset=Patient == "6", filter=n2gate.results, smooth=FALSE, xbin = 32 ,stats=TRUE ,abs=TRUE ,digits=3 ) ## displaying multiple filters in one panel with stats recGate1<-rectangleGate("FL3-H"=c(2.3,4.1),"FL2-H"=c(6.8,9)) recGate2<-rectangleGate("FL3-H"=c(1,3),"FL2-H"=c(4,6)) filters1<-filters(list(recGate1,recGate2)) trans<-transform("FL2-H"=asinh,"FL3-H"=asinh) trans_data<-transform(GvHD[1:2],trans) #replicate filters object across samples flist <- list(filters1 , filters1) names(flist) <- sampleNames(trans_data) xyplot(`FL2-H` ~ `FL3-H` ,data=trans_data ,filter= flist ,stats=TRUE ,margin=FALSE , xbin = 32 , smooth = FALSE ) #display recGate2 as a overlay overlay <- Subset(trans_data,recGate1) xyplot(`FL2-H` ~ `FL3-H` ,data=trans_data ,filter=recGate2 ,stats=TRUE ,margin=FALSE , smooth = FALSE , xbin = 32 ,overlay= list(rect2 = overlay) ,par.settings = list(overlay.symbol = list(cex = 0.1)) )