CATALYST 1.8.7
daFrame
classcluster
: FlowSOM clustering & ConsensusClusterPlus metaclusteringplotMedExprs
: Median marker-expressions by clusterplotClusterExprs
: Marker-densities by clustermergeClusters
: Manual cluster mergingplotClusterHeatmap
: Heatmap of (meta)clustering resultsplotAbundances
: Relative population abundancesrunDR
Most of the pipeline and visualizations presented herein were adapted from Nowicka et al. (2017)’s “CyTOF workflow: differential discovery in high-throughput high-dimensional cytometry datasets”. For the complete workflow, go here.
# load required packages
suppressPackageStartupMessages({
library(CATALYST)
library(flowCore)
library(diffcyt)
library(SummarizedExperiment)
})
PBMC_fs
:flowSet
holding PBMCs samples from 4 patients, each containing between 500 and 1000 cells. For each sample, the expression of 10 cell surface and 14 signaling markers was measured before (REF) and upon BCR/FcR-XL stimulation (BCRXL) with B cell receptor/Fc receptor crosslinking for 30’, resulting in a total of 8 samples.PBMC_panel
:fcs_colname
column), its targeted protein marker (antigen
column), and the marker_class
(“type” or “state”).PBMC_md
:file_name
, sample_id
, condition
, and patient_id
.# load example data
data(PBMC_fs, PBMC_panel, PBMC_md)
PBMC_fs
## A flowSet with 8 experiments.
##
## column names:
## CD3(110:114)Dd CD45(In115)Dd pNFkB(Nd142)Dd pp38(Nd144)Dd CD4(Nd145)Dd CD20(Sm147)Dd CD33(Nd148)Dd pStat5(Nd150)Dd CD123(Eu151)Dd pAkt(Sm152)Dd pStat1(Eu153)Dd pSHP2(Sm154)Dd pZap70(Gd156)Dd pStat3(Gd158)Dd CD14(Gd160)Dd pSlp76(Dy164)Dd pBtk(Er166)Dd pPlcg2(Er167)Dd pErk(Er168)Dd pLat(Er170)Dd IgM(Yb171)Dd pS6(Yb172)Dd HLA-DR(Yb174)Dd CD7(Yb176)Dd
head(PBMC_panel)
## # A tibble: 6 x 3
## fcs_colname antigen marker_class
## <chr> <chr> <chr>
## 1 CD3(110:114)Dd CD3 type
## 2 CD45(In115)Dd CD45 type
## 3 pNFkB(Nd142)Dd pNFkB state
## 4 pp38(Nd144)Dd pp38 state
## 5 CD4(Nd145)Dd CD4 type
## 6 CD20(Sm147)Dd CD20 type
head(PBMC_md)
## file_name sample_id condition patient_id
## 1 PBMC_patient1_BCRXL.fcs BCRXL1 BCRXL Patient1
## 2 PBMC_patient1_Ref.fcs Ref1 Ref Patient1
## 3 PBMC_patient2_BCRXL.fcs BCRXL2 BCRXL Patient2
## 4 PBMC_patient2_Ref.fcs Ref2 Ref Patient2
## 5 PBMC_patient3_BCRXL.fcs BCRXL3 BCRXL Patient3
## 6 PBMC_patient3_Ref.fcs Ref3 Ref Patient3
The code snippet below demonstrates how to construct a flowSet
from a set of FCS files. However, we also give the option to directly specify the path to a set of FCS files (see next section).
# download exemplary set of FCS files
url <- "http://imlspenticton.uzh.ch/robinson_lab/cytofWorkflow"
fcs_zip <- "PBMC8_fcs_files.zip"
download.file(paste0(url, "/", fcs_zip), destfile = fcs_zip, mode = "wb")
unzip(fcs_zip)
# read in FCS files as flowSet
fcs_files <- list.files(pattern = ".fcs$")
fs <- read.flowSet(fcs_files, transformation = FALSE, truncate_max_range = FALSE)
daFrame
classData used and returned throughout differential analysis are held in objects of the daFrame
class. Its constructor requires the following inputs:
x
: a flowSet
holding the raw measurement data, or a character string that specifies a path to a set of FCS files.panel
: a 2 column data.frame that contains for each marker of interest i) its column name in the raw input data, and ii) its targeted protein marker.md
: a data.frame with columns describing the experimental design.Optionally, cols_to_use
will specify which columns (channels) to keep from the input data. Here, we keep all measurement parameters (default value cols_to_use = NULL
).
# construct daFrame
(daf <- daFrame(PBMC_fs, PBMC_panel, PBMC_md))
## class: daFrame
## dim: 5428 24
## metadata(2): experiment_info cofactor
## assays(1): exprs
## rownames: NULL
## rowData names(3): sample_id condition patient_id
## colnames(24): CD3 CD45 ... HLA_DR CD7
## colData names(3): channel_name marker_name marker_class
## reducedDimNames(0):
## spikeNames(0):
We provide flexibility in the way the panel and metadata table can be set up. Specifically, column names are allowed to differ from the example above, and multiple factors (patient ID, conditions, batch etc.) can be specified. Arguments panel_cols
and md_cols
should then be used to specify which columns hold the required information. An example is given below:
# alter panel column names
panel2 <- PBMC_panel
colnames(panel2)[1:2] <- c("channel_name", "marker")
# alter metadata column names & add 2nd condition
md2 <- PBMC_md
colnames(md2) <- c("file", "sampleID", "cond1", "patientID")
md2$cond2 <- rep(c("A", "B"), 4)
# construct daFrame
daFrame(PBMC_fs, panel2, md2,
panel_cols = list(channel = "channel_name", antigen = "marker"),
md_cols = list(file = "file", id = "sampleID",
factors = c("cond1", "cond2", "patientID")))
Note that, independent of the input panel and metadata tables, the constructor will fix the names of mandatory slots for latter data accession (sample_id
in the rowData
, channel_name
and marker_name
in the colData
). The md
table will be stored under experiment_info
inside the metadata
.
plotCounts
: Number of cells measured per sampleThe number of cells measured per sample may be plotted with plotCounts
, or directly accessed via n_cells()
. This plot should be used as a guide together with other readouts to identify samples where not enough cells were assayed.
n_cells(daf)
## BCRXL1 Ref1 BCRXL2 Ref2 BCRXL3 Ref3 BCRXL4 Ref4
## 528 881 665 438 563 660 934 759
plotCounts(daf, color_by = "condition")
plotMDS
: Multi-dimensional scaling plotA multi-dimensional scaling (MDS) plot on median expresion values may be rendered with plotMDS
. Such a plot will give a sense of similarities between samples in an unsupervised way and of key difference in expression before conducting any formal testing. In our example, we can see a clear separation between reference (REF) and stimulation condition (BCRXL).
plotMDS(daf, color_by = "condition")
plotExprHeatmap
: Heatmap of (scaled) median marker expressionsplotExprHeatmap
will show a heatmap on median marker intensities with hierarchically clustered rows (samples) and columns (markers). This plot should give an idea of which markers will drive sample clustering, and how similiar samples are in their expression profile. We specify bin_anno = TRUE
to display expression values inside each bin, and row_anno = TRUE
to include row annotations for each factor in metadata(daf)
.
plotExprHeatmap(daf, bin_anno = TRUE, row_anno = TRUE)
cluster
: FlowSOM clustering & ConsensusClusterPlus metaclusteringCATALYST provides a simple wrapper to perform high resolution FlowSOM
clustering and lower resolution ConsensusClusterPlus
metaclustering. By default, the data will be initially clustered into xdim = 10
x ydim = 10
= 100 groups. Secondly, the function will metacluster populations into 2 through maxK
(default 20) clusters. To make analyses reproducible, the random seed may be set via seed
. By default, if the colData()$marker_class
column is specified, the set of markers with marker class “type” will be used for clustering. Alternatively, the markers that should be used for clustering can be specified with argument cols_to_use
.
# specify markers to use for clustering
lineage_markers <- c("CD3", "CD45", "CD4", "CD20",
"CD33", "CD123", "CD14", "IgM", "HLA_DR", "CD7")
daf <- cluster(daf, cols_to_use = lineage_markers,
xdim = 10, ydim = 10, maxK = 20, verbose = FALSE, seed = 1)
Let K = xdim
x ydim
be the number of FlowSOM clusters. cluster
will add information to the following slots of the input daFrame
:
rowData
:
cluster_id
: cluster ID as inferred by FlowSOM. One of 1, …, K.colData
:
marker_class
: factor "type"
or "state"
. Specifyies whether a marker has been used for clustering or not, respectively.metadata
:
SOM_codes
: a table with dimensions K x (# type markers). Contains the SOM codes.cluster_codes
: a table with dimensions K x (maxK
+ 1). Contains the cluster codes for all metaclusterings.delta_area
: a ggplot
object (see below for details).The delta area represents the amount of extra cluster stability gained when clustering into k groups as compared to k-1 groups. It can be expected that high stability of clusters can be reached when clustering into the number of groups that best fits the data. The “natural” number of clusters present in the data should thus corresponds to the value of k where there is no longer a considerable increase in stability (pleateau onset). For more details, the user can refer to the original description of the consensus clustering method (Monti et al. 2003).
# access & render delta area plot
metadata(daf)$delta_area
plotMedExprs
: Median marker-expressions by clusterA combined boxplot and jitter of median marker intensitied can be generated via plotMedExprs
. In order to compare medians for each cluster, and potentially identify changes across conditions early on, we specify facet = "cluster_id"
.
p <- plotMedExprs(daf, k = "meta8", facet = "cluster_id",
group_by = "condition", shape_by = "patient_id")
p$facet$params$ncol <- 4
p
Alternatively, we can facet the above plot by antigen
in order to compare marker expressions calculated over all cells across conditions:
p <- plotMedExprs(daf, facet = "antigen",
group_by = "condition", shape_by = "patient_id")
p$facet$params$ncol <- 6
p
plotClusterExprs
: Marker-densities by clusterDistributions of marker intensities (arcsinh-transformed) across cell populations of interest can be plotted with plotClusterExprs
. We specify markers = "type"
(equivalent to type_markers(daf)
), to include type-markers only. Here, blue densities are calculated over all cells and serve as a reference.
plotClusterExprs(daf, k = "meta8", markers = "type")
mergeClusters
: Manual cluster mergingProvided with a 2 column data.frame containing old_cluster
and new_cluster
IDs, mergeClusters
allows for manual cluster merging of any clustering available within the input daFrame
(i.e. the xdim
x ydim
FlowSOM clusters, and any of the 2-maxK
ConsensusClusterPlus metaclusters). For latter accession (visualization, differential testing), the function will assign a unique ID (specified with id
) to each merging, and add a column to the cluster_codes
inside the metadata
slot of the input daFrame
.
data(merging_table)
head(merging_table)
## # A tibble: 6 x 2
## old_cluster new_cluster
## <dbl> <chr>
## 1 1 B-cells IgM+
## 2 2 surface-
## 3 3 NK cells
## 4 4 CD8 T-cells
## 5 5 B-cells IgM-
## 6 6 monocytes
daf <- mergeClusters(daf, k = "meta20", table = merging_table, id = "merging1")
head(cluster_codes(daf))
## som100 meta2 meta3 meta4 meta5 meta6 meta7 meta8 meta9 meta10 meta11
## 1 1 1 1 1 1 1 1 1 1 1 1
## 2 2 1 1 1 1 1 1 1 1 1 1
## 3 3 1 1 1 1 1 1 1 1 2 2
## 4 4 1 1 1 1 1 1 1 1 2 2
## 5 5 1 1 1 1 1 1 1 1 2 2
## 6 6 1 2 2 2 2 2 2 2 3 3
## meta12 meta13 meta14 meta15 meta16 meta17 meta18 meta19 meta20
## 1 1 1 1 1 1 1 1 1 1
## 2 1 1 1 1 1 1 1 1 1
## 3 2 2 2 2 2 2 2 2 2
## 4 2 2 2 2 2 2 2 2 2
## 5 2 2 2 2 2 2 2 2 2
## 6 3 3 3 3 3 3 3 3 3
## merging1
## 1 B-cells IgM+
## 2 B-cells IgM+
## 3 surface-
## 4 surface-
## 5 surface-
## 6 NK cells
plotClusterHeatmap
: Heatmap of (meta)clustering resultsClusterings and metaclusters maybe be viewing with the plotClusterHeatmap
. In its 1st panel, the function will display median (arcsinh-transformed and optionally scaled) cell-type marker expressions (across all samples). Depending on argument hm2
, the 2nd panel will vary as follows:
"abundances"
: cluster frequencies by sample;"state_markers"
: median cell state marker expressions across clusters (analogous to the left-hand side heatmap);Argument scale
(default TRUE
) specifies whether scaled values should be plotted. These correspond to arcsinh-transformed expression values scaled between 0 and 1 using low (1%) and high (99%) percentiles as boundaries. Note that, in any case, hierarchical clustering is performed on the unscaled data.
While it is only used here for visualization, this additional transformation of the arcsinh-transformed data can sometimes give a better representation of relative differences in marker expression between cell populations.
# median pS6 expression by sample as 2nd heatmap
plotClusterHeatmap(daf, hm2 = "pS6", k = "meta12", m = "meta6")
# population frequencies by sample as 2nd heatmap
plotClusterHeatmap(daf, hm2 = "abundances",
draw_freqs = TRUE, cluster_anno = FALSE)
plotAbundances
: Relative population abundancesRelative population abundances for any clustering of interest can be plotted with plotAbundances
. Argument by
will specify whether to plot proportions for each sample or cluster.
If by = "sample_id"
, the function displays each sample’s cell type composition, and the size of a given stripe reflects the proportion of the corresponding cell type the given sample. Argument group
then specifies the facetting. If by = "cluster_id"
, argument group
then specifies the grouping and color coding.
plotAbundances(daf, k = "meta12", by = "sample_id", group_by = "condition")
plotAbundances(daf, k = "merging1", by = "cluster_id",
group_by = "condition", shape_by = "patient_id")
runDR
The set of dimension reduction (DR) methods available from the scater can be applied using runDR
. To make results reproducible, the random seed should be set via set.seed
prior to running runDR
. The subset of markers to use is specified via cols_to_use
, and the cells to use with argument rows_to_use
. When rows_to_use
is a single numeric value N, runDR
will draw a random subset of N cells per sample.
CATALYST
allows DRs to use different sets of cells. In the example below, we can run a PCA on the full dataset, but more time-consuming methods, e.g., t-SNE or UMAP, on only a subset of cells.
# run PCA on all cells
set.seed(1)
daf <- runDR(daf, "PCA")
# run UMAP on 200 cells per sample
set.seed(2)
daf <- runDR(daf, "UMAP", rows_to_use = 200)
As in a regular SingleCellExperiment, DRs available within the daFrame
can be viewed via reducedDimNames
and accessed with reducedDim(s)
.
# view & access DRs
reducedDimNames(daf)
## [1] "PCA" "UMAP"
reducedDims(daf)
## List of length 2
## names(2): PCA UMAP
head(reducedDim(daf, "PCA"))
## PCA1 PCA2
## [1,] -0.8842050 0.1580524
## [2,] -0.2409804 -1.3562566
## [3,] 3.5161452 -0.6373593
## [4,] -0.5935625 -1.0838852
## [5,] -1.3820122 -0.6742902
## [6,] 1.4908868 0.7474590
# all cells
nrow(reducedDim(daf, "PCA")) == nrow(daf)
## [1] TRUE
# 200 per sample
nrow(reducedDim(daf, "UMAP")) == 200 * nlevels(sample_ids(daf))
## [1] TRUE
plotDR
: Visualizing reduced dimensionsplotDR
is a flexbile wrapper to plot reduced dimensions stored within the reducedDims
slot of a daFrame
. The coloring of cell projections is specified with argument color_by
, and may correspond to a marker, clustering, or any factor available in the daFrame
’s rowData
. In order to facilitate comparing reduced dimensions across experimental groups or samples, embeddings can be split via facet
.
# color by marker expression & split by condition
plotDR(daf, "UMAP", color_by = "pS6", facet = "condition")
# color by 8 metaclusters & split by sample ID
plotDR(daf, "UMAP", color_by = "meta8", facet = "sample_id")
CATALYST has been designed to be compatible with the diffcyt package (Weber 2018), which implements statistical methods for differential discovery in high-dimensional cytometry (including flow cytometry, mass cytometry or CyTOF, and oligonucleotide-tagged cytometry) using high-resolution clustering and moderated tests. The input to the diffcyt pipeline can either be raw data, or a daFrame
object. We give an exmaple of the latter below.
Please refer to the diffcyt vignette and R documentation (??diffcyt
) for more detailed information.
# create design & constrast matrix
formula <- createFormula(
experiment_info = PBMC_md,
cols_fixed = "condition",
cols_random = "patient_id")
contrast <- createContrast(c(0, 1))
# test for
# - differential abundance (DA) of clusters
# - differential states (DS) within clusters
res_DA <- diffcyt(daf, contrast = contrast, formula = formula,
analysis_type = "DA", method_DA = "diffcyt-DA-GLMM", clustering_to_use = "meta10")
res_DS <- diffcyt(daf, contrast = contrast, formula = formula,
analysis_type = "DS", method_DS = "diffcyt-DS-LMM", clustering_to_use = "meta10")
plotDiffHeatmap
: Heatmap of differential testing resultsDifferential testing results returned by diffcyt can be displayed with the plotDiffHeatmap
function.
For differential abundance (DA) tests, plotDiffHeatmap
will display
For differential state (DS) tests, plotDiffHeatmap
will display
Clusters (DA) and cluster-marker combinations (DS), respectively, will be marked as significant if their adjusted p-value falls below the threshold value specified with th
(default 0.1), and will be ordered by significance if order = TRUE
(the default). The number of top findings to display can be specified with top_n
(default 20). When normalize = TRUE
, the right-hand side heatmap will display Z-score normalized values. For DA, cluster frequencies will be arcsine-square-root scaled prior to normalization.
plotDiffHeatmap(daf, res_DA, all = TRUE, th = 0.05)
plotDiffHeatmap(daf, res_DS, top_n = 20)
sessionInfo()
## R version 3.6.1 (2019-07-05)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.3 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] SummarizedExperiment_1.14.1 DelayedArray_0.10.0
## [3] BiocParallel_1.18.1 matrixStats_0.55.0
## [5] Biobase_2.44.0 GenomicRanges_1.36.1
## [7] GenomeInfoDb_1.20.0 IRanges_2.18.3
## [9] S4Vectors_0.22.1 BiocGenerics_0.30.0
## [11] diffcyt_1.4.4 flowCore_1.50.0
## [13] CATALYST_1.8.7 BiocStyle_2.12.0
##
## loaded via a namespace (and not attached):
## [1] readxl_1.3.1 backports_1.1.5
## [3] circlize_0.4.8 drc_3.0-1
## [5] plyr_1.8.4 igraph_1.2.4.1
## [7] ConsensusClusterPlus_1.48.0 lazyeval_0.2.2
## [9] shinydashboard_0.7.1 splines_3.6.1
## [11] ggplot2_3.2.1 scater_1.12.2
## [13] TH.data_1.0-10 digest_0.6.21
## [15] htmltools_0.4.0 viridis_0.5.1
## [17] fansi_0.4.0 magrittr_1.5
## [19] cluster_2.1.0 openxlsx_4.1.0.1
## [21] limma_3.40.6 ComplexHeatmap_2.0.0
## [23] RcppParallel_4.4.4 sandwich_2.5-1
## [25] colorspace_1.4-1 rrcov_1.4-7
## [27] ggrepel_0.8.1 haven_2.1.1
## [29] xfun_0.10 dplyr_0.8.3
## [31] crayon_1.3.4 RCurl_1.95-4.12
## [33] jsonlite_1.6 graph_1.62.0
## [35] lme4_1.1-21 zeallot_0.1.0
## [37] survival_2.44-1.1 zoo_1.8-6
## [39] glue_1.3.1 gtable_0.3.0
## [41] nnls_1.4 zlibbioc_1.30.0
## [43] XVector_0.24.0 GetoptLong_0.1.7
## [45] car_3.0-3 BiocSingular_1.0.0
## [47] shape_1.4.4 SingleCellExperiment_1.6.0
## [49] DEoptimR_1.0-8 abind_1.4-5
## [51] scales_1.0.0 mvtnorm_1.0-11
## [53] edgeR_3.26.8 Rcpp_1.0.2
## [55] plotrix_3.7-6 viridisLite_0.3.0
## [57] xtable_1.8-4 clue_0.3-57
## [59] foreign_0.8-72 rsvd_1.0.2
## [61] FlowSOM_1.16.0 tsne_0.1-3
## [63] DT_0.9 htmlwidgets_1.5.1
## [65] httr_1.4.1 FNN_1.1.3
## [67] RColorBrewer_1.1-2 ellipsis_0.3.0
## [69] pkgconfig_2.0.3 XML_3.98-1.20
## [71] uwot_0.1.4 utf8_1.1.4
## [73] locfit_1.5-9.1 labeling_0.3
## [75] tidyselect_0.2.5 rlang_0.4.0
## [77] reshape2_1.4.3 later_1.0.0
## [79] munsell_0.5.0 cellranger_1.1.0
## [81] tools_3.6.1 cli_1.1.0
## [83] ggridges_0.5.1 shinyBS_0.61
## [85] evaluate_0.14 stringr_1.4.0
## [87] fastmap_1.0.1 yaml_2.2.0
## [89] knitr_1.25 zip_2.0.4
## [91] robustbase_0.93-5 purrr_0.3.2
## [93] nlme_3.1-141 mime_0.7
## [95] compiler_3.6.1 beeswarm_0.2.3
## [97] plotly_4.9.0 curl_4.2
## [99] png_0.1-7 tibble_2.1.3
## [101] pcaPP_1.9-73 stringi_1.4.3
## [103] RSpectra_0.15-0 forcats_0.4.0
## [105] lattice_0.20-38 Matrix_1.2-17
## [107] nloptr_1.2.1 shinyjs_1.0
## [109] vctrs_0.2.0 pillar_1.4.2
## [111] lifecycle_0.1.0 BiocManager_1.30.7
## [113] GlobalOptions_0.1.1 BiocNeighbors_1.2.0
## [115] data.table_1.12.4 bitops_1.0-6
## [117] irlba_2.3.3 corpcor_1.6.9
## [119] httpuv_1.5.2 R6_2.4.0
## [121] bookdown_0.14 promises_1.1.0
## [123] gridExtra_2.3 rio_0.5.16
## [125] vipor_0.4.5 codetools_0.2-16
## [127] boot_1.3-23 MASS_7.3-51.4
## [129] gtools_3.8.1 assertthat_0.2.1
## [131] rjson_0.2.20 withr_2.1.2
## [133] multcomp_1.4-10 GenomeInfoDbData_1.2.1
## [135] hms_0.5.1 grid_3.6.1
## [137] minqa_1.2.4 tidyr_1.0.0
## [139] rmarkdown_1.16 DelayedMatrixStats_1.6.1
## [141] carData_3.0-2 Rtsne_0.15
## [143] shiny_1.4.0 ggbeeswarm_0.6.0
Bodenmiller, Bernd, Eli R Zunder, Rachel Finck, Tiffany J Chen, Erica S Savig, Robert V Bruggner, Erin F Simonds, et al. 2012. “Multiplexed Mass Cytometry Profiling of Cellular States Perturbed by Small-Molecule Regulators.” Nature Biotechnology 30 (9):858–67. https://doi.org/10.1038/nbt.2317.
Bruggner, Robert V, Bernd Bodenmiller, David L Dill, Robert J Tibshirani, and Garry P Nolan. 2014. “Automated Identification of Stratifying Signatures in Cellular Subpopulations.” Proceedings of the National Academy of Sciences of the United States of America 111 (26). National Academy of Sciences:E2770–E2777. https://doi.org/10.1073/pnas.1408792111.
Monti, Stefano, Pablo Tamayo, Jill Mesirov, and Todd Golub. 2003. “Consensus Clustering: A Resampling-Based Method for Class Discovery and Visualization of Gene Expression Microarray Data.” Machine Learning 52 (1):91–118. https://doi.org/10.1023/A:1023949509487.
Nowicka, M, C Krieg, LM Weber, FJ Hartmann, S Guglietta, B Becher, MP Levesque, and MD Robinson. 2017. “CyTOF Workflow: Differential Discovery in High-Throughput High-Dimensional Cytometry Datasets [Version 2; Referees: 2 Approved].” F1000Research 6 (748). https://doi.org/10.12688/f1000research.11622.2.
Weber, LM. 2018. “Diffcyt: Differential Discovery in High-Dimensional Cytometry via High-Resolution Clustering.” R Package Version 1.0.0. https://github.com/lmweber/diffcyt.