runDiffusionMap {scater} | R Documentation |
Produce a diffusion map for the cells, based on the data in a SingleCellExperiment object.
runDiffusionMap(object, ncomponents = 2, ntop = 500, feature_set = NULL, exprs_values = "logcounts", scale_features = TRUE, use_dimred = NULL, n_dimred = NULL, rand_seed = NULL, ...)
object |
A SingleCellExperiment object |
ncomponents |
Numeric scalar indicating the number of diffusion components to obtain. |
ntop |
Numeric scalar specifying the number of most variable features to use for constructing the diffusion map. |
feature_set |
Character vector of row names, a logical vector or a numeric vector of indices indicating a set of features to use to construct the diffusion map.
This will override any |
exprs_values |
Integer scalar or string indicating which assay of |
scale_features |
Logical scalar, should the expression values be standardised so that each feature has unit variance? |
use_dimred |
String or integer scalar specifying the entry of |
n_dimred |
Integer scalar, number of dimensions of the reduced dimension slot to use when |
rand_seed |
Deprecated, numeric scalar that can be passed to |
... |
Additional arguments to pass to |
The function DiffusionMap
is used internally to compute the diffusion map.
Setting use_dimred
allows users to easily construct a diffusion map from low-rank approximations of the original expression matrix (e.g., after PCA).
In such cases, arguments such as ntop
, feature_set
, exprs_values
and scale_features
will be ignored.
The behaviour of DiffusionMap
seems to be non-deterministic, in a manner that is not responsive to any set.seed
call.
The reason for this is unknown.
A SingleCellExperiment object containing the coordinates of the first ncomponent
diffusion map components for each cell.
This is stored in the "DiffusionMap"
entry of the reducedDims
slot.
Aaron Lun, based on code by Davis McCarthy
Haghverdi L, Buettner F, Theis FJ. Diffusion maps for high-dimensional single-cell analysis of differentiation data. Bioinformatics. 2015; doi:10.1093/bioinformatics/btv325
## Set up an example SingleCellExperiment data("sc_example_counts") data("sc_example_cell_info") example_sce <- SingleCellExperiment( assays = list(counts = sc_example_counts), colData = sc_example_cell_info ) example_sce <- normalize(example_sce) example_sce <- runDiffusionMap(example_sce) reducedDimNames(example_sce) head(reducedDim(example_sce))