1 Overview

This vignette demonstrates how to use the SpaceMarkersExperiment class, a SpatialExperiment subclass that serves as a unified container for all SpaceMarkers inputs and outputs. Instead of tracking separate matrices, data frames, and lists across pipeline steps, a single SpaceMarkersExperiment object holds everything: expression data, spatial coordinates, pattern features, optimal parameters, hotspots, interacting genes, and interaction scores.

We will reproduce the same breast cancer analysis from the main vignette using this integrated workflow.

2 Installation

if (!require("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("SpaceMarkers")
library(SpaceMarkers)

3 Setup

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] TRUE
untar(sp_cache, exdir = data_dir)

4 Loading Data into a SpaceMarkersExperiment

4.1 Using load10X

load10X() searches the Visium directory for expression and spatial files using flexible pattern matching — it works regardless of whether your files are in an outs/, data/, or flat directory layout. It log-transforms the counts and returns a SpaceMarkersExperiment directly:

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_5

Alternatively, you can load first and add features separately with add_features():

# Step 1: Load without features
sme <- load10X(data_dir)

# Step 2: Add features later
sme <- add_features(sme, cogaps_result, method = "CoGAPS")

If the features cover fewer spots than the expression data (or vice versa), add_features() will keep only the shared spots and emit a warning.

The object inherits from SpatialExperiment so all standard accessors work:

# Expression matrix
dim(assay(sme, "logcounts"))
## [1] 36601  4898
# Spatial coordinates
head(spatialCoords(sme))
##                        y     x
## AAACAACGAATAGTTC-1  4662 14376
## AAACAAGTATCTCCCA-1 16545 25987
## AAACAATCTACTAGCA-1  5392 18039
## AAACACCAATAACTGC-1 18606 14703
## AAACAGAGCGACTCCT-1  8032 24950
## AAACAGCTTTCAGAAG-1 14818 13367
# Spatial patterns stored in colData
head(spatial_patterns(sme))
## DataFrame with 6 rows and 5 columns
##                      Pattern_1   Pattern_2   Pattern_3   Pattern_4   Pattern_5
##                      <numeric>   <numeric>   <numeric>   <numeric>   <numeric>
## AAACAACGAATAGTTC-1 0.467625588 1.04939e-01 2.57606e-01 0.684806228 4.74709e-02
## AAACAAGTATCTCCCA-1 0.269075811 4.39442e-01 2.05647e-01 0.292133719 1.16758e-02
## AAACAATCTACTAGCA-1 0.110593386 1.14852e-02 2.30915e-01 0.411131471 9.50832e-02
## AAACACCAATAACTGC-1 0.000250838 1.68580e-01 1.22360e-01 0.000156279 8.04193e-01
## AAACAGAGCGACTCCT-1 0.284930885 1.10251e-01 9.05316e-08 0.242940620 3.43081e-08
## AAACAGCTTTCAGAAG-1 0.158373639 9.74108e-06 1.72347e-01 0.305995703 7.16760e-01

4.2 Filtering

For this demonstration we will use two patterns and a curated gene list:

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, ]

# Keep all patterns for now; we will subset per analysis below
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_5

5 Running SpaceMarkers with the SME Object

5.1 Undirected Analysis

The SpaceMarkers() function now accepts a SpaceMarkersExperiment as its first argument. All results are stored back in the returned object. We keep all five CoGAPS patterns on the SME (so the overlap heatmap shows every pairwise relationship) and restrict the gene-level interaction tests to Pattern_1 × Pattern_2 — our narrative focus — via pattern_pairs, which SpaceMarkers() forwards to get_pairwise_interacting_genes().

sme_undirected <- SpaceMarkers(
    x = sme,
    directed = FALSE,
    min.gene.expr = 3,
    pattern_pairs = rbind(c("Pattern_1", "Pattern_2"))
)

5.1.1 Inspecting Results

Summary scores are accessed via dedicated accessors:

# Analysis type
analysis_type(sme_undirected)
## [1] "undirected"
# IM scores (summary: genes x pattern pairs)
head(undirected_scores(sme_undirected))
##              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.480979

Detailed intermediate results are stored in metadata() and accessed via convenience functions:

