transformGamPoi: Variance Stabilizing Transformation for Gamma Poisson Data

View source: R/transformGamPoi.R

transformGamPoiR Documentation

Variance Stabilizing Transformation for Gamma Poisson Data

Description

Variance Stabilizing Transformation for Gamma Poisson Data

Usage

transformGamPoi(
  data,
  transformation = c("acosh", "shifted_log", "randomized_quantile_residuals",
    "pearson_residuals", "analytic_pearson_residuals"),
  overdispersion = 0.05,
  size_factors = TRUE,
  ...,
  on_disk = NULL,
  verbose = FALSE
)

Arguments

data

any matrix-like object (e.g. matrix, dgCMatrix, DelayedArray, HDF5Matrix) with one column per sample and row per gene. It can also be an object of type glmGamPoi, in which case it is directly used to calculate the variance-stabilized values.

transformation

one of c("acosh", "shifted_log", "randomized_quantile_residuals", "pearson_residuals", "analytic_pearson_residuals"). See acosh_transform, shifted_log_transform, or residual_transform for more information.

overdispersion

the simplest count model is the Poisson model. However, the Poisson model assumes that variance = mean. For many applications this is too rigid and the Gamma-Poisson allows a more flexible mean-variance relation (variance = mean + mean^2 * overdispersion).
overdispersion can either be

  • a single boolean that indicates if an overdispersion is estimated for each gene.

  • a numeric vector of length nrow(data) fixing the overdispersion to those values.

  • the string "global" to indicate that one dispersion is fit across all genes.

Note that overdispersion = 0 and overdispersion = FALSE are equivalent and both reduce the Gamma-Poisson to the classical Poisson model. Default: 0.05 which is roughly the overdispersion observed on ostensibly homogeneous cell lines.

size_factors

in large scale experiments, each sample is typically of different size (for example different sequencing depths). A size factor is an internal mechanism of GLMs to correct for this effect.
size_factors is either a numeric vector with positive entries that has the same lengths as columns in the data that specifies the size factors that are used. Or it can be a string that species the method that is used to estimate the size factors (one of c("normed_sum", "deconvolution", "poscounts")). Note that "normed_sum" and "poscounts" are fairly simple methods and can lead to suboptimal results. For the best performance, I recommend to use size_factors = "deconvolution" which calls scran::calculateSumFactors(). However, you need to separately install the scran package from Bioconductor for this method to work. Also note that size_factors = 1 and size_factors = FALSE are equivalent. If only a single gene is given, no size factor is estimated (ie. size_factors = 1). Default: "normed_sum".

...

additional parameters passed to acosh_transform, shifted_log_transform, or residual_transform

on_disk

a boolean that indicates if the dataset is loaded into memory or if it is kept on disk to reduce the memory usage. Processing in memory can be significantly faster than on disk. Default: NULL which means that the data is only processed in memory if data is an in-memory data structure.

verbose

boolean that decides if information about the individual steps are printed. Default: FALSE

Value

a matrix (or a vector if the input is a vector) with the transformed values.

References

Ahlmann-Eltze, Constantin, and Wolfgang Huber. "Transformation and Preprocessing of Single-Cell RNA-Seq Data." bioRxiv (2021).

Ahlmann-Eltze, Constantin, and Wolfgang Huber. "glmGamPoi: Fitting Gamma-Poisson Generalized Linear Models on Single Cell Count Data." Bioinformatics (2020)

Dunn, Peter K., and Gordon K. Smyth. "Randomized quantile residuals." Journal of Computational and Graphical Statistics 5.3 (1996): 236-244.

Hafemeister, Christoph, and Rahul Satija. "Normalization and variance stabilization of single-cell RNA-seq data using regularized negative binomial regression." Genome biology 20.1 (2019): 1-15.

Hafemeister, Christoph, and Rahul Satija. "Analyzing scRNA-seq data with the sctransform and offset models." (2020)

Lause, Jan, Philipp Berens, and Dmitry Kobak. "Analytic Pearson residuals for normalization of single-cell RNA-seq UMI data." Genome Biology (2021).

See Also

acosh_transform, shifted_log_transform, and residual_transform

Examples

  # Load a single cell dataset
  sce <- TENxPBMCData::TENxPBMCData("pbmc4k")
  # Reduce size for this example
  set.seed(1)
  sce_red <- sce[sample(which(rowSums2(counts(sce)) > 0), 1000),
                 sample(ncol(sce), 200)]

  assay(sce_red, "acosh") <- transformGamPoi(sce_red, "acosh")
  assay(sce_red, "shifted_log") <- transformGamPoi(sce_red, "shifted_log")

  # Residual Based Variance Stabilizing Transformation
  rq <- transformGamPoi(sce_red, transformation = "randomized_quantile", on_disk = FALSE,
                        verbose = TRUE)
  pearson <- transformGamPoi(sce_red, transformation = "pearson", on_disk = FALSE, verbose = TRUE)

  plot(rowMeans2(counts(sce_red)), rowVars(assay(sce_red, "acosh")), log = "x")
  points(rowMeans2(counts(sce_red)), rowVars(assay(sce_red, "shifted_log")), col = "red")
  points(rowMeans2(counts(sce_red)), rowVars(rq), col = "blue")


  # Plot first two principal components
  acosh_pca <- prcomp(t(assay(sce_red, "acosh")), rank. = 2)
  rq_pca <- prcomp(t(rq), rank. = 2)
  pearson_pca <- prcomp(t(pearson), rank. = 2)

  plot(acosh_pca$x, asp = 1)
  points(rq_pca$x, col = "blue")
  points(pearson_pca$x, col = "green")


const-ae/transformGamPoi documentation built on April 14, 2023, 11:33 p.m.