normalize {scater}R Documentation

Normalise a SingleCellExperiment object using pre-computed size factors

Description

Compute normalised expression values from a SingleCellExperiment object using the size factors stored in the object. Return the object with the normalised expression values added.

Usage

normalizeSCE(object, exprs_values = "counts", return_log = TRUE,
  log_exprs_offset = NULL, centre_size_factors = TRUE,
  return_norm_as_exprs = TRUE)

## S4 method for signature 'SingleCellExperiment'
normalize(object, exprs_values = "counts",
  return_log = TRUE, log_exprs_offset = NULL, centre_size_factors = TRUE,
  return_norm_as_exprs = TRUE)

normalise(...)

Arguments

object

a SingleCellExperiment object.

exprs_values

character string indicating which slot of the assayData from the SingleCellExperiment object should be used to compute log-transformed expression values. Valid options are 'counts', 'tpm', 'cpm' and 'fpkm'. Defaults to the first available value of the options in the order shown.

return_log

logical(1), should normalized values be returned on the log scale? Default is TRUE. If TRUE, output is stored as "logcounts" in the returned object; if FALSE output is stored as "normcounts"

log_exprs_offset

scalar numeric value giving the offset to add when taking log2 of normalised values to return as expression values. If NULL, value is taken from metadata(object)$log.exprs.offset if defined, otherwise 1.

centre_size_factors

logical, should size factors centred at unity be stored in the returned object if exprs_values="counts"? Defaults to TRUE. Regardless, centred size factors will always be used to calculate exprs from count data. This argument is ignored for other exprs_values, where no size factors are used/modified.

return_norm_as_exprs

logical, should the normalised expression values be returned to the exprs slot of the object? Default is TRUE. If FALSE, values in the exprs slot will be left untouched. Regardless, normalised expression values will be returned in the norm_exprs(object) slot.

...

arguments passed to normalize when calling normalise.

Details

normalize is exactly the same as normalise, the option provided for those who have a preference for North American or British/Australian spelling.

Value

an SingleCellExperiment object

Warning about centred size factors

Centring the size factors ensures that the computed exprs can be interpreted as being on the same scale as log-counts. This does not affect relative comparisons between cells in the same object, as all size factors are scaled by the same amount. However, if two different SingleCellExperiment objects are run separately through normalize, the size factors in each object will be rescaled differently. This means that the size factors and exprs will not be comparable between objects.

This lack of comparability is not always obvious. For example, if we subsetted an existing SingleCellExperiment, and ran normalize separately on each subset, the resulting exprs in each subsetted object would not be comparable to each other. This is despite the fact that all cells were originally derived from a single SingleCellExperiment object.

In general, it is advisable to only compare size factors and exprs between cells in one SingleCellExperiment object. If objects are to be combined, new size factors should be computed using all cells in the combined object, followed by running normalize.

Author(s)

Davis McCarthy and Aaron Lun

Examples

data("sc_example_counts")
data("sc_example_cell_info")
example_sce <- SingleCellExperiment(
assays = list(counts = sc_example_counts), colData = sc_example_cell_info)
keep_gene <- rowSums(counts(example_sce)) > 0
example_sce <- example_sce[keep_gene,]

## Apply TMM normalisation taking into account all genes
example_sce <- normaliseExprs(example_sce, method = "TMM")
## Scale counts relative to a set of control features (here the first 100 features)
example_sce <- normaliseExprs(example_sce, method = "none",
feature_set = 1:100)

## normalize the object using the saved size factors
example_sce <- normalize(example_sce)


[Package scater version 1.6.3 Index]