testPseudotime: Test for differences along pseudotime

testPseudotimeR Documentation

Test for differences along pseudotime

Description

Implements a simple method of testing for significant differences with respect to pseudotime, based on fitting linear models with a spline basis matrix.

Usage

testPseudotime(x, ...)

## S4 method for signature 'ANY'
testPseudotime(
  x,
  pseudotime,
  df = 5,
  get.lfc = TRUE,
  get.spline.coef = FALSE,
  trend.only = TRUE,
  block = NULL,
  row.data = NULL,
  BPPARAM = NULL
)

## S4 method for signature 'SummarizedExperiment'
testPseudotime(x, ..., assay.type = "logcounts")

Arguments

x

A numeric matrix-like object containing log-expression values for cells (columns) and genes (rows). Alternatively, a SummarizedExperiment containing such a matrix.

...

For the generic, further arguments to pass to specific method.

For the SummarizedExperiment method, further arguments to pass to the ANY method.

pseudotime

A numeric vector of length equal to ncol(x), containing the pseudotime orderings along a single lineage. Alternatively, a numeric matrix with number of rows equal to ncol(x), where each column contains an ordering across one of multiple lineages. Alternatively, a PseudotimeOrdering object containing such a matrix in pathStat(pseudotime).

df

Integer scalar specifying the degrees of freedom for the splines.

get.lfc

Logical scalar indicating whether to return an overall log-fold change along each path.

get.spline.coef

Logical scalar indicating whether to return the estimates of the spline coefficients.

trend.only

Deprecated and ignored.

block

Factor of length equal to the number of cells in x, specifying the blocking factor.

row.data

A DataFrame with the same number and order of rows in x, containing per-gene annotations to be cbinded to the output DataFrame(s).

BPPARAM

A BiocParallelParam object from the BiocParallel package, used to control parallelization.

assay.type

String or integer scalar specifying the assay containing the log-expression matrix.

Details

This function fits a natural spline to the expression of each gene with respect to pseudotime. It then does an ANOVA to test whether any of the spline coefficients are non-zero. In this manner, genes exhibiting a significant (and potentially non-linear) trend with respect to the pseudotime can be detected as those with low p-values.

Branched trajectories with multiple paths are represented by a 2-dimensional pseudotime. In this case, only one path is tested at a time by only using one column of pseudotime to form the spline basis matrix. Cells with NA values in any given pseudotime column are assumed to be assigned to a different path and are ignored when fitting the corresponding model.

By default, estimates of the spline coefficients are not returned as they are difficult to interpret. Rather, a log-fold change of expression along each path is estimated to provide some indication of the overall magnitude and direction of any change.

block can be used to fit a separate linear model to each of multiple batches, after which the statistics are combined across batches as described in testLinearModel. This avoids potential confounding effects from batch-specific differences in the distribution of cells across pseudotime.

Value

If pseudotime is a vector, a DataFrame is returned containing the statistics for each gene (row), including the p-value and its BH-adjusted equivalent. If get.lfc=TRUE, an overall log-fold change is returned for each path.

If get.spline.coef=TRUE, the estimated spline coefficients are also returned (single path) or the differences in the spline fits to the first path are returned (multiple paths).

If pseudotime is a 2-dimensional object, a list of DataFrames is instead returned. Each DataFrame has the same format as described above and contains test statistics for each column (i.e., lineage) in pseudotime.

Author(s)

Aaron Lun

See Also

orderCells, to generate the pseudotime matrix.

testLinearModel, which performs the tests under the hood.

Examples

y <- matrix(rnorm(10000), ncol=100)
u <- runif(100)
testPseudotime(y, u)

# Handling a blocking factor.
b <- gl(2, 50)
testPseudotime(y, u, block=b)


LTLA/TSCAN documentation built on Aug. 16, 2024, 12:40 p.m.