fastRUVIII: A fast version of the ruv::RUVIII algorithm

Description Usage Arguments Value Author(s) Examples

View source: R/fastRUVIII.R

Description

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

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
fastRUVIII(
  Y,
  M,
  ctl,
  k = NULL,
  eta = NULL,
  svd_k = 50,
  include.intercept = TRUE,
  average = FALSE,
  BPPARAM = SerialParam(),
  BSPARAM = ExactParam(),
  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.

svd_k

If BSPARAM is set to RandomParam or IrlbaParam class from BiocSingular package, then svd_k will be used to used to reduce the computational cost of singular value decomposition. Default to 50.

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.

BPPARAM

A BiocParallelParam class object from the BiocParallel package is used. Default is SerialParam().

BSPARAM

A BiocSingularParam class object from the BiocSingular package is used. Default is ExactParam().

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

1
2
3
4
5
6
7
8
9
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, BSPARAM = BiocSingular::ExactParam())
improved2 = scMerge::fastRUVIII(Y = Y, M = M, ctl = ctl, 
k = 20, BSPARAM = BiocSingular::RandomParam(), svd_k = 50)
old = ruv::RUVIII(Y = Y, M = M, ctl = ctl, k = 20)
all.equal(improved1, old)
all.equal(improved2, old)

scMerge documentation built on Nov. 8, 2020, 7:04 p.m.