Perform ClusterDE on a bulk microarray dataset
Pan Liu
pliu3@fredhutch.orgYiwen Yang
yyang5@fredhutch.org27 August 2026
Source:vignettes/ClusterDE-bulk-microarray.Rmd
ClusterDE-bulk-microarray.RmdDownload data
Load glioblastoma (GBM) dataset processed in Seurat object, originally from Verhaak et al. and can be downloaded at https://gdc.cancer.gov/about-data/publications/gbm_exp. Here we subset the Classical and Mesenchymal subtypes as an example.
Generate synthetic null data
Generate 20 null replicates using ClusterDE stable to eliminate randomness of null generation. Note that for normalized bulk microarray data, the distribution family should be “gaussian”.
RNGkind("L'Ecuyer-CMRG")
seed <- 123
null_data <- ClusterDE::constructNull(
gbm_sub,
data_type = "bulk_microarray",
family = "gaussian",
n_rep = 20,
seed = seed
)Sample clustering
Use Consensus Clustering to divide the null samples into two clusters (following the same clustering procedure as the real data).
set.seed(seed)
null_clusters <- lapply(null_data, function(mat) {
mat <- as.matrix(mat)
mat_centered <- mat - rowMeans(mat, na.rm = T)
clusters <- ConsensusClusterPlus::ConsensusClusterPlus(
mat_centered,
maxK = 3,
reps = 100,
pItem = 0.8,
pFeature = 1,
clusterAlg = "hc",
distance = "pearson",
seed = seed,
plot = F,
title = tempdir()
)[[2]]$consensusClass
})DE analysis
Use any established bulk methods to find DE genes in both real and null data. We use Wilcoxon as an example here.
original_pval <- apply(Seurat::GetAssayData(gbm_sub, layer = "counts"), 1, function(x) {
stats <- wilcox.test(x ~ gbm_sub$subtype)
stats$p.value
})
original_pval <- p.adjust(original_pval, method = "BH")
print(paste0("Wilcox found ", sum(original_pval < 0.05), " DEG"))
#> [1] "Wilcox found 667 DEG"
null_pvals <- lapply(seq_along(null_data), function(i) {
null_pval <- apply(null_data[[i]], 1, function(x) {
stats <- wilcox.test(x ~ null_clusters[[i]])
stats$p.value
})
p.adjust(null_pval, method = "BH")
})FDR control
Finally, compare original p-values and null p-values using ClusterDE FDR control, we observed much fewer DEGs comparing to Wilcoxon:
Session information
sessionInfo()
#> R version 4.5.1 (2025-06-13)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.4 LTS
#>
#> Matrix products: default
#> BLAS/LAPACK: FlexiBLAS OPENBLAS; LAPACK version 3.12.0
#>
#> Random number generation:
#> RNG: L'Ecuyer-CMRG
#> Normal: Inversion
#> Sample: Rejection
#>
#> locale:
#> [1] LC_CTYPE=C.UTF-8 LC_NUMERIC=C LC_TIME=C.UTF-8
#> [4] LC_COLLATE=C.UTF-8 LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
#> [7] LC_PAPER=C.UTF-8 LC_NAME=C LC_ADDRESS=C
#> [10] LC_TELEPHONE=C LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C
#>
#> time zone: America/Los_Angeles
#> tzcode source: system (glibc)
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] BiocStyle_2.36.0
#>
#> loaded via a namespace (and not attached):
#> [1] RColorBrewer_1.1-3 jsonlite_2.0.0
#> [3] magrittr_2.0.3 spatstat.utils_3.1-4
#> [5] farver_2.1.2 rmarkdown_2.29
#> [7] fs_1.6.6 ragg_1.4.0
#> [9] vctrs_0.6.5 ROCR_1.0-11
#> [11] memoise_2.0.1 spatstat.explore_3.4-3
#> [13] S4Arrays_1.8.1 htmltools_0.5.8.1
#> [15] SparseArray_1.8.0 sass_0.4.10
#> [17] sctransform_0.4.2 parallelly_1.45.0
#> [19] KernSmooth_2.23-26 bslib_0.9.0
#> [21] htmlwidgets_1.6.4 desc_1.4.3
#> [23] ica_1.0-3 plyr_1.8.9
#> [25] plotly_4.10.4 zoo_1.8-14
#> [27] cachem_1.1.0 igraph_2.1.4
#> [29] mime_0.13 lifecycle_1.0.4
#> [31] pkgconfig_2.0.3 Matrix_1.7-5
#> [33] R6_2.6.1 fastmap_1.2.0
#> [35] GenomeInfoDbData_1.2.14 MatrixGenerics_1.20.0
#> [37] fitdistrplus_1.2-2 future_1.58.0
#> [39] shiny_1.10.0 digest_0.6.37
#> [41] S4Vectors_0.46.0 patchwork_1.3.0
#> [43] Seurat_5.3.0 tensor_1.5
#> [45] RSpectra_0.16-2 irlba_2.3.5.1
#> [47] GenomicRanges_1.60.0 textshaping_1.0.1
#> [49] progressr_0.15.1 spatstat.sparse_3.1-0
#> [51] mgcv_1.9-3 httr_1.4.7
#> [53] polyclip_1.10-7 abind_1.4-8
#> [55] compiler_4.5.1 withr_3.0.2
#> [57] backports_1.5.0 ConsensusClusterPlus_1.72.0
#> [59] fastDummies_1.7.5 MASS_7.3-65
#> [61] DelayedArray_0.34.1 tools_4.5.1
#> [63] lmtest_0.9-40 httpuv_1.6.16
#> [65] future.apply_1.20.0 goftest_1.2-3
#> [67] glue_1.8.0 nlme_3.1-168
#> [69] promises_1.3.3 grid_4.5.1
#> [71] checkmate_2.3.2 Rtsne_0.17
#> [73] cluster_2.1.8.1 reshape2_1.4.4
#> [75] generics_0.1.4 gtable_0.3.6
#> [77] spatstat.data_3.1-6 tidyr_1.3.1
#> [79] data.table_1.17.4 XVector_0.48.0
#> [81] sp_2.2-0 BiocGenerics_0.54.0
#> [83] spatstat.geom_3.4-1 RcppAnnoy_0.0.22
#> [85] ggrepel_0.9.6 RANN_2.6.2
#> [87] pillar_1.10.2 stringr_1.5.1
#> [89] spam_2.11-1 RcppHNSW_0.6.0
#> [91] later_1.4.2 splines_4.5.1
#> [93] dplyr_1.1.4 lattice_0.22-9
#> [95] survival_3.8-3 deldir_2.0-4
#> [97] gamlss.data_6.0-6 tidyselect_1.2.1
#> [99] SingleCellExperiment_1.30.1 miniUI_0.1.2
#> [101] pbapply_1.7-2 knitr_1.50
#> [103] gridExtra_2.3 bookdown_0.43
#> [105] IRanges_2.42.0 SummarizedExperiment_1.38.1
#> [107] scattermore_1.2 stats4_4.5.1
#> [109] xfun_0.52 Biobase_2.68.0
#> [111] matrixStats_1.5.0 UCSC.utils_1.4.0
#> [113] stringi_1.8.7 ClusterDE_0.99.4
#> [115] lazyeval_0.2.2 yaml_2.3.10
#> [117] evaluate_1.0.3 codetools_0.2-20
#> [119] tibble_3.3.0 BiocManager_1.30.27
#> [121] cli_3.6.5 uwot_0.2.3
#> [123] xtable_1.8-4 reticulate_1.42.0
#> [125] systemfonts_1.2.3 jquerylib_0.1.4
#> [127] GenomeInfoDb_1.44.0 dichromat_2.0-0.1
#> [129] Rcpp_1.0.14 globals_0.18.0
#> [131] spatstat.random_3.4-1 scDesign3_1.7.6
#> [133] png_0.1-8 spatstat.univar_3.1-3
#> [135] parallel_4.5.1 pkgdown_2.0.9
#> [137] ggplot2_3.5.2 mclust_6.1.1
#> [139] dotCall64_1.2 bettermc_1.2.2.9000
#> [141] listenv_0.9.1 gamlss.dist_6.1-1
#> [143] mvtnorm_1.3-3 viridisLite_0.4.2
#> [145] gamlss_5.4-22 scales_1.4.0
#> [147] ggridges_0.5.6 crayon_1.5.3
#> [149] SeuratObject_5.1.0 purrr_1.0.4
#> [151] rlang_1.1.6 cowplot_1.1.3