plotVolcano {OUTRIDER} | R Documentation |
The OUTRIDER package provides mutliple functions to visualize the data and the results of a full data set analysis.
This is the list of all plotting function provided by OUTRIDER:
plotAberrantPerSample()
plotVolcano()
plotExpressionRank()
plotQQ()
plotCountCorHeatmap()
plotFPKM()
plotDispEsts()
plotPowerAnalysis()
plotEncDimSearch()
For a detailed description of each plot function please see the details. Most of the functions share the same parameters.
plotVolcano(ods, sampleID, main, padjCutoff = 0.05, zScoreCutoff = 0, pch = 16, basePlot = FALSE, col = c("gray", "firebrick")) plotQQ(ods, geneID, main, global = FALSE, padjCutoff = 0.05, zScoreCutoff = 0, samplePoints = TRUE, legendPos = "topleft", outlierRatio = 0.001, conf.alpha = 0.05, pch = 16, xlim = NULL, ylim = NULL, col = NULL) plotExpressionRank(ods, geneID, main, padjCutoff = 0.05, zScoreCutoff = 0, normalized = TRUE, basePlot = FALSE, log = TRUE, col = c("gray", "firebrick"), ...) plotCountCorHeatmap(ods, normalized = TRUE, rowCentered = TRUE, rowCoFactor = NULL, rowColSet = "Set1", colCoFactor = NULL, colColSet = "Set2", nCluster = 4, main = "Count correlation heatmap", dendrogram = "both", basePlot = TRUE, names = c("both", "row", "col", "none"), ...) plotAberrantPerSample(ods, main, padjCutoff = 0.05, zScoreCutoff = 0, outlierRatio = 0.001, col = brewer.pal(3, "Dark2")[c(1, 2)], yadjust = c(1.2, 1.2), labLine = c(3.5, 3), ylab = "#Aberrantly expressed genes", labCex = par()$cex, ymax = NULL, ...) plotFPKM(ods) ## S4 method for signature 'OutriderDataSet' plotDispEsts(object, compareDisp, xlim, ylim, main = "Dispersion estimates versus mean expression", ...) plotPowerAnalysis(ods) plotEncDimSearch(ods)
ods, object |
An OutriderDataSet object. |
sampleID, geneID |
A sample or gene ID, which should be plotted. Can also be a vector. Integers are treated as indices. |
main |
Title for the plot, if missing a default title will be used. |
padjCutoff, zScoreCutoff |
Significance or Z-score cutoff to mark outliers |
pch |
Integer or character to be used for plotting the points |
basePlot |
if |
col |
Set color for the points. If set, it must be a character vector of length 2. (1. normal point; 2. outlier point) |
global |
Flag to plot a global Q-Q plot, default FALSE |
samplePoints |
Sample points for Q-Q plot, defaults to max 30k points |
legendPos |
Set legendpos, by default topleft. |
outlierRatio |
The fraction to be used for the outlier sample filtering |
conf.alpha |
If set, a confidence interval is plotted, defaults to 0.05 |
xlim, ylim |
The x/y limits for the plot or NULL to use the full data range |
normalized |
If TRUE, the normalized counts are used, the default, otherwise the raw counts |
log |
If TRUE, the default, counts are plotted in log10. |
... |
Additional parameters passed to plot() or plot_ly() if not stated otherwise in the details for each plot function |
rowCentered |
If TRUE, the counts are row-wise (gene-wise) centered |
rowCoFactor |
A vector of co-factors for color coding the rows |
rowColSet |
A vector of colors or a color set from RColorBrewer |
colCoFactor |
A vector of co-factors for color coding the columns |
colColSet |
A vector of colors or a color set from RColorBrewer |
nCluster |
An integer to be used for cutting the dendrogram into groups. If this argument is set the resulting clusters are saved in the returned OutriderDataSet. |
dendrogram |
A character string indicating whether to draw 'none', 'row', 'column' or 'both' dendrograms. |
names |
character string indicating whether to draw 'none', 'row', 'col', or 'both' names. |
yadjust |
Option to adjust position of Median and 90 percentile labels. |
labLine |
Option to move axis labels |
ylab |
The y axis label |
labCex |
The label cex parameter |
ymax |
If set, ymax is the upper bound for the plot range on the y axis. |
compareDisp |
If TRUE, the default, and if the autoCorrect normalization was used it computes the dispersion without autoCorrect and plots it for comparison. |
plotAberrantPerSample
: The number of aberrant events per sample are
plotted sorted by rank. The ... parameters are passed on to the
aberrant
function.
plotVolcano
: the volcano plot is sample-centric. It plots for a given
sample the negative log10 nominal P-values against the Z-scores for all
genes.
plotExpressionRank
: This function plots for a given gene the
expression level against the expression rank for all samples. This can
be used with normalized and unnormalized expression values.
plotQQ
: the quantile-quantile plot for a given gene or if
global
is set to TRUE
over the full data set. Here the
observed P-values are plotted against the expected ones in the negative
log10 space.
plotCountCorHeatmap
: The correlation heatmap of the count data
of the full data set. Default the values are log transformed and
row centered. This function returns an OutriderDataSet with annotated
clusters if requested. The ... arguments are passed to the
heatmap.2
function.
plotFPKM
: The distribution of FPKM values. If the OutriderDataSet
object contains the passedFilter
column, it will plot both FPKM
distributions for the expressed genes and for the filtered genes.
plotDispEsts
: Plots the dispersion of the OutriderDataSet
model against the normalized mean count. If autoCorrect is used it will also
estimate the dispersion without normalization for comparison.
plotPowerAnalysis
: The power analysis plot should give the user a
ruff estimate of the events one can be detected with OUTRIDER. Based on
the dispersion of the provided OUTRIDER data set the theoretical P-value
over the mean expression is plotted. This is done for different expression
levels. The curves are smooths to make the reading of the plot easier.
If base R graphics are used nothing is returned else the plotly or the gplot object is returned.
ods <- makeExampleOutriderDataSet(dataset="Kremer") implementation <- 'autoencoder' ods <- filterExpression(ods, minCounts=TRUE) ods <- OUTRIDER(ods, implementation=implementation) plotAberrantPerSample(ods) plotVolcano(ods, 1) plotVolcano(ods, 'MUC1404', basePlot=TRUE) plotExpressionRank(ods, 1) plotExpressionRank(ods, "FAAH", normalized=FALSE, log=FALSE, main="Over expression outlier", basePlot=TRUE) plotQQ(ods, 1) plotQQ(ods, global=TRUE, outlierRatio=0.001) sex <- sample(c("female", "male"), dim(ods)[2], replace=TRUE) colData(ods)$sex <- sex ods <- plotCountCorHeatmap(ods, colCoFactor="sex", normalized=FALSE) ods <- plotCountCorHeatmap(ods, nCluster=4) head(colData(ods)$clusterNumber) mcols(ods)$basepairs <- 1 mcols(ods)$passedFilter <- rowMeans(counts(ods)) > 10 plotFPKM(ods) plotDispEsts(ods, compareDisp=FALSE) plotPowerAnalysis(ods) ## Not run: ods <- findEncodingDim(ods) plotEncDimSearch(ods) ## End(Not run)