# Undirected hotspots
hs <- hotspots(sme_undirected, type = "undirected")
head(hs)
##                               barcode     y     x Pattern_1 Pattern_2 Pattern_3
## AAACAACGAATAGTTC-1 AAACAACGAATAGTTC-1  4662 14376 Pattern_1      <NA> Pattern_3
## AAACAAGTATCTCCCA-1 AAACAAGTATCTCCCA-1 16545 25987 Pattern_1 Pattern_2 Pattern_3
## AAACAATCTACTAGCA-1 AAACAATCTACTAGCA-1  5392 18039      <NA>      <NA> Pattern_3
## AAACACCAATAACTGC-1 AAACACCAATAACTGC-1 18606 14703      <NA>      <NA>      <NA>
## AAACAGAGCGACTCCT-1 AAACAGAGCGACTCCT-1  8032 24950 Pattern_1      <NA>      <NA>
## AAACAGCTTTCAGAAG-1 AAACAGCTTTCAGAAG-1 14818 13367      <NA>      <NA>      <NA>
##                    Pattern_4 Pattern_5
## AAACAACGAATAGTTC-1 Pattern_4      <NA>
## AAACAAGTATCTCCCA-1 Pattern_4      <NA>
## AAACAATCTACTAGCA-1 Pattern_4      <NA>
## AAACACCAATAACTGC-1      <NA> Pattern_5
## AAACAGAGCGACTCCT-1      <NA>      <NA>
## AAACAGCTTTCAGAAG-1 Pattern_4 Pattern_5
# Per-pair interaction details
pair_names <- names(interactions(sme_undirected))
pair_names
## [1] "Pattern_1_Pattern_2"
# Detailed stats for the first pair (if any interactions were found)
if (length(pair_names) > 0) {
    pair1 <- interactions(sme_undirected, pair = pair_names[1])
    if (length(pair1$interacting_genes) > 0) {
        head(pair1$interacting_genes[[1]])
    } else {
        message("No interacting genes found for this pair.")
    }
}
##              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.480979

5.2 Directed Analysis

The directed analysis keeps every pattern on the SME so the directed overlap heatmap covers all source → target pairs, and restricts gene-level scoring to Pattern_1 → Pattern_5 (immune → cancer) via pattern_pairs, which SpaceMarkers() forwards to calculate_gene_scores_directed().

sme_directed <- SpaceMarkers(
    x = sme,
    directed = TRUE,
    min.gene.expr = 3,
    pattern_pairs = rbind(c("Pattern_1", "Pattern_5"))
)
## number of iterations= 82 
## number of iterations= 259 
## number of iterations= 64 
## number of iterations= 134 
## number of iterations= 60 
## number of iterations= 79 
## number of iterations= 180 
## number of iterations= 38 
## number of iterations= 121 
## number of iterations= 89
analysis_type(sme_directed)
## [1] "directed"
head(directed_scores(sme_directed))
##          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-03

Directed analysis also stores hotspot and influence intermediates:

# Pattern hotspots
head(hotspots(sme_directed, type = "pattern"))
##                               barcode     y     x Pattern_1 Pattern_2 Pattern_3
## AAACAACGAATAGTTC-1 AAACAACGAATAGTTC-1  4662 14376 Pattern_1      <NA> Pattern_3
## AAACAAGTATCTCCCA-1 AAACAAGTATCTCCCA-1 16545 25987 Pattern_1      <NA> Pattern_3
## AAACAATCTACTAGCA-1 AAACAATCTACTAGCA-1  5392 18039      <NA>      <NA> Pattern_3
## AAACACCAATAACTGC-1 AAACACCAATAACTGC-1 18606 14703      <NA>      <NA>      <NA>
## AAACAGAGCGACTCCT-1 AAACAGAGCGACTCCT-1  8032 24950 Pattern_1      <NA>      <NA>
## AAACAGCTTTCAGAAG-1 AAACAGCTTTCAGAAG-1 14818 13367      <NA>      <NA> Pattern_3
##                    Pattern_4 Pattern_5
## AAACAACGAATAGTTC-1 Pattern_4      <NA>
## AAACAAGTATCTCCCA-1      <NA>      <NA>
## AAACAATCTACTAGCA-1      <NA>      <NA>
## AAACACCAATAACTGC-1      <NA> Pattern_5
## AAACAGAGCGACTCCT-1      <NA>      <NA>
## AAACAGCTTTCAGAAG-1      <NA> Pattern_5
# Influence hotspots
head(hotspots(sme_directed, type = "influence"))
##                               barcode     y     x Pattern_1 Pattern_2 Pattern_3
## AAACAACGAATAGTTC-1 AAACAACGAATAGTTC-1  4662 14376 Pattern_1      <NA> Pattern_3
## AAACAAGTATCTCCCA-1 AAACAAGTATCTCCCA-1 16545 25987      <NA>      <NA> Pattern_3
## AAACAATCTACTAGCA-1 AAACAATCTACTAGCA-1  5392 18039      <NA>      <NA> Pattern_3
## AAACACCAATAACTGC-1 AAACACCAATAACTGC-1 18606 14703      <NA>      <NA>      <NA>
## AAACAGAGCGACTCCT-1 AAACAGAGCGACTCCT-1  8032 24950 Pattern_1      <NA>      <NA>
## AAACAGCTTTCAGAAG-1 AAACAGCTTTCAGAAG-1 14818 13367      <NA>      <NA>      <NA>
##                    Pattern_4 Pattern_5
## AAACAACGAATAGTTC-1 Pattern_4      <NA>
## AAACAAGTATCTCCCA-1      <NA>      <NA>
## AAACAATCTACTAGCA-1      <NA>      <NA>
## AAACACCAATAACTGC-1      <NA> Pattern_5
## AAACAGAGCGACTCCT-1      <NA>      <NA>
## AAACAGCTTTCAGAAG-1      <NA> Pattern_5
# Influence map
head(influence_map(sme_directed))
##                               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.522630852

