runTSNE {scater} | R Documentation |
Perform t-stochastic neighbour embedding (t-SNE) for the cells, based on the data in a SingleCellExperiment object.
runTSNE(object, ncomponents = 2, ntop = 500, feature_set = NULL, exprs_values = "logcounts", scale_features = TRUE, use_dimred = NULL, n_dimred = NULL, rand_seed = NULL, perplexity = min(50, floor(ncol(object)/5)), pca = TRUE, initial_dims = 50, ...)
object |
A SingleCellExperiment object. |
ncomponents |
Numeric scalar indicating the number of t-SNE dimensions to obtain. |
ntop |
Numeric scalar specifying the number of most variable features to use for t-SNE. |
feature_set |
Character vector of row names, a logical vector or a numeric vector of indices indicating a set of features to use for t-SNE.
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 |
Numeric scalar that can be passed to |
perplexity |
Numeric scalar defining the perplexity parameter, see |
pca |
Logical scalar passed to |
initial_dims |
Integer scalar passed to |
... |
Additional arguments to pass to |
The function Rtsne
is used internally to compute the t-SNE.
Note that the algorithm is not deterministic, so different runs of the function will produce differing results.
Users are advised to test multiple random seed, and then use rand_seed
to set a random seed for replicable results.
The value of the perplexity
parameter can have a large effect on the results.
By default, the function will try to provide a reasonable setting, by scaling the perplexity with the number of cells until it reaches a maximum of 50.
However, it is often worthwhile to manually try multiple values to ensure that the conclusions are robust.
Setting use_dimred
allows users to easily perform t-SNE on 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.
A SingleCellExperiment object containing the coordinates of the first ncomponent
t-SNE dimensions for each cell.
This is stored in the "TSNE"
entry of the reducedDims
slot.
Aaron Lun, based on code by Davis McCarthy
L.J.P. van der Maaten. Barnes-Hut-SNE. In Proceedings of the International Conference on Learning Representations, 2013.
## 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 <- runTSNE(example_sce) reducedDimNames(example_sce) head(reducedDim(example_sce))