Region based spatial transcriptomics analysis with scider
Mengbo Li, Ning Liu, Quoc Hoang Nguyen and Yunshun Chen
Source:vignettes/getting-started.Rmd
getting-started.RmdIntroduction
scider is a user-friendly R package providing functions
to model the global density of cells in a slide of spatial
transcriptomics data. All functions in the package are built on the
SpatialExperiment object, allowing integration into various
spatial transcriptomics-related packages from Bioconductor. After
modelling density, scider supports several downstream
analyses, including colocalization analysis, boundary detection and
differential density analysis.
The plotting functions return ggplot objects, so figures
can be further customised with any ggplot2 layer.
This vignette walks through a typical workflow on a Xenium breast cancer dataset.
Load Example Data
The example dataset is a small SpatialExperiment object
(one Xenium breast cancer section) hosted as a GitHub Release asset.
Download it once and cache locally:
data_url <- "https://github.com/ChenLaboratory/example_data/releases/download/v1.0.0/spe_xenium_bc_s1rep1.rds"
cache_dir <- tools::R_user_dir("scider", "cache")
data_file <- file.path(cache_dir, "spe_xenium_bc_s1rep1.rds")
if (!file.exists(data_file)) {
dir.create(cache_dir, recursive = TRUE, showWarnings = FALSE)
download.file(data_url, data_file, mode = "wb")
}
spe <- readRDS(data_file)
spe
#> class: SpatialExperiment
#> dim: 313 163797
#> metadata(0):
#> assays(1): counts
#> rownames(313): ABCC11 ACTA2 ... ZEB2 ZNF562
#> rowData names(3): ID Symbol Type
#> colnames(163797): cell_1 cell_2 ... cell_167779 cell_167780
#> colData names(12): cell_id transcript_counts ... gene_counts cell_type
#> reducedDimNames(0):
#> mainExpName: NULL
#> altExpNames(0):
#> spatialCoords names(2) : x_centroid y_centroid
#> imgData names(1): sample_idData exploration
The Xenium data was annotated by 10X Genomics. The cell type annotation is summarized below.
table(colData(spe)$cell_type)
#>
#> Stromal Invasive_Tumor DCIS Macrophages Myoepi Endothelial CD4+_T_Cells
#> 41422 38149 24604 12798 9938 8931 8453
#> CD8+_T_Cells Others B_Cells Perivascular-Like DCs Mast_Cells
#> 6940 5769 4987 847 792 167We can visualize the spatial distribution of all the cell types across the tissue as follows.
plotSpatial(spe, group.by = "cell_type")
Or visualize the spatial distribution of specific cell types (e.g., DCIS and Invasive Tumor).
COI <- c("DCIS", "Invasive_Tumor")
sel <- spe$cell_type %in% COI
plotSpatial(spe[, sel], group.by = "cell_type")
We can also examine the expression level of each individual gene
spatially. Because the returned object is a ggplot, we can
add a title with ggtitle().
plotSpatial(spe, feature = "FOXA1", type = "log") +
ggtitle("FOXA1 expression (log-scale)")
Cell type density
scider implements grid-based density analysis for
spatial transcriptomics data. We first calculate density for each cell
type using function gridDensity. The calculated density and
grid information are saved in the metadata of the
SpatialExperiment object. The size of the grid (bin) is set
by the grid.length.x parameter, and the level of smoothing
is set by bandwidth.
spe <- gridDensity(spe, grid.length.x = 50, bandwidth = 50)We can visualise the overall (all cell types) density of the whole
tissue using function plotDensity. The probs
parameter indicates the density percentile below which the grids are
filtered (not plotted). Setting probs to 0 will plot
everything including the background.
plotDensity(spe, probs = 0)
Alternatively, user can specify cutoff to filter grids
with low cell density. The density values represent the expected numbers
of cells at each grid.
plotDensity(spe, cutoff = 5)
We can also visualise the density of specific cell type(s).
plotDensity(spe, coi = COI)
Find Regions-of-interest (ROIs)
After obtaining grid-based density for each cell type, we can then
detect regions-of-interest (ROIs) based on density of certain cell
type(s). This is done by using the findROI function.
Here we identify ROIs based on the cell density of DCIS and invasive
tumor. The min.density parameter filters grids with low
cell density, whereas the ngrid.min gives the minimum
number of grids that each ROI should have.
spe <- findROI(spe, coi = COI, min.density = 0.5, ngrid.min = 20)
#> Overwriting the probs argument. Grids are filtered by the min.density value.The ROIs can then be visualized using plotROI.
plotROI(spe, roi = COI)
Testing relationship between cell types
After defining ROIs, we can test whether the densities of any two
cell types are spatially correlated, within each ROI or overall while
accounting for ROI variation. This is done with function
corDensity, which uses a modified t-test (Dutilleul’s
correction) that accounts for the spatial autocorrelation between
neighbouring grids.
results <- corDensity(spe, roi = COI)We can examine the correlation between each pair of cell types in each ROI.
results$ROI
#> DataFrame with 3276 rows and 9 columns
#> celltype1 celltype2 ROI ngrid cor.coef t df p.Pos p.Neg
#> <character> <character> <character> <numeric> <numeric> <numeric> <numeric> <numeric> <numeric>
#> 1 Stromal Invasive tumor 1 202 -0.151757 -1.16577 57.6512 0.875746 0.124253901
#> 2 Stromal Invasive tumor 2 128 -0.777733 -4.30982 12.1338 0.999505 0.000494573
#> 3 Stromal Invasive tumor 3 248 -0.700439 -3.20855 10.6886 0.995690 0.004310468
#> 4 Stromal Invasive tumor 4 171 -0.457845 -2.54249 24.3734 0.991120 0.008880219
#> 5 Stromal Invasive tumor 5 227 -0.488500 -3.92329 49.1096 0.999864 0.000135806
#> ... ... ... ... ... ... ... ... ... ...
#> 3272 D cs Mast cells 38 54 0.3382772 1.827257 25.8390 0.0396169 0.960383
#> 3273 D cs Mast cells 39 108 0.2618742 1.662013 37.5172 0.0524223 0.947578
#> 3274 D cs Mast cells 40 132 0.0235890 0.157033 44.2912 0.4379664 0.562034
#> 3275 D cs Mast cells 41 103 0.3117195 1.585540 23.3578 0.0631452 0.936855
#> 3276 D cs Mast cells 42 62 -0.0201452 -0.146636 52.9616 0.5580119 0.441988Or the correlation between each pair of cell types across all ROIs:
results$overall
#> DataFrame with 78 rows and 5 columns
#> celltype1 celltype2 cor.coef p.Pos p.Neg
#> <character> <character> <numeric> <numeric> <numeric>
#> 1 Stromal Invasive tumor -0.4230995 1.00000e+00 4.29019e-53
#> 2 Stromal Dcis -0.5108656 1.28106e-01 8.94811e-62
#> 3 Stromal Macrophages 0.3641173 5.08375e-32 9.99999e-01
#> 4 Stromal Myoepi -0.0423066 1.48122e-19 4.37604e-03
#> 5 Stromal Endothelial 0.4622224 1.09692e-55 1.00000e+00
#> ... ... ... ... ... ...
#> 74 B cells D cs 0.3357256 4.33879e-61 1.000000
#> 75 B cells Mast cells 0.1900816 2.09341e-27 0.999998
#> 76 Perivascular like D cs 0.1259144 2.18311e-11 0.999801
#> 77 Perivascular like Mast cells 0.0859226 9.10731e-07 0.952813
#> 78 D cs Mast cells 0.0992215 3.65284e-29 0.998556We can visualise the statistics between each pair of cell types using
function plotCorHeatmap in the ROIs:
plotCorHeatmap(results$ROI)Or the correlation between cell type pairs across the whole slide:
plotCorHeatmap(results$overall)Cell density contour
Another strategy to partition a tissue spatially is to construct cell
density contours. This can be done using the getContour
function.
The equal.cell parameter ensures roughly equal cell
numbers across all density levels.
spe <- getContour(spe, coi = COI, equal.cell = TRUE)
#> Finding contour using total density of input COIs.
#> Using bins = 10 to draw contours with equal cell numbers.The constructed contours can be visualised using
plotContour.
plotContour(spe, coi = COI)
Session Information
sessionInfo()
#> R version 4.6.0 (2026-04-24 ucrt)
#> Platform: x86_64-w64-mingw32/x64
#> Running under: Windows 11 x64 (build 26200)
#>
#> Matrix products: default
#> LAPACK version 3.12.1
#>
#> locale:
#> [1] LC_COLLATE=English_Australia.utf8 LC_CTYPE=English_Australia.utf8 LC_MONETARY=English_Australia.utf8
#> [4] LC_NUMERIC=C LC_TIME=English_Australia.utf8
#>
#> time zone: Australia/Sydney
#> tzcode source: internal
#>
#> attached base packages:
#> [1] stats4 stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] ggplot2_4.0.3 SpatialExperiment_1.21.0 SingleCellExperiment_1.33.2 SummarizedExperiment_1.41.1
#> [5] Biobase_2.71.0 GenomicRanges_1.63.2 Seqinfo_1.1.0 IRanges_2.45.0
#> [9] S4Vectors_0.49.2 BiocGenerics_0.57.1 generics_0.1.4 MatrixGenerics_1.23.0
#> [13] matrixStats_1.5.0 scider_1.9.6
#>
#> loaded via a namespace (and not attached):
#> [1] DBI_1.3.0 deldir_2.0-4 rlang_1.2.0 magrittr_2.0.5 snakecase_0.11.1
#> [6] otel_0.2.0 e1071_1.7-17 compiler_4.6.0 spatstat.geom_3.7-3 mgcv_1.9-4
#> [11] fftwtools_0.9-11 vctrs_0.7.3 stringr_1.6.0 pkgconfig_2.0.3 fastmap_1.2.0
#> [16] magick_2.9.1 XVector_0.51.0 labeling_0.4.3 promises_1.5.0 purrr_1.2.2
#> [21] xfun_0.57 jsonlite_2.0.0 goftest_1.2-3 later_1.4.8 DelayedArray_0.37.1
#> [26] spatstat.utils_3.2-2 R6_2.6.1 stringi_1.8.7 RColorBrewer_1.1-3 spatstat.data_3.1-9
#> [31] spatstat.univar_3.1-7 lubridate_1.9.5 Rcpp_1.1.1-1.1 knitr_1.51 tensor_1.5.1
#> [36] httpuv_1.6.17 Matrix_1.7-5 splines_4.6.0 igraph_2.3.0 timechange_0.4.0
#> [41] tidyselect_1.2.1 abind_1.4-8 codetools_0.2-20 spatstat.random_3.4-5 spatstat.explore_3.8-0
#> [46] lattice_0.22-9 tibble_3.3.1 shiny_1.13.0 withr_3.0.3 S7_0.2.2
#> [51] evaluate_1.0.5 sf_1.1-1 units_1.0-1 proxy_0.4-29 polyclip_1.10-7
#> [56] pillar_1.11.1 KernSmooth_2.23-26 plotly_4.12.0 dbscan_1.2.4 scales_1.4.0
#> [61] xtable_1.8-8 class_7.3-23 glue_1.8.1 janitor_2.2.1 pheatmap_1.0.13
#> [66] lazyeval_0.2.3 tools_4.6.0 hexDensity_1.4.10 hexbin_1.28.5 data.table_1.18.2.1
#> [71] grid_4.6.0 tidyr_1.3.2 nlme_3.1-169 fastmatrix_0.6-6 cli_3.6.6
#> [76] spatstat.sparse_3.1-0 S4Arrays_1.11.1 viridisLite_0.4.3 dplyr_1.2.1 gtable_0.3.6
#> [81] SpatialPack_0.4-1 digest_0.6.39 classInt_0.4-11 SparseArray_1.11.13 rjson_0.2.23
#> [86] htmlwidgets_1.6.4 farver_2.1.2 htmltools_0.5.9 lifecycle_1.0.5 httr_1.4.8
#> [91] mime_0.13