6 Visualization

Since the SME holds all results, every visualization function accepts the object directly and reads from the stored slots.

6.1 Patterns on the tissue

plot_spatial() automatically extracts coordinates from the object, finds the tissue image stored in the object’s visiumDir, and falls back to a blank background if none is available. The source argument selects where feature_col is resolved: colData (patterns, continuous features), assay (gene expression), hotspots (hotspot labels), influence_map (per-pattern influence values), or interaction (a three-way categorical overlay).

plot_spatial(sme, feature_col = "Pattern_1", title = "Pattern_1 (immune)")

plot_spatial(sme, feature_col = "Pattern_5", title = "Pattern_5 (cancer)")

6.2 Undirected hotspots

Pass source = "hotspots" so plot_spatial() reads hotspot labels from metadata(sme)$hotspots$undirected rather than the continuous pattern column in colData (which shares the same name). colors maps the “hot” level to a single fill.

plot_spatial(sme_undirected, feature_col = "Pattern_1",
             source = "hotspots", hotspot_type = "undirected",
             colors = "steelblue", title = "Pattern_1 hotspots")

…and Pattern_2, from the same find_all_hotspots() pass:

plot_spatial(sme_undirected, feature_col = "Pattern_2",
             source = "hotspots", hotspot_type = "undirected",
             colors = "firebrick", title = "Pattern_2 hotspots")

6.3 Undirected overlap scores

With hotspots in hand we can already answer the broader spatial question — which pattern pairs actually share territory. plot_overlap_scores() reads overlap_scores(sme) directly and auto-detects the undirected shape (symmetric pattern1 × pattern2 heatmap) across every pair:

plot_overlap_scores(sme_undirected, title = "Undirected overlap")

6.4 Undirected interaction overlay

For a chosen pair, a three-way overlay classifies each hotspot spot as Pattern_1 only, interacting (hot in both), or Pattern_2 only. This is the spatial region the per-gene tests below operate on.

plot_spatial(sme_undirected,
             source = "interaction", hotspot_type = "undirected",
             interaction_patterns = c("Pattern_1", "Pattern_2"),
             colors = c(Pattern_1   = "steelblue",
                        interacting = "purple",
                        Pattern_2   = "firebrick"))

6.5 IM scores bar plot and top genes

The column names of undirected_scores(sme) are the pattern-pair keys — use the first pair key robustly instead of hard-coding:

pair_col <- setdiff(colnames(undirected_scores(sme_undirected)), "Gene")[1]
plot_im_scores(sme_undirected, interaction = pair_col, nGenes = 10)

plot_spatial() accepts gene names directly — it falls through to the expression assay via source = "auto" (the default):

ims <- undirected_scores(sme_undirected)
top_genes <- head(ims[order(ims[[pair_col]], decreasing = TRUE), "Gene"], 5)

for (g in top_genes) {
    print(plot_spatial(sme_undirected, feature_col = g, title = g))
}

6.6 Directed hotspots

The directed workflow splits hotspots into two flavors: pattern-value GMM hotspots (hotspot_type = "pattern") and influence-map GMM hotspots (hotspot_type = "influence"). Here we focus on the pattern-level hotspots that define the source/target regions of each directed pair.

plot_spatial(sme_directed, feature_col = "Pattern_1",
             source = "hotspots", hotspot_type = "pattern",
             colors = "steelblue", title = "Pattern_1 GMM hotspots")

