Perform non-metric multi-dimensional scaling (nMDS) on samples, based on the
data in a SingleCellExperiment
object.
calculateNMDS(x, ...)
# S4 method for ANY
calculateNMDS(
x,
FUN = vegdist,
nmdsFUN = c("isoMDS", "monoMDS"),
ncomponents = 2,
ntop = 500,
subset_row = NULL,
scale = FALSE,
transposed = FALSE,
keep_dist = FALSE,
...
)
# S4 method for SummarizedExperiment
calculateNMDS(
x,
...,
assay.type = assay_name,
assay_name = exprs_values,
exprs_values = "counts",
FUN = vegdist
)
# S4 method for SingleCellExperiment
calculateNMDS(
x,
...,
assay.type = assay_name,
assay_name = exprs_values,
exprs_values = "counts",
dimred = NULL,
n_dimred = NULL,
FUN = vegdist
)
runNMDS(x, ..., altexp = NULL, name = "NMDS")
For calculateNMDS
, a numeric matrix of expression values
where rows are features and columns are cells.
Alternatively, a TreeSummarizedExperiment
containing such a matrix.
For runNMDS
a SingleCellExperiment
additional arguments to pass to FUN
and
nmdsFUN
.
a function
or character
value with a function
name returning a dist
object
a character
value to choose the scaling
implementation, either “isoMDS” for
MASS::isoMDS
or “monoMDS” for
vegan::monoMDS
Numeric scalar indicating the number of NMDS dimensions to obtain.
Numeric scalar specifying the number of features with the highest variances to use for dimensionality reduction.
Vector specifying the subset of features to use for dimensionality reduction. This can be a character vector of row names, an integer vector of row indices or a logical vector.
Logical scalar, should the expression values be standardized?
Logical scalar, is x transposed with cells in rows?
Logical scalar indicating whether the dist
object
calculated by FUN
should be stored as ‘dist’ attribute of
the matrix returned/stored by calculateNMDS
/ runNMDS
.
a single character
value for specifying which
assay to use for calculation.
a single character
value for specifying which
assay to use for calculation.
(Please use assay.type
instead. At some point assay_name
will be disabled.)
a single character
value for specifying which
assay to use for calculation.
(Please use assay.type
instead.)
String or integer scalar specifying the existing dimensionality reduction results to use.
Integer scalar or vector specifying the dimensions to use if dimred is specified.
String or integer scalar specifying an alternative experiment containing the input data.
String specifying the name to be used to store the result in the reducedDims of the output.
For calculateNMDS
, a matrix is returned containing the MDS
coordinates for each sample (row) and dimension (column).
Either MASS::isoMDS
or
vegan::monoMDS
are used internally to compute
the NMDS components. If you supply a custom FUN
, make sure that
the arguments of FUN
and nmdsFUN
do not collide.
MASS::isoMDS
,
vegan::monoMDS
for NMDS component calculation.
plotMDS
, to quickly visualize the
results.
# generate some example data
mat <- matrix(1:60, nrow = 6)
df <- DataFrame(n = c(1:6))
tse <- TreeSummarizedExperiment(assays = list(counts = mat),
rowData = df)
#
calculateNMDS(tse)
#> initial value 0.383462
#> iter 5 value 0.161655
#> iter 10 value 0.113278
#> final value 0.003270
#> converged
#> [,1] [,2]
#> [1,] -0.60647893 0.282119175
#> [2,] -0.38097446 -0.130374816
#> [3,] -0.17579928 -0.184259003
#> [4,] -0.02959516 -0.158516980
#> [5,] 0.09082607 -0.111879947
#> [6,] 0.14427170 -0.047910688
#> [7,] 0.20118089 -0.001837781
#> [8,] 0.22817639 0.061022767
#> [9,] 0.25717122 0.116809029
#> [10,] 0.27122156 0.174828245
#> attr(,"Stress")
#> [1] 0.003269857
#
data(esophagus)
esophagus <- runNMDS(esophagus, FUN = vegan::vegdist, name = "BC")
#> initial value 0.000000
#> final value 0.000000
#> converged
esophagus <- runNMDS(esophagus, FUN = vegan::vegdist, name = "euclidean",
method = "euclidean")
#> initial value 0.000000
#> final value 0.000000
#> converged
reducedDims(esophagus)
#> List of length 2
#> names(2): BC euclidean