SpaceMarkers Functions Step-by-Step
2026-04-29
1 Overview
This vignette mirrors the SpaceMarkers with SpatialExperiment Objects vignette
but replaces the single SpaceMarkers() call with each constituent function
called explicitly. Every step dispatches on the SpaceMarkersExperiment (SME)
and returns an updated SME, so the analysis is a chain of SME-in / SME-out
calls — no manual slot assignment or pre-extracted matrices required.
The same breast cancer Visium dataset is used throughout.
2 Installation
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("SpaceMarkers")library(SpaceMarkers)3 Setup
3.1 Links to Data
data_dir <- "visiumBrCA"
unlink(file.path(data_dir), recursive = TRUE)
dir.create(data_dir, showWarnings = FALSE)
main_10xlink <- "https://cf.10xgenomics.com/samples/spatial-exp/1.3.0"
counts_folder <- "Visium_Human_Breast_Cancer"
counts_file <- "Visium_Human_Breast_Cancer_filtered_feature_bc_matrix.h5"
counts_url <- paste(c(main_10xlink, counts_folder, counts_file), collapse = "/")
sp_folder <- "Visium_Human_Breast_Cancer"
sp_file <- "Visium_Human_Breast_Cancer_spatial.tar.gz"
sp_url <- paste(c(main_10xlink, sp_folder, sp_file), collapse = "/")3.2 Downloading Data
Downloads are routed through BiocFileCache so the ~260 MB dataset is
fetched once and reused across vignette builds and across the other
SpaceMarkers vignettes.
bfc <- BiocFileCache::BiocFileCache(ask = FALSE)
counts_cache <- BiocFileCache::bfcrpath(bfc, counts_url)
sp_cache <- BiocFileCache::bfcrpath(bfc, sp_url)
file.copy(counts_cache, file.path(data_dir, counts_file), overwrite = TRUE)## [1] TRUEuntar(sp_cache, exdir = data_dir)4 Loading Data into a SpaceMarkersExperiment
load10X() searches for expression and spatial files with flexible regex
matching, log-transforms counts, and returns an SME with pre-computed kernel
parameters in spatial_params().
cogaps_result <- readRDS(system.file("extdata", "CoGAPS_result.rds",
package = "SpaceMarkers", mustWork = TRUE))
sme <- load10X(data_dir, features = cogaps_result, method = "CoGAPS")
sme## class: SpaceMarkersExperiment
## dim: 36601 4898
## metadata(0):
## assays(1): logcounts
## rownames(36601): MIR1302-2HG FAM138A ... AC007325.4 AC007325.2
## rowData names(0):
## colnames(4898): AAACAACGAATAGTTC-1 AAACAAGTATCTCCCA-1 ...
## TTGTTTGTATTACACG-1 TTGTTTGTGTAAATTC-1
## colData names(6): Pattern_1 Pattern_2 ... Pattern_5 sample_id
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):
## spatialCoords names(2) : y x
## imgData names(0):
## spacemarkers analysis: none
## patterns: Pattern_1, Pattern_2, Pattern_3, Pattern_4, Pattern_54.1 Filtering
Retain a curated gene list and spots expressed in at least 3 cells:
data("curated_genes")
good_gene_threshold <- 3
keep_genes <- rownames(sme)[
apply(assay(sme, "logcounts"), 1, function(x) sum(x > 0) >= good_gene_threshold)
]
keep_genes <- intersect(keep_genes, curated_genes)
sme <- sme[keep_genes, ]
sme## class: SpaceMarkersExperiment
## dim: 114 4898
## metadata(0):
## assays(1): logcounts
## rownames(114): C1orf127 CDC42 ... ATRX UBE2A
## rowData names(0):
## colnames(4898): AAACAACGAATAGTTC-1 AAACAAGTATCTCCCA-1 ...
## TTGTTTGTATTACACG-1 TTGTTTGTGTAAATTC-1
## colData names(6): Pattern_1 Pattern_2 ... Pattern_5 sample_id
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):
## spatialCoords names(2) : y x
## imgData names(0):
## spacemarkers analysis: none
## patterns: Pattern_1, Pattern_2, Pattern_3, Pattern_4, Pattern_55 Accessors
The SME holds everything the pipeline needs. A few accessors are worth knowing before walking through the steps:
| Accessor | Returns |
|---|---|
assay(sme, "logcounts") |
Genes × spots expression matrix |
spatialCoords(sme) |
Spot coordinates |
spatial_patterns(sme) |
Pattern values per spot |
spatial_params(sme) |
Pre-computed kernel parameters (2 × n) |
The SME pipeline methods read these internally — you do not need to extract them by hand before each step.
spatial_params(sme)## Pattern_1 Pattern_2 Pattern_3 Pattern_4 Pattern_5
## sigmaOpt 176.5483 176.5483 176.5483 176.5483 176.5483
## threshOpt 4.0000 4.0000 4.0000 4.0000 4.00006 Undirected Analysis
The undirected analysis identifies genes that are differentially expressed in the overlap zone between two patterns. We focus on Pattern_1 (immune) and Pattern_2, which show substantial spatial co-occurrence in this dataset.
Every step of the pipeline is an SME-in / SME-out method: the SME is passed
through a chain of step functions that each read their inputs from the object
and write their outputs back to the appropriate slot. Writes happen via the
setter generics (hotspots<-, interactions<-, undirected_scores<-,
overlap_scores<-, analysis_type<-); no manual @spacemarkers or
metadata() assignment is needed.
We’ll walk through the five pipeline steps one at a time, visualizing what each step produces.
6.1 Step 0: Visualize the starting patterns
Before running any pipeline step, inspect the CoGAPS patterns on the tissue.
plot_spatial() resolves a feature name against colData (for patterns),
rownames(sme) (for gene expression), or the hotspot metadata (via the
hotspot_type argument).
plot_spatial(sme, feature_col = "Pattern_1", title = "Pattern_1 (immune)")
And the second pattern we’ll use for the undirected analysis:
plot_spatial(sme, feature_col = "Pattern_2", title = "Pattern_2")
6.2 Step 1: Find hotspots
find_all_hotspots() labels each spot as a hotspot (or NA) for every
pattern, using a Gaussian-kernel null distribution with the pre-computed
spatial_params(sme). The result is stored in
metadata(sme)$hotspots$undirected via the hotspots<- setter.
sme_ud <- find_all_hotspots(sme)# Fraction of spots classified as hotspots per pattern
hs <- hotspots(sme_ud, "undirected")
pat_cols <- setdiff(colnames(hs), c("barcode", "y", "x"))
vapply(pat_cols, function(p) mean(!is.na(hs[[p]])), numeric(1))## Pattern_1 Pattern_2 Pattern_3 Pattern_4 Pattern_5
## 0.2592895 0.3938342 0.2719477 0.4013883 0.3084933Plot the hotspot labels for each pattern — pass source = "hotspots" so
plot_spatial() reads from the hotspots metadata rather than the continuous
pattern column in colData. The colors argument controls fill: for
categorical hotspot labels ("hot" / NA) plot_spatial() passes it into
ggplot2::scale_fill_manual(), so a single color string maps the "hot"
level to that color. Other aesthetic knobs: point_size, stroke, alpha,
bg_color, text_size, and crop.
plot_spatial(sme_ud, feature_col = "Pattern_1",
source = "hotspots", hotspot_type = "undirected",
colors = "steelblue", title = "Pattern_1 hotspots")
…and the Pattern_2 hotspots picked up by the same call (a single
find_all_hotspots() populates hotspots for every pattern column):
plot_spatial(sme_ud, feature_col = "Pattern_2",
source = "hotspots", hotspot_type = "undirected",
colors = "firebrick", title = "Pattern_2 hotspots")
A combined interaction overlay classifies each hotspot spot as Pattern_1
only, interacting (hot in both), or Pattern_2 only. Pass
source = "interaction" with interaction_patterns = c("Pattern_1", "Pattern_2") and a three-color vector to distinguish the three classes.
plot_spatial(sme_ud,
source = "interaction", hotspot_type = "undirected",
interaction_patterns = c("Pattern_1", "Pattern_2"),
colors = c(Pattern_1 = "steelblue",
interacting = "purple",
Pattern_2 = "firebrick"))
6.3 Step 2: Calculate overlap scores
Before scoring genes we already know where every pair of hotspots sits
in space, so we can compute the pairwise spatial overlap.
calculate_overlap_undirected() measures how much each pair’s hotspots
overlap (Szymkiewicz–Simpson by default); the resulting table is stored
via overlap_scores<- and analysis_type(sme) is set to
"undirected".
sme_ud <- calculate_overlap_undirected(sme_ud)
overlap_scores(sme_ud)## pattern1 pattern2 overlapScore
## 2 Pattern_1 Pattern_2 0.50944882
## 3 Pattern_1 Pattern_3 0.17795276
## 4 Pattern_1 Pattern_4 0.31811024
## 5 Pattern_1 Pattern_5 0.15984252
## 8 Pattern_2 Pattern_3 0.42792793
## 9 Pattern_2 Pattern_4 0.33903577
## 10 Pattern_2 Pattern_5 0.22964924
## 14 Pattern_3 Pattern_4 0.74774775
## 15 Pattern_3 Pattern_5 0.04654655
## 20 Pattern_4 Pattern_5 0.12045003analysis_type(sme_ud)## [1] "undirected"plot_overlap_scores(sme_ud)
The heatmap highlights which pattern pairs share territory; pairs with substantial overlap are the interesting candidates for the gene-level tests in the next step.
6.4 Step 3: Find pairwise interacting genes
get_pairwise_interacting_genes() tests each gene for differential
expression between a pair’s overlap zone and each pattern’s exclusive
zone. The pattern_pairs argument scopes the test to a single pair so
the narrative stays focused; omit it to compute all pairs. The full
per-gene result list is written to metadata(sme)$interactions.
sme_ud <- get_pairwise_interacting_genes(
sme_ud,
mode = "DE", analysis = "enrichment",
minOverlap = 10, workers = 1,
pattern_pairs = rbind(c("Pattern_1", "Pattern_2"))
)# Detailed per-gene statistics for the Pattern_1 x Pattern_2 pair
names(interactions(sme_ud))## [1] "Pattern_1_Pattern_2"pair1 <- interactions(sme_ud)[[1]]
if (length(pair1$interacting_genes) > 0) {
head(pair1$interacting_genes[[1]])
}## Gene Pattern_1 x Pattern_2 KW.obs.tot KW.obs.groups KW.df
## COMP COMP vsBoth 2552 3 2
## MAP2K2 MAP2K2 vsBoth 2552 3 2
## HLA-A HLA-A vsBoth 2552 3 2
## CST1 CST1 vsBoth 2552 3 2
## ARHGEF40 ARHGEF40 vsBoth 2552 3 2
## SLC12A9 SLC12A9 vsBoth 2552 3 2
## KW.statistic KW.pvalue KW.p.adj Dunn.zP1_Int Dunn.zP2_Int
## COMP 46.91225 6.503296e-11 2.471253e-10 -5.034959 -6.686948
## MAP2K2 34.41648 3.361679e-08 1.094947e-07 -5.136201 -5.259909
## HLA-A 181.00503 4.957447e-40 3.532181e-39 -7.453216 -13.450969
## CST1 48.70199 2.657646e-11 1.044730e-10 -3.785701 -6.974570
## ARHGEF40 13.93738 9.408859e-04 1.986315e-03 -2.813471 -3.621020
## SLC12A9 17.76551 1.387617e-04 3.101732e-04 -2.655819 -4.203384
## Dunn.zP2_P1 Dunn.pval_1_Int Dunn.pval_2_Int Dunn.pval_2_1
## COMP -0.8161070 4.779513e-07 2.278734e-11 4.144389e-01
## MAP2K2 0.7093684 2.803481e-07 1.441270e-07 4.780959e-01
## HLA-A -4.7157840 9.109187e-14 3.038637e-41 2.407818e-06
## CST1 -2.5359309 1.532757e-04 3.068081e-12 1.121489e-02
## ARHGEF40 -0.3419204 4.900978e-03 2.934437e-04 7.324108e-01
## SLC12A9 -1.0981652 7.911614e-03 2.629536e-05 2.721324e-01
## Dunn.pval_1_Int.adj Dunn.pval_2_Int.adj Dunn.pval_2_1.adj
## COMP 1.537756e-06 1.053915e-10 5.444109e-01
## MAP2K2 9.429891e-07 5.078761e-07 6.170773e-01
## HLA-A 5.185229e-13 3.212274e-40 7.322404e-06
## CST1 3.739254e-04 1.513586e-11 1.945082e-02
## ARHGEF40 9.628470e-03 6.715928e-04 9.083531e-01
## SLC12A9 1.439654e-02 6.949488e-05 3.683743e-01
## SpaceMarkersMetric
## COMP 5.115554
## MAP2K2 4.808176
## HLA-A 4.476319
## CST1 3.512456
## ARHGEF40 3.483833
## SLC12A9 3.4809796.5 Step 4: Calculate IM scores
get_im_scores() aggregates the per-gene Kruskal–Wallis + Dunn’s-test
statistics from the previous step into a single score matrix (genes × pattern
pairs), stored via undirected_scores<-.
sme_ud <- get_im_scores(sme_ud)
head(undirected_scores(sme_ud))## Gene Pattern_1_Pattern_2
## COMP COMP 5.115554
## MAP2K2 MAP2K2 4.808176
## HLA-A HLA-A 4.476319
## CST1 CST1 3.512456
## ARHGEF40 ARHGEF40 3.483833
## SLC12A9 SLC12A9 3.480979The top-ranked genes for the Pattern_1 × Pattern_2 interaction:
pair_col <- setdiff(colnames(undirected_scores(sme_ud)), "Gene")[1]
plot_im_scores(sme_ud, interaction = pair_col, nGenes = 10)
Plot the top gene’s spatial expression — the IM score is a summary of where a gene is differentially expressed in the overlap zone.
ims <- undirected_scores(sme_ud)
top_gene <- ims[order(ims[[pair_col]], decreasing = TRUE), "Gene"][1]
plot_spatial(sme_ud, feature_col = top_gene,
title = sprintf("%s (top %s gene)", top_gene, pair_col))
7 Directed Analysis
The directed analysis quantifies how much each pattern’s spatial influence drives gene expression in the hotspot of the opposing pattern. We focus on Pattern_1 (immune) → Pattern_5 (cancer), matching the main vignette.
The directed workflow uses Gaussian mixture model thresholds rather than the
kernel null distribution used in the undirected workflow, so the hotspot
step differs — you supply type = "pattern" (on the pattern values) or
type = "influence" (on the smoothed influence map).
7.1 Step 0: Visualize the starting patterns
plot_spatial(sme, feature_col = "Pattern_1", title = "Pattern_1 (immune)")
The directed analysis focuses on how Pattern_1 (immune) radiates into the Pattern_5 (cancer) region, so we’ll track the cancer pattern in the same way:
plot_spatial(sme, feature_col = "Pattern_5", title = "Pattern_5 (cancer)")
7.2 Step 1: Calculate influence
calculate_influence() applies kernel smoothing to each pattern column
independently, producing per-spot values that describe how far each pattern
“radiates” into the surrounding tissue. The result is a data.frame stored in
metadata(sme)$influence via influence_map<-.
sme_dir <- calculate_influence(sme)
head(influence_map(sme_dir))## barcode x y Pattern_1 Pattern_2
## AAACAACGAATAGTTC-1 AAACAACGAATAGTTC-1 14376 4662 0.34859369 0.10298807
## AAACAAGTATCTCCCA-1 AAACAAGTATCTCCCA-1 25987 16545 0.16456297 0.31619677
## AAACAATCTACTAGCA-1 AAACAATCTACTAGCA-1 18039 5392 0.06432877 0.11121815
## AAACACCAATAACTGC-1 AAACACCAATAACTGC-1 14703 18606 0.05964437 0.15881943
## AAACAGAGCGACTCCT-1 AAACAGAGCGACTCCT-1 24950 8032 0.34172043 0.26491974
## AAACAGCTTTCAGAAG-1 AAACAGCTTTCAGAAG-1 13367 14818 0.08296669 0.06691587
## Pattern_3 Pattern_4 Pattern_5
## AAACAACGAATAGTTC-1 0.17965177 0.62482126 0.022594385
## AAACAAGTATCTCCCA-1 0.23027045 0.35776379 0.039713306
## AAACAATCTACTAGCA-1 0.24245830 0.40456609 0.059994366
## AAACACCAATAACTGC-1 0.06351464 0.07489106 0.566755550
## AAACAGAGCGACTCCT-1 0.02251637 0.17334740 0.009465141
## AAACAGCTTTCAGAAG-1 0.07077331 0.27682185 0.522630852Plot the smoothed influence map for each pattern (source = "influence_map"
makes plot_spatial read from metadata(sme)$influence instead of the
continuous pattern column in colData):
plot_spatial(sme_dir, feature_col = "Pattern_1",
source = "influence_map",
title = "Pattern_1 influence")
Pattern_5’s influence map, on the same tissue:
plot_spatial(sme_dir, feature_col = "Pattern_5",
source = "influence_map",
title = "Pattern_5 influence")
7.3 Step 2: Find pattern hotspots
find_hotspots_gmm(type = "pattern") fits a two-component Gaussian mixture
model to each pattern’s value distribution and labels spots above the
boundary as hotspots. Stored under metadata(sme)$hotspots$pattern.
sme_dir <- find_hotspots_gmm(sme_dir, type = "pattern")## number of iterations= 83
## number of iterations= 300
## number of iterations= 70
## number of iterations= 132
## number of iterations= 65plot_spatial(sme_dir, feature_col = "Pattern_1",
source = "hotspots", hotspot_type = "pattern",
colors = "steelblue", title = "Pattern_1 GMM hotspots")
…and Pattern_5’s GMM hotspot, which defines the target region the directed analysis will test:
plot_spatial(sme_dir, feature_col = "Pattern_5",
source = "hotspots", hotspot_type = "pattern",
colors = "firebrick", title = "Pattern_5 GMM hotspots")
7.4 Step 3: Find influence hotspots
Running find_hotspots_gmm(type = "influence") applies the same GMM
thresholding but to the smoothed influence map from Step 1 — spots under
strong influence of each pattern. Stored under
metadata(sme)$hotspots$influence.
sme_dir <- find_hotspots_gmm(sme_dir, type = "influence")## number of iterations= 87
## number of iterations= 224
## number of iterations= 37
## number of iterations= 114
## number of iterations= 70plot_spatial(sme_dir, feature_col = "Pattern_1",
source = "hotspots", hotspot_type = "influence",
colors = "steelblue", title = "Pattern_1 influence hotspots")
Pattern_5’s influence hotspot — the complementary source region used in the reverse-direction test (Pattern_5 influencing Pattern_1):
plot_spatial(sme_dir, feature_col = "Pattern_5",
source = "hotspots", hotspot_type = "influence",
colors = "firebrick", title = "Pattern_5 influence hotspots")
With both pattern and influence hotspots in place, the directed interaction zone is well-defined. A directed analysis runs two independent t-tests per pair — one for each direction — so the per-spot classification is best read one direction at a time:
- Forward (Pattern_1 → Pattern_5): inside Pattern_1’s GMM
hotspot, split spots into
"Pattern_1"(not in Pattern_5’s influence) and"Pattern_1 near Pattern_5"(in Pattern_5’s influence — the t-test’s interaction zone). The “control” vs “treatment” groups the gene-level test compares. - Reverse (Pattern_5 → Pattern_1): the mirror —
"Pattern_5"vs"Pattern_5 near Pattern_1".
Each plot also carries a "<target> influence" context layer so the
reach of the target’s influence is visible beyond the source’s
hotspot.
labs_fwd <- overlap_map(sme_dir, c("Pattern_1", "Pattern_5"),
direction = "forward")
table(labs_fwd, useNA = "ifany")## labs_fwd
## Pattern_1 Pattern_1 near Pattern_5 Pattern_5 influence
## 922 246 1451
## <NA>
## 2279plot_spatial(sme_dir,
source = "interaction",
interaction_patterns = c("Pattern_1", "Pattern_5"),
direction = "forward",
colors = c(Pattern_1 = "steelblue",
"Pattern_1 near Pattern_5" = "purple",
"Pattern_5 influence" = "gray60"))
labs_rev <- overlap_map(sme_dir, c("Pattern_1", "Pattern_5"),
direction = "reverse")
table(labs_rev, useNA = "ifany")## labs_rev
## Pattern_5 Pattern_5 near Pattern_1 Pattern_1 influence
## 1586 246 767
## <NA>
## 2299plot_spatial(sme_dir,
source = "interaction",
interaction_patterns = c("Pattern_1", "Pattern_5"),
direction = "reverse",
colors = c(Pattern_5 = "firebrick",
"Pattern_5 near Pattern_1" = "goldenrod",
"Pattern_1 influence" = "gray60"))
7.5 Step 4: Directed overlap scores
Before scoring genes we look at the pairwise directional overlap across
every source → target pair. calculate_overlap_directed() quantifies
how much each pattern’s GMM hotspot overlaps with the opposing pattern’s
influence zone, stores the table via overlap_scores<-, and sets
analysis_type(sme) to "directed".
sme_dir <- calculate_overlap_directed(sme_dir)
overlap_scores(sme_dir)## pattern influence relAbundance
## 2 Pattern_2 near.Pattern_1 1.2466920
## 3 Pattern_3 near.Pattern_1 0.6873046
## 4 Pattern_4 near.Pattern_1 0.6675868
## 5 Pattern_5 near.Pattern_1 0.6492605
## 6 Pattern_1 near.Pattern_2 1.0483733
## 8 Pattern_3 near.Pattern_2 1.0719174
## 9 Pattern_4 near.Pattern_2 0.1683509
## 10 Pattern_5 near.Pattern_2 0.0000000
## 11 Pattern_1 near.Pattern_3 0.8482350
## 12 Pattern_2 near.Pattern_3 1.3481380
## 14 Pattern_4 near.Pattern_3 2.0392518
## 15 Pattern_5 near.Pattern_3 0.5023964
## 16 Pattern_1 near.Pattern_4 0.7357006
## 17 Pattern_2 near.Pattern_4 0.1710897
## 18 Pattern_3 near.Pattern_4 2.0268984
## 20 Pattern_5 near.Pattern_4 0.0837588
## 21 Pattern_1 near.Pattern_5 0.6078959
## 22 Pattern_2 near.Pattern_5 0.2715306
## 23 Pattern_3 near.Pattern_5 0.4156515
## 24 Pattern_4 near.Pattern_5 0.1083318analysis_type(sme_dir)## [1] "directed"plot_overlap_scores() auto-detects whether the stored scores are
undirected (columns pattern1, pattern2, overlapScore) or directed
(columns pattern, influence, relAbundance) and picks the right
axes and fill scale accordingly.
plot_overlap_scores(sme_dir, title = "Directed overlap (relative abundance)")
The heatmap flags which source → target pairs have the most directional interaction — good candidates for the per-gene tests in the next step.
7.6 Step 5: Directed gene scores
calculate_gene_scores_directed() scores each gene for each source → target
pair, using the intersection of the target’s GMM hotspot and the source’s
influence hotspot as the interaction zone. pattern_pairs restricts to the
Pattern_1 × Pattern_5 pair. Results written via directed_scores<-.
sme_dir <- calculate_gene_scores_directed(
sme_dir,
pattern_pairs = rbind(c("Pattern_1", "Pattern_5"))
)head(directed_scores(sme_dir))## statistic p.value n1 n2 effect_size gene
## C1orf127 1.026674 3.053435e-01 246 922 0.08665055 C1orf127
## CDC42 5.566032 4.968939e-08 246 922 0.40810687 CDC42
## C1QC 6.081864 2.580757e-09 246 922 0.39911510 C1QC
## C1QB 2.676324 7.733878e-03 246 922 0.18059121 C1QB
## ZNF593 7.414034 8.697873e-13 246 922 0.55808811 ZNF593
## YTHDF2 3.475487 5.741073e-04 246 922 0.27064062 YTHDF2
## cell_interaction p.adj
## C1orf127 Pattern_1_near_Pattern_5 3.588573e-01
## CDC42 Pattern_1_near_Pattern_5 1.452459e-07
## C1QC Pattern_1_near_Pattern_5 8.287501e-09
## C1QB Pattern_1_near_Pattern_5 1.366918e-02
## ZNF593 Pattern_1_near_Pattern_5 3.361212e-12
## YTHDF2 Pattern_1_near_Pattern_5 1.223331e-03The top gene for the strongest directed interaction, plotted on the tissue:
ds <- directed_scores(sme_dir)
ci <- sort(unique(ds$cell_interaction))[1]
top_gene <- ds[ds$cell_interaction == ci, ]
top_gene <- top_gene[order(-top_gene$effect_size), "gene"][1]
plot_spatial(sme_dir, feature_col = top_gene,
title = sprintf("%s (top directed gene, %s)", top_gene, ci))
8 Putting it together
Once the individual steps make sense, the whole analysis is a single pipe chain:
sme_ud <- sme |>
find_all_hotspots() |>
calculate_overlap_undirected() |>
get_pairwise_interacting_genes(
mode = "DE", analysis = "enrichment",
minOverlap = 10, workers = 1,
pattern_pairs = rbind(c("Pattern_1", "Pattern_2"))
) |>
get_im_scores()
sme_dir <- sme |>
calculate_influence() |>
find_hotspots_gmm(type = "pattern") |>
find_hotspots_gmm(type = "influence") |>
calculate_overlap_directed() |>
calculate_gene_scores_directed(
pattern_pairs = rbind(c("Pattern_1", "Pattern_5"))
)9 Cleanup
unlink(file.path(data_dir), recursive = TRUE)10 Session Info
sessionInfo()## R version 4.6.0 RC (2026-04-17 r89917)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.4 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] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] SpaceMarkers_2.3.1 BiocStyle_2.41.0
##
## loaded via a namespace (and not attached):
## [1] RColorBrewer_1.1-3 jsonlite_2.0.0
## [3] shape_1.4.6.1 magrittr_2.0.5
## [5] spatstat.utils_3.2-2 magick_2.9.1
## [7] farver_2.1.2 rmarkdown_2.31
## [9] GlobalOptions_0.1.4 vctrs_0.7.3
## [11] memoise_2.0.1 spatstat.explore_3.8-0
## [13] matrixTests_0.2.3.1 tinytex_0.59
## [15] rstatix_0.7.3 bmp_0.3.1
## [17] htmltools_0.5.9 S4Arrays_1.13.0
## [19] readbitmap_0.1.5 curl_7.1.0
## [21] broom_1.0.12 SparseArray_1.13.0
## [23] Formula_1.2-5 sass_0.4.10
## [25] bslib_0.10.0 htmlwidgets_1.6.4
## [27] plyr_1.8.9 httr2_1.2.2
## [29] plotly_4.12.0 cachem_1.1.0
## [31] lifecycle_1.0.5 pkgconfig_2.0.3
## [33] Matrix_1.7-5 R6_2.6.1
## [35] fastmap_1.2.0 MatrixGenerics_1.25.0
## [37] digest_0.6.39 colorspace_2.1-2
## [39] S4Vectors_0.51.1 tensor_1.5.1
## [41] GenomicRanges_1.65.0 RSQLite_2.4.6
## [43] labeling_0.4.3 filelock_1.0.3
## [45] spatstat.sparse_3.1-0 httr_1.4.8
## [47] polyclip_1.10-7 abind_1.4-8
## [49] compiler_4.6.0 bit64_4.8.0
## [51] withr_3.0.2 S7_0.2.2
## [53] backports_1.5.1 tiff_0.1-12
## [55] BiocParallel_1.47.0 carData_3.0-6
## [57] viridis_0.6.5 DBI_1.3.0
## [59] MASS_7.3-65 rappdirs_0.3.4
## [61] DelayedArray_0.39.0 rjson_0.2.23
## [63] tools_4.6.0 otel_0.2.0
## [65] ape_5.8-1 goftest_1.2-3
## [67] glue_1.8.1 nanoparquet_0.5.1
## [69] nlme_3.1-169 grid_4.6.0
## [71] reshape2_1.4.5 generics_0.1.4
## [73] hdf5r_1.3.12 gtable_0.3.6
## [75] spatstat.data_3.1-9 tidyr_1.3.2
## [77] data.table_1.18.2.1 car_3.1-5
## [79] XVector_0.53.0 BiocGenerics_0.59.0
## [81] spatstat.geom_3.7-3 pillar_1.11.1
## [83] stringr_1.6.0 circlize_0.4.18
## [85] splines_4.6.0 dplyr_1.2.1
## [87] BiocFileCache_3.3.0 lattice_0.22-9
## [89] survival_3.8-6 bit_4.6.0
## [91] deldir_2.0-4 tidyselect_1.2.1
## [93] SingleCellExperiment_1.35.0 knitr_1.51
## [95] gridExtra_2.3 bookdown_0.46
## [97] IRanges_2.47.0 Seqinfo_1.3.0
## [99] SummarizedExperiment_1.43.0 stats4_4.6.0
## [101] xfun_0.57 Biobase_2.73.1
## [103] mixtools_2.0.0.1 matrixStats_1.5.0
## [105] stringi_1.8.7 lazyeval_0.2.3
## [107] yaml_2.3.12 codetools_0.2-20
## [109] evaluate_1.0.5 kernlab_0.9-33
## [111] effsize_0.8.1 tibble_3.3.1
## [113] qvalue_2.45.0 BiocManager_1.30.27
## [115] cli_3.6.6 segmented_2.2-1
## [117] jquerylib_0.1.4 dichromat_2.0-0.1
## [119] Rcpp_1.1.1-1.1 spatstat.random_3.4-5
## [121] dbplyr_2.5.2 png_0.1-9
## [123] spatstat.univar_3.1-7 parallel_4.6.0
## [125] ggplot2_4.0.3 blob_1.3.0
## [127] jpeg_0.1-11 SpatialExperiment_1.23.0
## [129] viridisLite_0.4.3 scales_1.4.0
## [131] purrr_1.2.2 rlang_1.2.0