plotRLE {scater} | R Documentation |
Produce a relative log expression (RLE) plot of one or more transformations of cell expression values.
plotRLE(object, exprs_mats = list(logcounts = "logcounts"), exprs_logged = c(TRUE), colour_by = NULL, style = "minimal", legend = "auto", order_by_colour = TRUE, ncol = 1, ...)
object |
an |
exprs_mats |
named list of expression matrices. Entries can either be a
character string, in which case the corresponding expression matrix will be
extracted from the SingleCellExperiment |
exprs_logged |
logical vector of same length as |
colour_by |
character string defining the column of |
style |
character(1), either |
legend |
character, specifying how the legend(s) be shown? Default is
|
order_by_colour |
logical, should cells be ordered (grouped) by the
|
ncol |
integer, number of columns for the facetting of the plot. Default is 1. |
... |
further arguments passed to |
Unwanted variation can be highly problematic and so its detection is often crucial. Relative log expression (RLE) plots are a powerful tool for visualising such variation in high dimensional data. RLE plots are particularly useful for assessing whether a procedure aimed at removing unwanted variation, i.e. a normalisation procedure, has been successful. These plots, while originally devised for gene expression data from microarrays, can also be used to reveal unwanted variation in single-cell expression data, where such variation can be problematic.
If style is "full", as usual with boxplots, the box shows the inter-quartile range and whiskers extend no more than 1.5 * IQR from the hinge (the 25th or 75th percentile). Data beyond the whiskers are called outliers and are plotted individually. The median (50th percentile) is shown with a white bar.
If style is "minimal", then median is shown with a circle, the IQR in a grey line, and "whiskers" (as defined above) for the plots are shown with coloured lines. No outliers are shown for this plot style.
a ggplot plot object
Davis McCarthy
Gandolfo LC, Speed TP. RLE Plots: Visualising Unwanted Variation in High Dimensional Data. arXiv [stat.ME]. 2017. Available: http://arxiv.org/abs/1704.03590
data("sc_example_counts") data("sc_example_cell_info") example_sce <- SingleCellExperiment( assays = list(counts = sc_example_counts), colData = sc_example_cell_info) example_sce <- normalize(example_sce) drop_genes <- apply(logcounts(example_sce), 1, function(x) {var(x) == 0}) example_sce <- example_sce[!drop_genes, ] plotRLE(example_sce, list(logcounts= "logcounts", counts = "counts"), c(TRUE, FALSE), colour_by = "Mutation_Status", style = "minimal") plotRLE(example_sce, list(logcounts = "logcounts", counts = "counts"), c(TRUE, FALSE), colour_by = "Mutation_Status", style = "full", outlier.alpha = 0.1, outlier.shape = 3, outlier.size = 0)