DuoClustering2018 1.2.0
This vignette describes how each of the included clustering methods was applied to the collection of data sets in order to generate the clustering result summaries provided with the package. It also shows how to apply a new clustering method to the included data sets, to generate results that can be compared to those already included.
The code below describes how we applied each of the included clustering methods
to the data sets for our paper (Duò, Robinson, and Soneson 2018). The apply_*()
functions,
describing how the respective clustering methods were run, are available from
the GitHub
repository
corresponding to the publication. In order to apply a new clustering algorithm
to one of the data sets using the same framework, it is necessary to generate a
function with the same format. The input arguments to this function should be:
SingleCellExperiment
objectThe function should return a list with three elements:
st
- a vector with the timing information. Should have five elements, named
user.self
, sys.self
, user.child
, sys.child
and elapsed
.cluster
- a named vector of cluster assignments for all cells.est_k
- the number of clusters estimated by the method (if available,
otherwise NA
).If the method does not allow specification of the desired number of clusters,
but has another parameter affecting the resolution, this can be accommodated as
well (see the solution for Seurat
in the code below).
First, load the package and define the data set and clustering method to use
(note that in order to apply a method named <method>
, there has to be a
function named apply_<method>()
, with the above specifications, available in
the workspace).
suppressPackageStartupMessages({
library(DuoClustering2018)
})
scename <- "sce_filteredExpr10_Koh"
sce <- sce_filteredExpr10_Koh()
## snapshotDate(): 2019-04-29
## see ?DuoClustering2018 and browseVignettes('DuoClustering2018') for documentation
## downloading 0 resources
## loading from cache
## 'EH1501 : 1501'
method <- "PCAHC"
Next, define the list of hyperparameter values. The package contains the hyperparameter values for the methods included in our paper.
## Load parameter files. General dataset and method parameters as well as
## dataset/method-specific parameters
params <- duo_clustering_all_parameter_settings_v2()[[paste0(scename, "_",
method)]]
## snapshotDate(): 2019-04-29
## see ?DuoClustering2018 and browseVignettes('DuoClustering2018') for documentation
## downloading 0 resources
## loading from cache
## 'EH1619 : 1619'
params
## $nPC
## [1] 30
##
## $range_clusters
## [1] 2 3 4 5 6 7 8 9 10 11 12 13 14 15
Finally, define the number of times to apply the clustering method (for each value of the number of clusters), and run the clustering across a range of imposed numbers of clusters (defined in the parameter list).
## Set number of times to run clustering for each k
n_rep <- 5
## Run clustering
set.seed(1234)
L <- lapply(seq_len(n_rep), function(i) { ## For each run
cat(paste0("run = ", i, "\n"))
if (method == "Seurat") {
tmp <- lapply(params$range_resolutions, function(resolution) {
## For each resolution
cat(paste0("resolution = ", resolution, "\n"))
## Run clustering
res <- get(paste0("apply_", method))(sce = sce, params = params,
resolution = resolution)
## Put output in data frame
df <- data.frame(dataset = scename,
method = method,
cell = names(res$cluster),
run = i,
k = length(unique(res$cluster)),
resolution = resolution,
cluster = res$cluster,
stringsAsFactors = FALSE, row.names = NULL)
tm <- data.frame(dataset = scename,
method = method,
run = i,
k = length(unique(res$cluster)),
resolution = resolution,
user.self = res$st[["user.self"]],
sys.self = res$st[["sys.self"]],
user.child = res$st[["user.child"]],
sys.child = res$st[["sys.child"]],
elapsed = res$st[["elapsed"]],
stringsAsFactors = FALSE, row.names = NULL)
kest <- data.frame(dataset = scename,
method = method,
run = i,
k = length(unique(res$cluster)),
resolution = resolution,
est_k = res$est_k,
stringsAsFactors = FALSE, row.names = NULL)
list(clusters = df, timing = tm, kest = kest)
}) ## End for each resolution
} else {
tmp <- lapply(params$range_clusters, function(k) { ## For each k
cat(paste0("k = ", k, "\n"))
## Run clustering
res <- get(paste0("apply_", method))(sce = sce, params = params, k = k)
## Put output in data frame
df <- data.frame(dataset = scename,
method = method,
cell = names(res$cluster),
run = i,
k = k,
resolution = NA,
cluster = res$cluster,
stringsAsFactors = FALSE, row.names = NULL)
tm <- data.frame(dataset = scename,
method = method,
run = i,
k = k,
resolution = NA,
user.self = res$st[["user.self"]],
sys.self = res$st[["sys.self"]],
user.child = res$st[["user.child"]],
sys.child = res$st[["sys.child"]],
elapsed = res$st[["elapsed"]],
stringsAsFactors = FALSE, row.names = NULL)
kest <- data.frame(dataset = scename,
method = method,
run = i,
k = k,
resolution = NA,
est_k = res$est_k,
stringsAsFactors = FALSE, row.names = NULL)
list(clusters = df, timing = tm, kest = kest)
}) ## End for each k
}
## Summarize across different values of k
assignments <- do.call(rbind, lapply(tmp, function(w) w$clusters))
timings <- do.call(rbind, lapply(tmp, function(w) w$timing))
k_estimates <- do.call(rbind, lapply(tmp, function(w) w$kest))
list(assignments = assignments, timings = timings, k_estimates = k_estimates)
}) ## End for each run
## Summarize across different runs
assignments <- do.call(rbind, lapply(L, function(w) w$assignments))
timings <- do.call(rbind, lapply(L, function(w) w$timings))
k_estimates <- do.call(rbind, lapply(L, function(w) w$k_estimates))
## Add true group for each cell
truth <- data.frame(cell = as.character(rownames(colData(sce))),
trueclass = as.character(colData(sce)$phenoid),
stringsAsFactors = FALSE)
assignments$trueclass <- truth$trueclass[match(assignments$cell, truth$cell)]
## Combine results
res <- list(assignments = assignments, timings = timings,
k_estimates = k_estimates)
df <- dplyr::full_join(res$assignments %>%
dplyr::select(dataset, method, cell, run, k,
resolution, cluster, trueclass),
res$k_estimates %>%
dplyr::select(dataset, method, run, k,
resolution, est_k)
) %>% dplyr::full_join(res$timings %>% dplyr::select(dataset, method, run, k,
resolution, elapsed))
The resulting df
data frames can then be combined across data sets, filterings
and methods and used as input to the provided plotting functions.
sessionInfo()
## R version 3.6.0 (2019-04-26)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.2 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.9-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.9-bioc/R/lib/libRlapack.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=C
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] plyr_1.8.4 ExperimentHub_1.10.0
## [3] AnnotationHub_2.16.0 BiocFileCache_1.8.0
## [5] dbplyr_1.4.0 tidyr_0.8.3
## [7] dplyr_0.8.0.1 DuoClustering2018_1.2.0
## [9] SingleCellExperiment_1.6.0 SummarizedExperiment_1.14.0
## [11] DelayedArray_0.10.0 BiocParallel_1.18.0
## [13] matrixStats_0.54.0 Biobase_2.44.0
## [15] GenomicRanges_1.36.0 GenomeInfoDb_1.20.0
## [17] IRanges_2.18.0 S4Vectors_0.22.0
## [19] BiocGenerics_0.30.0 BiocStyle_2.12.0
##
## loaded via a namespace (and not attached):
## [1] viridis_0.5.1 httr_1.4.0
## [3] bit64_0.9-7 viridisLite_0.3.0
## [5] shiny_1.3.2 assertthat_0.2.1
## [7] interactiveDisplayBase_1.22.0 BiocManager_1.30.4
## [9] blob_1.1.1 GenomeInfoDbData_1.2.1
## [11] yaml_2.2.0 pillar_1.3.1
## [13] RSQLite_2.1.1 lattice_0.20-38
## [15] glue_1.3.1 digest_0.6.18
## [17] promises_1.0.1 XVector_0.24.0
## [19] colorspace_1.4-1 htmltools_0.3.6
## [21] httpuv_1.5.1 Matrix_1.2-17
## [23] pkgconfig_2.0.2 bookdown_0.9
## [25] zlibbioc_1.30.0 xtable_1.8-4
## [27] purrr_0.3.2 scales_1.0.0
## [29] later_0.8.0 tibble_2.1.1
## [31] ggplot2_3.1.1 lazyeval_0.2.2
## [33] mime_0.6 magrittr_1.5
## [35] crayon_1.3.4 mclust_5.4.3
## [37] memoise_1.1.0 evaluate_0.13
## [39] ggthemes_4.1.1 tools_3.6.0
## [41] stringr_1.4.0 munsell_0.5.0
## [43] AnnotationDbi_1.46.0 compiler_3.6.0
## [45] rlang_0.3.4 grid_3.6.0
## [47] RCurl_1.95-4.12 rappdirs_0.3.1
## [49] labeling_0.3 bitops_1.0-6
## [51] rmarkdown_1.12 gtable_0.3.0
## [53] DBI_1.0.0 curl_3.3
## [55] reshape2_1.4.3 R6_2.4.0
## [57] gridExtra_2.3 knitr_1.22
## [59] bit_1.1-14 stringi_1.4.3
## [61] Rcpp_1.0.1 tidyselect_0.2.5
## [63] xfun_0.6
Duò, A, MD Robinson, and D Soneson. 2018. “A Systematic Performance Evaluation of Clustering Methods for Single-Cell RNA-seq Data.” F1000Research 7:1141.