fastRUVIII {scMerge}R Documentation

A fast version of the ruv::RUVIII algorithm

Description

Perform a fast version of the ruv::RUVIII algorithm for scRNA-Seq data noise estimation

Usage

fastRUVIII(Y, M, ctl, k = NULL, eta = NULL, fast_svd = FALSE,
  rsvd_prop = 0.1, include.intercept = TRUE, average = FALSE,
  fullalpha = NULL, return.info = FALSE, inputcheck = TRUE)

Arguments

Y

The unnormalised scRNA-Seq data matrix. A m by n matrix, where m is the number of observations and n is the number of features.

M

The replicate mapping matrix. The mapping matrix has m rows (one for each observation), and each column represents a set of replicates. The (i, j)-th entry of the mapping matrix is 1 if the i-th observation is in replicate set j, and 0 otherwise. See ruv::RUVIII for more details.

ctl

An index vector to specify the negative controls. Either a logical vector of length n or a vector of integers.

k

The number of unwanted factors to remove. This is inherited from the ruvK argument from the scMerge::scMerge function.

eta

Gene-wise (as opposed to sample-wise) covariates. See ruv::RUVIII for details.

fast_svd

If TRUE, fast algorithms will be used for singular value decomposition calculation via the irlba and rsvd packages. We recommend using this option when the number of cells is large (e.g. more than 1000 cells).

rsvd_prop

If fast_svd = TRUE, then rsvd_prop will be used to used to reduce the computational cost of randomised singular value decomposition. We recommend setting this number to less than 0.25 to achieve a balance between numerical accuracy and computational costs.

include.intercept

When eta is specified (not NULL) but does not already include an intercept term, this will automatically include one. See ruv::RUVIII for details.

average

Average replicates after adjustment. See ruv::RUVIII for details.

fullalpha

Not used. Please ignore. See ruv::RUVIII for details.

return.info

Additional information relating to the computation of normalised matrix. We recommend setting this to true.

inputcheck

We recommend setting this to true.

Value

A normalised matrix of the same dimensions as the input matrix Y.

Author(s)

Yingxin Lin, John Ormerod, Kevin Wang

Examples

L = ruvSimulate(m = 200, n = 500, nc = 400, nCelltypes = 3, nBatch = 2, lambda = 0.1, sce = FALSE)
Y = L$Y; M = L$M; ctl = L$ctl
improved1 = scMerge::fastRUVIII(Y = Y, M = M, ctl = ctl, k = 20, fast_svd = FALSE)
improved2 = scMerge::fastRUVIII(Y = Y, M = M, ctl = ctl, k = 20, fast_svd = TRUE, rsvd_prop = 0.1)
old = ruv::RUVIII(Y = Y, M = M, ctl = ctl, k = 20)
all.equal(improved1, old)
all.equal(improved2, old)

[Package scMerge version 1.0.0 Index]