…and Pattern_5, the target region the immune pattern radiates into:

plot_spatial(sme_directed, feature_col = "Pattern_5",
             source = "hotspots", hotspot_type = "pattern",
             colors = "firebrick", title = "Pattern_5 GMM hotspots")

6.7 Directed overlap scores

plot_overlap_scores() detects the directed shape (pattern × influence with a relative-abundance fill) automatically, again across every source → target pair:

plot_overlap_scores(sme_directed, title = "Directed overlap (relative abundance)")

6.8 Directed interaction overlays

Because a directed analysis runs two independent t-tests per pair — one per source → target direction — the per-spot classification is best rendered one direction at a time. The forward view shows Pattern_1’s hotspot split by whether each spot is also in Pattern_5’s influence zone (the test’s interaction group):

labs_fwd <- overlap_map(sme_directed, 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> 
##                     2279
plot_spatial(sme_directed,
             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"))

And the reverse (Pattern_5 → Pattern_1):

labs_rev <- overlap_map(sme_directed, 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> 
##                     2299
plot_spatial(sme_directed,
             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 Constructing a SpaceMarkersExperiment Manually

If you already have expression and coordinate data in R, you can construct a SpaceMarkersExperiment directly:

mat <- matrix(rnorm(100), nrow = 10, ncol = 10)
rownames(mat) <- paste0("gene_", 1:10)
colnames(mat) <- paste0("spot_", 1:10)

coords <- matrix(runif(20) * 100, ncol = 2)
colnames(coords) <- c("y", "x")
rownames(coords) <- colnames(mat)

cd <- S4Vectors::DataFrame(
    CellType_A = runif(10),
    CellType_B = runif(10),
    row.names = colnames(mat)
)

my_sme <- SpaceMarkersExperiment(
    assays = list(logcounts = mat),
    colData = cd,
    spatialCoords = coords,
    spaceMarkers = list(
        params = list(pattern_names = c("CellType_A", "CellType_B"))
    )
)
my_sme
## class: SpaceMarkersExperiment 
## dim: 10 10 
## metadata(0):
## assays(1): logcounts
## rownames(10): gene_1 gene_2 ... gene_9 gene_10
## rowData names(0):
## colnames(10): spot_1 spot_2 ... spot_9 spot_10
## colData names(3): CellType_A CellType_B sample_id
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):
## spatialCoords names(2) : y x
## imgData names(0):
## spacemarkers analysis: none 
##   patterns: CellType_A, CellType_B

7.1 From an Existing SpatialExperiment

You can also extract features from any SpatialExperiment object using get_spatial_features():

spe <- SpatialExperiment::SpatialExperiment(
    assays = list(logcounts = mat),
    colData = cd,
    spatialCoords = coords
)

# Extract numeric colData columns as spatial features
features <- get_spatial_features(spe)
head(features)
##        CellType_A  CellType_B
## spot_1  0.3905260 0.285896066
## spot_2  0.9292278 0.003628604
## spot_3  0.8214613 0.156325010
## spot_4  0.1798048 0.825986990
## spot_5  0.3074150 0.682848027
## spot_6  0.3597079 0.962238463

8 Backward Compatibility

The legacy file-path interface still works. Set returnSME = FALSE to get the old-style output:

# Legacy usage (returns data.frame, not SME)
result <- SpaceMarkers(
    features = "path/to/features.rds",
    data = "path/to/visium_dir",
    directed = FALSE,
    returnSME = FALSE
)

9 Summary

Feature Legacy Workflow SME Workflow
Input File paths SpaceMarkersExperiment or file paths
Output Scattered data.frames/lists Single SpaceMarkersExperiment
Expression Separate matrix assay(sme, "logcounts")
Coordinates Separate data.frame spatialCoords(sme)
Patterns Merged in spPatterns spatial_patterns(sme)
Kernel Parameters Separate matrix spatial_params(sme)
All Hyperparameters Scattered arguments params(sme)
Undirected Scores Separate data.frame undirected_scores(sme)
Directed Scores Separate data.frame directed_scores(sme)
LR Scores Separate matrix lr_scores(sme)
Overlap Scores Separate data.frame overlap_scores(sme)
Hotspots Separate data.frame hotspots(sme, type = ...)
Interactions Nested list interactions(sme, pair = ...)
Influence Separate data.frame influence_map(sme)
Spatial Plotting Manual coord extraction + plot_spatial_data_over_image plot_spatial(sme, feature_col = ...)

10 Cleanup

unlink(file.path(data_dir), recursive = TRUE)

11 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