forked from microbiome/mia
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathrunNMDS.R
209 lines (196 loc) · 6.19 KB
/
runNMDS.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
#' Perform non-metric MDS on sample-level data
#'
#' Perform non-metric multi-dimensional scaling (nMDS) on samples, based on the
#' data in a \code{SingleCellExperiment} object.
#'
#' @inheritParams getDominant
#' @inheritParams runDPCoA
#'
#' @details
#' For \code{addNMDS} a
#' \code{\link[TreeSummarizedExperiment]{TreeSummarizedExperiment}}
#'
#' @param keep.dist \code{Logical scalar}. Indicates whether the \code{dist}
#' object calculated by \code{FUN} should be stored as \sQuote{dist} attribute
#' of the matrix returned/stored by \code{getNMDS}/ \code{addNMDS}. (Default:
#' \code{FALSE})
#'
#' @param keep_dist Deprecated. Use \code{keep.dist} instead.
#'
#' @param FUN \code{Function} or \code{Character scalar}. A value with a
#' function name returning a \code{\link[stats:dist]{dist}} object
#'
#' @param nmds.fun \code{Character scalar}. A value to choose the scaling
#' implementation, either \dQuote{isoMDS} for
#' \code{\link[MASS:isoMDS]{MASS::isoMDS}} or \dQuote{monoMDS} for
#' \code{\link[vegan:monoMDS]{vegan::monoMDS}}
#'
#' @param nmdsFUN Deprecated. Use \code{nmds.fun} instead.
#'
#' @param name \code{Character scalar}. A name for the column of the
#' \code{colData} where results will be stored. (Default: \code{"NMDS"})
#'
#' @param ... additional arguments to pass to \code{FUN} and
#' \code{nmds.fun}.
#'
#' @param dimred \code{Character scalar} or \code{integer scalar}. Specifies
#' the existing dimensionality reduction results to use.
#'
#' @param ndimred \code{integer vector}. Specifies the dimensions to use if
#' dimred is specified.
#'
#' @param n_dimred Deprecated. Use \code{ndimred} instead.
#'
#' @param exprs_values Deprecated. Use \code{assay.type} instead.
#'
#' @return For \code{getNMDS}, a matrix is returned containing the MDS
#' coordinates for each sample (row) and dimension (column).
#'
#' @details
#' Either \code{\link[MASS:isoMDS]{MASS::isoMDS}} or
#' \code{\link[vegan:monoMDS]{vegan::monoMDS}} are used internally to compute
#' the NMDS components. If you supply a custom \code{FUN}, make sure that
#' the arguments of \code{FUN} and \code{nmds.fun} do not collide.
#'
#' @name runNMDS
#'
#' @seealso
#' \code{\link[MASS:isoMDS]{MASS::isoMDS}},
#' \code{\link[vegan:monoMDS]{vegan::monoMDS}}
#' for NMDS component calculation.
#'
#' \code{\link[scater:plotReducedDim]{plotMDS}}, to quickly visualize the
#' results.
#'
#' @examples
#' # generate some example data
#' mat <- matrix(1:60, nrow = 6)
#' df <- DataFrame(n = c(1:6))
#' tse <- TreeSummarizedExperiment(assays = list(counts = mat),
#' rowData = df)
#' #
#' getNMDS(tse)
#'
#' #
#' data(esophagus)
#' esophagus <- addNMDS(esophagus, FUN = vegan::vegdist, name = "BC")
#' esophagus <- addNMDS(esophagus, FUN = vegan::vegdist, name = "euclidean",
#' method = "euclidean")
#' reducedDims(esophagus)
NULL
.format_nmds_isoMDS <- function(nmds){
ans <- nmds$points
attr(ans, "Stress") <- nmds$stress
ans
}
#' @importFrom vegan scores
.format_nmds_monoMDS <- function(nmds){
ans <- scores(nmds)
colnames(ans) <- NULL
attr(ans, "Stress") <- nmds$stress
ans
}
.format_nmds <- function(nmds, nmds.fun, sample_names){
ans <- switch(
nmds.fun,
"isoMDS" = .format_nmds_isoMDS(nmds),
"monoMDS" = .format_nmds_monoMDS(nmds))
rownames(ans) <- sample_names
ans
}
.get_nmds_args_isoMDS <- function(args){
args[c("maxit","trace","tol","p")]
}
.get_nmds_args_monoMDS <- function(args){
args[c(
"model","threshold","maxit","weakties","stress","scaling","pc",
"smin","sfgrmin","sratmax")]
}
.get_nmds_args <- function(nmds.fun, ...){
args <- list(...)
args <- switch(
nmds.fun,
"isoMDS" = .get_nmds_args_isoMDS(args),
"monoMDS" = .get_nmds_args_monoMDS(args))
args <- args[!vapply(args,is.null,logical(1))]
args
}
#' @importFrom MASS isoMDS
#' @importFrom stats cmdscale
#' @importFrom vegan vegdist monoMDS
.calculate_nmds <- function(
x, FUN = vegdist,
nmds.fun = nmdsFUN,
nmdsFUN = c("isoMDS","monoMDS"),
ncomponents = 2, ntop = 500, subset.row = subset_row,
subset_row = NULL, scale = FALSE, transposed = FALSE,
keep.dist = keep_dist,
keep_dist = FALSE, ...){
nmds.fun <- match.arg(nmds.fun)
nmdsArgs <- .get_nmds_args(nmds.fun, ...)
if(!transposed) {
x <- .get_mat_for_reddim(
x, subset_row = subset.row, ntop = ntop, scale = scale)
}
x <- as.matrix(x)
sample_names <- rownames(x)
sample_dist <- do.call(
FUN, c(list(x), list(...)))
attributes(sample_dist) <- attributes(sample_dist)[c("class","Size")]
y <- cmdscale(sample_dist, k = ncomponents)
ans <- do.call(
nmds.fun,
c(list(sample_dist, y = y, k = ncomponents), nmdsArgs))
ans <- .format_nmds(ans, nmds.fun, sample_names)
if (keep.dist) {
attr(ans,"dist") <- sample_dist
}
ans
}
#' @rdname runNMDS
#' @export
setMethod("getNMDS", "ANY", .calculate_nmds)
#' @rdname runNMDS
#' @importFrom SummarizedExperiment assay
#' @export
setMethod("getNMDS", "SummarizedExperiment",
function(x, ..., assay.type = assay_name, assay_name = exprs_values,
exprs_values = "counts", FUN = vegdist) {
.calculate_nmds(assay(x, assay.type), FUN = FUN, ...)
}
)
#' @rdname runNMDS
#' @export
setMethod("getNMDS", "SingleCellExperiment",
function(x, ..., assay.type = assay_name, assay_name = exprs_values,
exprs_values = "counts", dimred = NULL, ndimred = n_dimred,
n_dimred = NULL, FUN = vegdist){
mat <- .get_mat_from_sce(
x, exprs_values = assay.type, dimred = dimred, n_dimred = ndimred)
getNMDS(mat, transposed = !is.null(dimred), FUN = FUN,...)
}
)
#' @rdname runNMDS
#' @export
#' @aliases getNMDS
calculateNMDS <- function(x,...){
getNMDS(x,...)
}
#' @rdname runNMDS
#' @export
#' @importFrom SingleCellExperiment reducedDim<- altExp
addNMDS <- function(x, ..., altexp = NULL, name = "NMDS") {
if (!is.null(altexp)) {
y <- altExp(x, altexp)
} else {
y <- x
}
reducedDim(x, name) <- getNMDS(y, ...)
x
}
#' @rdname runNMDS
#' @export
#' @aliases addNMDS
runNMDS <- function(x,...){
addNMDS(x,...)
}