smoothclust 1.9.6
smoothclust is a method for identification of spatial domains and spatially-aware clustering in spatial transcriptomics data. The method generates spatial domains with smooth boundaries by smoothing gene expression profiles across neighboring spatial locations, followed by unsupervised clustering. Spatial domains consisting of consistent mixtures of cell types may then be further investigated by applying cell type compositional analyses or differential analyses.
The smoothclust package can be installed from Bioconductor as follows (using R version 4.4 onwards). This is the recommended installation for most users. Additional details are shown on the Bioconductor package page.
install.packages("BiocManager")
BiocManager::install("smoothclust")
The latest development version of the package can also be installed from the devel version of Bioconductor or from GitHub.
Input data can be provided either as a SpatialExperiment object within the Bioconductor framework, or as numeric matrices of expression values and spatial coordinates. See help file (?smoothclust) for details.
For matrix input, coordinate rows must correspond to expression columns in the same order.
In the example workflow below, we assume the input is in SpatialExperiment format.
The example workflow in this section demonstrates how to run smoothclust to generate spatial domains with smooth boundaries for a dataset from the 10x Genomics Visium platform.
library(smoothclust)
library(STexampleData)
library(scran)
library(scater)
library(ggspavis)
Load dataset from STexampleData package:
# load data
spe <- Visium_humanDLPFC()
## see ?STexampleData and browseVignettes('STexampleData') for documentation
## loading from cache
# keep spots over tissue
spe <- spe[, colData(spe)$in_tissue == 1]
dim(spe)
## [1] 33538 3639
assayNames(spe)
## [1] "counts"
Filtering is performed outside smoothclust, using the raw counts. First, we
identify mitochondrial genes and calculate spot-level quality-control (QC)
metrics from the complete count matrix. QC thresholds should be selected for
each dataset or sample; for this compact example, we retain tissue spots with
positive total and non-mitochondrial counts.
# calculate QC metrics from complete raw counts
raw_counts <- assay(spe, "counts")
is_mito <- grepl("(^mt-)", rowData(spe)$gene_name, ignore.case = TRUE)
total_counts <- Matrix::colSums(raw_counts)
detected_genes <- Matrix::colSums(raw_counts > 0)
mito_counts <- Matrix::colSums(raw_counts[is_mito, , drop = FALSE])
non_mito_counts <- total_counts - mito_counts
colData(spe)$sum <- total_counts
colData(spe)$detected <- detected_genes
colData(spe)$percent_mito <- 100 * mito_counts / total_counts
summary(colData(spe)[, c("sum", "detected", "percent_mito")])
## [1] "DataFrame object of length 3 with 0 metadata columns"
table(is_mito)
## is_mito
## FALSE TRUE
## 33525 13
# illustrative spot-level QC
keep_spot <- total_counts > 0 & non_mito_counts > 0
spe <- spe[, keep_spot]
raw_counts <- raw_counts[, keep_spot, drop = FALSE]
We calculate library size factors from all non-mitochondrial genes before low-expression filtering. This keeps normalization independent of the chosen expression threshold while preventing mitochondrial composition from affecting normalization of the retained genes. The factors are rescaled to have mean one.
# calculate non-mitochondrial size factors before gene filtering
non_mito_library_sizes <- Matrix::colSums(
raw_counts[!is_mito, , drop = FALSE]
)
size_factors <- non_mito_library_sizes / mean(non_mito_library_sizes)
summary(size_factors)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.04361 0.62235 0.89378 1.00000 1.28152 4.14171
Next, we apply an illustrative low-expression filter and remove mitochondrial genes. Filtering rules are study-specific; here, we retain genes detected in at least five spots. A spot may have zero counts across the retained genes after this feature filtering. Such spots can remain in the input because we supply positive size factors calculated from the broader, unfiltered gene set.
# filter low-expressed and mitochondrial genes
keep_expressed <- Matrix::rowSums(raw_counts > 0) >= 5
keep_gene <- keep_expressed & !is_mito
spe <- spe[keep_gene, ]
dim(spe)
## [1] 17998 3639
Run smoothclust using default parameter settings, which have been selected to be appropriate for Visium data from human tissue. By default, smoothclust normalizes the raw counts using library size factors, smooths the normalized values, and then applies a log-transformation. Supplying the factors calculated above ensures that normalization is defined using all non-mitochondrial genes rather than only the filtered analysis features.
The function returns processed expression values. Feature selection, dimensionality reduction, and clustering are performed outside smoothclust().
The method for smoothing can be specified by providing the method argument, with available options uniform, kernel, and knn. All methods include the focal location. For method = "knn", k counts the other nearest neighbors, giving k + 1 locations in each average. Additional arguments can be used to set parameter values, including bandwidth, k, and truncate, depending on the choice of method. For more details, see the function documentation (?smoothclust).
# run smoothclust using default parameters
spe <- smoothclust(spe, size_factors = size_factors)
The normalized, smoothed, and log-transformed expression values are stored in a
new assay named logcounts_smooth:
# check output object
assayNames(spe)
## [1] "counts" "logcounts_smooth"
Run clustering. We use a clustering workflow adapted from standard single-cell workflows. We select highly variable genes (HVGs) from smoothed log-expression by ranking residual variance above a fitted mean–variance trend using modelGeneVar(). These residuals provide empirical feature-ranking scores for clustering. We then perform k-means clustering on the top principal components (PCs) calculated from the selected genes.
We use a relatively small number of clusters for demonstration purposes in this example:
# preprocessing steps for clustering
# select top highly variable genes (HVGs)
dec <- modelGeneVar(spe, assay.type = "logcounts_smooth")
## Warning in .model_gene_var(x = assay(x, i = assay.type), ...): 'modelGeneVar' is deprecated.
## Use 'scrapper::modelGeneVariances' instead.
## See help("Deprecated")
## Warning in fitTrendVar(fm, fv, ...): 'fitTrendVar' is deprecated.
## Use 'scrapper::fitVarianceTrend' instead.
## See help("Deprecated")
## Warning in combineBlocks(collected, method = method, equiweight = equiweight, : 'combineBlocks' is deprecated.
## See help("Deprecated")
top_hvgs <- getTopHVGs(dec, prop = 0.1)
## Warning in getTopHVGs(dec, prop = 0.1): 'getTopHVGs' is deprecated.
## Use 'scrapper::chooseHighlyVariableGenes' instead.
## See help("Deprecated")
length(top_hvgs)
## [1] 877
# dimensionality reduction
# compute PCA on top HVGs
set.seed(123)
spe <- runPCA(
spe,
subset_row = top_hvgs,
exprs_values = "logcounts_smooth"
)
# run k-means clustering
set.seed(123)
k <- 5
clust <- kmeans(reducedDim(spe, "PCA"), centers = k)$cluster
table(clust)
## clust
## 1 2 3 4 5
## 999 1391 512 373 364
colLabels(spe) <- factor(clust)
Plot clusters / spatial domains generated by the smoothclust workflow:
# color palettes
pal8 <- "libd_layer_colors"
pal36 <- unname(palette.colors(36, "Polychrome 36"))
# plot clusters / spatial domains
plotCoords(spe, annotate = "label", pal = pal8)
Plot manually annotated reference labels, which can be used to evaluate the performance of the clustering in this dataset:
# plot reference labels
plotCoords(spe, annotate = "ground_truth", pal = pal8)
Labels supplied to boundary_density() must correspond to coordinate rows in the same order.
boundary_density() uses its own k argument for the number of neighbors, excluding the focal location and chosen independently of smoothing. For example:
boundary_k <- c(2L, 6L, 18L)
vapply(boundary_k, function(k_i) {
boundary_density(spatialCoords(spe), colLabels(spe), k = k_i)$boundary_density
}, numeric(1))
## [1] 0.05770816 0.05683796 0.08955452
sessionInfo()
## R version 4.6.1 (2026-06-24)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.5 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.24-bioc/R/lib/libRblas.so
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0 LAPACK version 3.12.0
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB 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
##
## time zone: America/New_York
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] ggspavis_1.19.1 scater_1.41.2
## [3] ggplot2_4.0.3 scran_1.41.1
## [5] scuttle_1.23.2 STexampleData_1.21.1
## [7] SpatialExperiment_1.23.0 SingleCellExperiment_1.35.2
## [9] SummarizedExperiment_1.43.0 Biobase_2.73.2
## [11] GenomicRanges_1.65.4 Seqinfo_1.3.2
## [13] IRanges_2.47.5 S4Vectors_0.51.10
## [15] MatrixGenerics_1.25.0 matrixStats_1.5.0
## [17] ExperimentHub_3.3.2 AnnotationHub_4.3.2
## [19] BiocFileCache_3.3.0 dbplyr_2.6.0
## [21] BiocGenerics_0.59.12 generics_0.1.4
## [23] smoothclust_1.9.6 BiocStyle_2.41.0
##
## loaded via a namespace (and not attached):
## [1] DBI_1.3.0 gridExtra_2.3.1 httr2_1.3.0
## [4] rlang_1.3.0 magrittr_2.0.5 otel_0.2.0
## [7] compiler_4.6.1 RSQLite_3.53.3 png_0.1-9
## [10] vctrs_0.7.3 pkgconfig_2.0.3 crayon_1.5.3
## [13] fastmap_1.2.0 magick_2.9.1 XVector_0.53.0
## [16] labeling_0.4.3 rmarkdown_2.32 ggbeeswarm_0.7.3
## [19] tinytex_0.61 purrr_1.2.2 bit_4.6.0
## [22] xfun_0.62 bluster_1.23.1 cachem_1.1.0
## [25] beachmat_2.29.3 jsonlite_2.0.0 blob_1.3.0
## [28] DelayedArray_0.39.8 BiocParallel_1.47.0 irlba_2.4.1
## [31] parallel_4.6.1 cluster_2.1.8.3 R6_2.6.1
## [34] bslib_0.12.0 RColorBrewer_1.1-3 limma_3.99.2
## [37] jquerylib_0.1.4 Rcpp_1.1.2 bookdown_0.48
## [40] knitr_1.52 BiocBaseUtils_1.15.1 Matrix_1.7-6
## [43] igraph_2.3.4 tidyselect_1.2.1 viridis_0.6.5
## [46] dichromat_2.0-1 abind_1.4-8 yaml_2.3.12
## [49] codetools_0.2-20 curl_8.1.0 lattice_0.23-1
## [52] tibble_3.3.1 withr_3.0.3 S7_0.2.2
## [55] KEGGREST_1.53.6 evaluate_1.0.5 Biostrings_2.81.9
## [58] pillar_1.11.1 BiocManager_1.30.27 filelock_1.0.3
## [61] BiocVersion_3.24.0 scales_1.4.0 glue_1.8.1
## [64] metapod_1.21.0 tools_4.6.1 BiocNeighbors_2.7.3
## [67] ScaledMatrix_1.21.0 locfit_1.5-9.12 ggside_0.4.1
## [70] grid_4.6.1 AnnotationDbi_1.75.2 edgeR_4.99.10
## [73] beeswarm_0.4.0 BiocSingular_1.29.2 vipor_0.4.7
## [76] cli_3.6.6 rsvd_1.0.6 rappdirs_0.3.4
## [79] viridisLite_0.4.3 S4Arrays_1.13.2 dplyr_1.2.1
## [82] gtable_0.3.6 sass_0.4.10 digest_0.6.39
## [85] ggrepel_0.9.8 SparseArray_1.13.4 dqrng_0.4.1
## [88] rjson_0.2.23 farver_2.1.2 memoise_2.0.1
## [91] htmltools_0.5.9 lifecycle_1.0.5 httr_1.4.9
## [94] statmod_1.5.2 bit64_4.8.6