MultiAssaySpatialExperiment use cases
Compiled: 2 September 2026; Modified: 6 July 2026
Package
MultiAssaySpatialExperiment 0.99.12
1 Introduction
This vignette demonstrates realistic workflows using MultiAssaySpatialExperiment
(MASE) for multi-assay spatial transcriptomics analysis. It is divided into two
parts:
Part 1: Data import: Reading spatial omics data from vendor outputs (Xenium, Visium, Visium HD, CosMx, MERSCOPE) and converting from other Bioconductor spatial formats.
Part 2: Real-world workflows: Complete analysis examples integrating multiple technologies, spatial transformations, and cross-platform label transfer.
Scientific question: How do cell-level and tissue-level spatial patterns complement each other, and can high-resolution Xenium data help interpret lower-resolution Visium spots?
library(MultiAssaySpatialExperiment)
library(SummarizedExperiment)
library(SingleCellExperiment)
library(S4Vectors)
library(sf)
# For visualization
library(ggplot2)
set.seed(123)2 Part 1: Data import
MASE provides readers for major spatial omics platforms and coercion methods for interoperability with other Bioconductor spatial classes.
2.1 Reading Xenium data
10x Genomics Xenium provides single-cell resolution subcellular imaging. The Xenium Explorer output directory contains:
cell_feature_matrix.h5: Cell × gene countscells.parquetorcells.csv.gz: Cell coordinates and metadatanucleus_boundaries.parquet: Cell segmentation polygons (optional)transcripts.parquet: Transcript-level coordinates (optional, large)
2.1.1 Basic usage
Two arguments matter more than the rest. segmentations chooses whether to load
the cell boundaries, the nucleus boundaries, or neither, and it is what decides
whether the object gets a shapes layer at all. add_transcripts controls the
transcript table, which is by far the largest part of a Xenium run and is worth
skipping until you need molecule-level detail.
mase_xenium <- readXeniumMASE(
data_dir = "path/to/xenium_output",
sample_id = "sample1",
segmentations = "cell",
add_transcripts = FALSE # Skip transcripts for faster loading
)The result holds the cell-by-gene matrix as a SingleCellExperiment in
experiments(mase)$rna, the cell centroids in spatialPoints(), and, when
segmentations is set, the cell or nucleus boundaries in spatialShapes().
spatialMap() records which spatial element each assay column sits at.
2.1.2 Options
| Parameter | Default | Purpose |
|---|---|---|
data_dir |
required | Path to Xenium output directory |
sample_id |
NULL |
Sample identifier (uses directory name if NULL) |
segmentations |
"cell" |
Segmentation to load: "cell", "nucleus", or "both" |
add_transcripts |
FALSE |
Load transcript-level coordinates (large!) |
images |
TRUE |
Load image metadata |
load_images |
FALSE |
Load image raster data |
2.1.3 When to load transcripts
Transcript-level data (millions of points) is useful for visualizing subcellular localization, quality control (transcript counts per cell), or custom aggregation strategies.
Skip transcripts (add_transcripts = FALSE) for faster loading and smaller
memory footprint when working with cell-level data only.
2.2 Reading Visium data
10x Genomics Visium provides spot-based spatial transcriptomics. Space Ranger output contains:
filtered_feature_bc_matrix.h5: Spot × gene countsspatial/tissue_positions.csv: Spot coordinatesspatial/scalefactors_json.json: Image scale factorsspatial/tissue_hires_image.png: H&E image (optional)
2.2.1 Basic usage
data selects the filtered matrix (spots under tissue) or the raw one (every
spot on the slide); filtered is what most analyses want. images = TRUE loads
the H&E image, which is what makes imgData() and the image-aware plotting
useful later, at the cost of reading the image file.
mase_visium <- readVisiumMASE(
data_dir = "path/to/spaceranger_output",
sample_id = "sample1",
images = TRUE,
data = "filtered"
)The result holds the spot-by-gene matrix in experiments() and the spot
coordinates in spatialPoints(). With images = TRUE the H&E image is
available through spatialImages(), and imgData() records which image belongs
to which specimen along with its scale factor.
2.2.2 Options
| Parameter | Default | Purpose |
|---|---|---|
data_dir |
required | Path to Space Ranger output directory |
sample_id |
NULL |
Sample identifier (uses directory name if NULL) |
type |
"HDF5" |
Matrix format: "HDF5" or "sparse" |
data |
"filtered" |
Matrix to read: "filtered" (spots under tissue) or "raw" |
images |
TRUE |
Load image metadata |
load_images |
FALSE |
Load image raster data |
unit |
"pixel" |
Coordinate units: "pixel" or "micron" |
min_area |
NULL |
Minimum polygon area for filtering |
block_size |
100 MB | Block size for chunked reading |
2.2.3 Visium HD support
Visium HD provides higher resolution with multiple binning levels (2 µm, 8 µm,
16 µm). Use the dedicated HD reader and specify bin size codes "002", "008",
or "016":
# Read Visium HD at 8 µm resolution
mase_visium_hd <- readVisiumHDMASE(
data_dir = "path/to/spaceranger_hd_output",
sample_id = "sample_hd",
bin_size = "008",
images = TRUE
)2.3 Reading CosMx data
NanoString CosMx provides subcellular resolution with compartment-specific segmentation (nucleus, membrane, cytoplasm). Output files include:
*_exprMat_file.csv: Expression matrix*_metadata_file.csv: Cell metadata*_tx_file.csv: Transcript coordinatesCellComposite/: Composite images
2.3.1 Basic usage
CosMx runs are organized into fields of view (FOVs), and a whole slide can be
more than you want in memory. fov_ids restricts the read to named FOVs; leaving
it NULL loads all of them.
mase_cosmx <- readCosMxMASE(
data_dir = "path/to/cosmx_output",
sample_id = "cosmx_sample",
fov_ids = NULL # NULL loads all FOVs; or c("1", "3", "5")
)As with the other readers, the counts land in experiments(), the cell
coordinates in spatialPoints(), and the segmentation in spatialShapes(),
which for CosMx includes subcellular compartments when the run provides them.
2.3.2 Options
| Parameter | Default | Purpose |
|---|---|---|
data_dir |
required | Path to CosMx output directory |
sample_id |
NULL |
Sample identifier (uses directory name if NULL) |
fov_ids |
NULL |
FOVs to load (NULL loads all) |
load_transcripts |
FALSE |
Load transcript-level coordinates |
images |
TRUE |
Load image metadata |
load_images |
FALSE |
Load image raster data |
min_area |
NULL |
Minimum polygon area for filtering |
2.3.3 FOV (Field of View) handling
CosMx data is organized by FOVs. Options:
fov_ids = NULL: Load all FOVs (combines into single MASE)fov_ids = c("1", "3", "5"): Load specific FOVsfov_ids = "1": Load a single FOV
2.4 Reading MERSCOPE data (MERFISH)
Vizgen MERSCOPE provides high-throughput spatial transcriptomics (often referred to as MERFISH). Output directory contains:
cell_by_gene.csv: Cell × gene matrixcell_metadata.csv: Cell coordinates and metadatadetected_transcripts.csv: Transcript coordinatescell_boundaries/: Cell segmentation polygons (Parquet)
2.4.1 Basic usage
segmentation names which segmentation run to load, since a MERSCOPE output
directory can hold more than one; "cellpose" is the usual default. As with
Xenium, load_transcripts gates the largest table in the run.
mase_merscope <- readMERSCOPEMASE(
data_dir = "path/to/vizgen_output",
sample_id = "merscope_sample",
segmentation = "cellpose",
load_transcripts = FALSE
)MERSCOPE runs are often three-dimensional, so spatialPoints() carries a z
column when the data provide one; the segmentation polygons are in
spatialShapes().
2.4.2 Options
| Parameter | Default | Purpose |
|---|---|---|
data_dir |
required | Path to MERSCOPE output directory |
sample_id |
NULL |
Sample identifier (uses directory name if NULL) |
fov_ids |
NULL |
FOVs to load (NULL loads all) |
segmentation |
"cellpose" |
Segmentation to load: "cellpose", "watershed", or "both" |
load_transcripts |
FALSE |
Load transcript-level coordinates |
images |
TRUE |
Load image metadata |
load_images |
FALSE |
Load image raster data |
min_area |
NULL |
Minimum polygon area for filtering |
min_qv |
20 | Minimum quality value for transcripts |
2.5 Choosing a reader
Select the reader that matches your platform output directory:
| Platform | Reader |
|---|---|
| 10x Xenium | readXeniumMASE() |
| 10x Visium | readVisiumMASE() |
| 10x Visium HD | readVisiumHDMASE() |
| NanoString CosMx | readCosMxMASE() |
| Vizgen MERSCOPE | readMERSCOPEMASE() |
2.6 Converting from other Bioconductor classes
If the data are already in a single-assay spatial class, there is no need to go back to the vendor output. Coercion is the usual way into MASE when you want to add a second assay to an analysis that started with just one.
2.6.1 From SpatialExperiment
library(SpatialExperiment)
# spe <- read10xVisium("path/to/visium")
mase_from_spe <- as(spe, "MultiAssaySpatialExperiment")The assays become the single element of the ExperimentList, colData carries
over as specimen metadata, and spatialCoords becomes a points layer. The
sampleMap and spatialMap are generated for you, so the result is immediately
a valid MASE with one assay: the point of the exercise is that you can now add a
second one alongside it.
2.6.2 From SpatialFeatureExperiment
library(SpatialFeatureExperiment)
# sfe <- read10xVisiumSFE("path/to/visium")
mase_from_sfe <- as(sfe, "MultiAssaySpatialExperiment")SpatialFeatureExperiment already carries geometry, so more survives the trip:
every assay is preserved, each colGeometry becomes a points or shapes layer
depending on whether it holds centroids or polygons, and the annotGeometries
are kept as further shapes layers. The generated spatialMap is what ties those
geometries back to the assay columns they describe.
2.7 Reading multiple samples
For multi-sample experiments, read each sample separately and combine:
mase1 <- readXeniumMASE(data_dir = "path/to/sample1", sample_id = "P001_xenium")
mase2 <- readXeniumMASE(data_dir = "path/to/sample2", sample_id = "P002_xenium")
mase_multi <- c(mase1, mase2)Each sample_id becomes a specimen in the combined object’s colData, and both
sets of spatial layers are carried across under their own names. c() also
rebuilds sampleMap and spatialMap over the union, which is what makes the
result subsettable by specimen afterwards: asking for one patient’s data returns
their assay columns and their geometry, and nothing from the other sample.
3 Part 2: Real-world workflows
Part 1 got data into the object. The workflows below are what you do once it is there, and each is chosen to exercise a different reason for using a multi-assay container: combining two platforms over one tissue, reconciling coordinate systems that do not agree, keeping many specimens straight, and the shape a typical analysis session takes end to end.
Unlike Part 1, these examples run on simulated data, so they execute when the vignette is built and you can modify them without vendor output to hand.
4 Workflow 1: Xenium + Visium integration
4.1 Background
A common experimental design combines Xenium, which gives single-cell resolution and subcellular localization but covers a limited gene panel (~300-500 genes), with Visium, which gives spot-level resolution (55 µm spots) across the whole transcriptome, at the cost of multiple cells per spot. Integrating both makes it possible to: 1. Use Xenium to identify cell types and their spatial organization 2. Use Visium to capture broader transcriptomic programs 3. Transfer cell type labels from Xenium to Visium spots 4. Deconvolve Visium spots using Xenium cell proportions
4.2 Data preparation
This example uses synthetic data that mimics a real experiment. In practice, you would load your own Xenium cell-by-gene matrix and cell centroids, and your own Visium spot-by-gene matrix and spot centroids.
# Xenium data: single-cell resolution
# Simulate 500 cells, 50 genes (subset of full panel)
n_cells <- 500
n_genes_xenium <- 50
xenium_counts <- matrix(
rpois(n_cells * n_genes_xenium, lambda = 5),
nrow = n_genes_xenium,
ncol = n_cells,
dimnames = list(
paste0("Gene", 1:n_genes_xenium),
paste0("Cell", 1:n_cells)
)
)
# Cell centroids (tissue coordinates in micrometers)
xenium_coords <- DataFrame(
x = runif(n_cells, 0, 1000),
y = runif(n_cells, 0, 1000),
instance_id = paste0("Cell", 1:n_cells)
)
# Cell type annotations (from clustering or known markers)
cell_types <- sample(
c("Neuron", "Astrocyte", "Microglia", "Oligodendrocyte"),
n_cells,
replace = TRUE,
prob = c(0.5, 0.2, 0.15, 0.15)
)
# Create SingleCellExperiment
xenium_sce <- SingleCellExperiment(
assays = list(counts = xenium_counts),
colData = DataFrame(
cell_type = cell_types,
tech = "xenium"
)
)# Visium data: spot-level resolution
# Simulate 100 spots, 1000 genes (whole transcriptome, downsampled)
n_spots <- 100
n_genes_visium <- 1000
# Create gene names: first 50 overlap with Xenium, rest are Visium-only
visium_gene_names <- c(
paste0("Gene", 1:n_genes_xenium), # Overlapping genes (50)
paste0("GeneV", 1:(n_genes_visium - n_genes_xenium)) # Visium-only (950)
)
visium_counts <- matrix(
rpois(n_spots * n_genes_visium, lambda = 10),
nrow = n_genes_visium,
ncol = n_spots,
dimnames = list(
visium_gene_names,
paste0("Spot", 1:n_spots)
)
)
# Spot centroids (same coordinate system as Xenium)
visium_coords <- DataFrame(
x = runif(n_spots, 0, 1000),
y = runif(n_spots, 0, 1000),
instance_id = paste0("Spot", 1:n_spots)
)
# Create SummarizedExperiment
visium_se <- SummarizedExperiment(
assays = list(counts = visium_counts),
colData = DataFrame(
tech = rep("visium", n_spots),
row.names = paste0("Spot", 1:n_spots)
)
)4.3 Build MASE object
Combine Xenium and Visium from the same tissue section into one MASE object.
Both assays share a single specimen (colData row); individual cells and spots
are observations linked through sampleMap and spatialMap.
# Specimen metadata: one row per biological unit (here, one tissue section)
specimen_id <- "section_01"
specimens <- DataFrame(
patient_id = "P001",
tissue = "mouse_brain",
section_id = specimen_id,
row.names = specimen_id
)
n_total <- n_cells + n_spots
# Sample map: link each assay column (cell or spot) to the shared specimen
sample_map <- DataFrame(
assay = factor(
c(rep("xenium", n_cells), rep("visium", n_spots)),
c("xenium", "visium")
),
primary = rep(specimen_id, n_total),
colname = c(paste0("Cell", 1:n_cells), paste0("Spot", 1:n_spots))
)
# Spatial map: link assay columns to spatial coordinates
spatial_map <- DataFrame(
assay = factor(
c(rep("xenium", n_cells), rep("visium", n_spots)),
c("xenium", "visium")
),
colname = c(paste0("Cell", 1:n_cells), paste0("Spot", 1:n_spots)),
element_type = rep("points", n_total),
region = factor(
c(rep("xenium_coords", n_cells), rep("visium_coords", n_spots)),
c("xenium_coords", "visium_coords")
),
instance_id = c(paste0("Cell", 1:n_cells), paste0("Spot", 1:n_spots))
)
# Construct MASE
mase <- MultiAssaySpatialExperiment(
experiments = ExperimentList(
xenium = xenium_sce,
visium = visium_se
),
colData = specimens,
sampleMap = sample_map,
points = PointsLayerList(
xenium_coords = xenium_coords,
visium_coords = visium_coords
),
spatialMap = spatial_map
)
mase
#> A MultiAssaySpatialExperiment object of 2 listed
#> experiments with user-defined names and respective classes.
#> Containing an ExperimentList class object of length 2:
#> [1] xenium: SingleCellExperiment with 50 rows and 500 columns
#> [2] visium: SummarizedExperiment with 1000 rows and 100 columns
#> Functionality:
#> experiments() - obtain the ExperimentList instance
#> colData() - the primary/phenotype DataFrame
#> sampleMap() - the sample coordination DataFrame
#> `$`, `[`, `[[` - extract colData columns, subset, or experiment
#> *Format() - convert into a long or wide DataFrame
#> assays() - convert ExperimentList to a SimpleList of matrices
#> exportClass() - save data to flat files
#> Spatial elements:
#> spatialImages: 0 elements
#> spatialLabels: 0 elements
#> spatialPoints: 2 elements
#> spatialShapes: 0 elements
#> imgData: NULL
#> spatialMap: present4.4 Exploratory analysis
With both assays in one object, the first thing worth checking is that they describe the same piece of tissue. Everything downstream, transferring labels, deconvolving spots, comparing expression, assumes the two coordinate systems agree, so it is cheap to confirm that before relying on it.
4.4.1 Do the two assays overlap in space?
Plotting the two side by side on the same axes shows both the resolution difference, thousands of cells against a coarse grid of spots, and whether they cover the same region:
# Extract coordinates
xenium_pts <- spatialPoint(mase, "xenium_coords")
visium_pts <- spatialPoint(mase, "visium_coords")
# Plot both assays
par(mfrow = c(1, 2))
# Xenium cells colored by type
# widen the x range so the legend sits beside the points rather than on top of them
plot(xenium_pts$x, xenium_pts$y,
col = as.factor(cell_types),
pch = 16, cex = 0.5,
main = "Xenium: Single cells",
xlab = "x (µm)", ylab = "y (µm)",
xlim = c(min(xenium_pts$x), max(xenium_pts$x) * 1.35))
legend("topright", legend = levels(as.factor(cell_types)),
col = seq_along(levels(as.factor(cell_types))), pch = 16,
cex = 0.7, bg = "white", box.col = "grey70")
# Visium spots
plot(visium_pts$x, visium_pts$y,
pch = 21, bg = "lightblue", cex = 2,
main = "Visium: Spots (55 µm)",
xlab = "x (µm)", ylab = "y (µm)")
4.5 Spatial annotation: Cell types in Visium spots
Use Xenium cell types to annotate Visium spots. For each spot, find which cells
fall within it and calculate cell type proportions. The loop below is explicit
for teaching; when spot polygons are stored in spatialShapes(), you can also
use annotateWithRegions() (see Working with MultiAssaySpatialExperiment).
# Create circular regions for Visium spots (55 µm diameter = 27.5 µm radius)
visium_circles <- st_buffer(
st_sf(visium_pts,
geometry = st_as_sfc(paste0("POINT(", visium_pts$x, " ", visium_pts$y, ")"))),
dist = 27.5
)
visium_shapes <- DataFrame(
geometry = st_geometry(visium_circles),
instance_id = visium_pts$instance_id
)
# Add shapes to MASE
spatialShapes(mase) <- ShapesLayerList(visium_spots = visium_shapes)
# For each Visium spot, find which Xenium cells are inside
xenium_sf <- st_sf(xenium_pts,
geometry = st_as_sfc(paste0("POINT(", xenium_pts$x, " ", xenium_pts$y, ")")))
visium_sf <- st_sf(
instance_id = visium_shapes$instance_id,
geometry = visium_shapes$geometry
)
# Spatial intersection
intersections <- st_intersects(visium_sf, xenium_sf)
# Calculate cell type proportions per spot
spot_cell_types <- lapply(seq_along(intersections), function(i) {
cell_indices <- intersections[[i]]
if (length(cell_indices) == 0) {
return(data.frame(
spot = visium_pts$instance_id[i],
Neuron = 0, Astrocyte = 0, Microglia = 0, Oligodendrocyte = 0,
n_cells = 0
))
}
types <- cell_types[cell_indices]
props <- base::table(types) / length(types)
data.frame(
spot = visium_pts$instance_id[i],
Neuron = as.numeric(props["Neuron"]),
Astrocyte = as.numeric(props["Astrocyte"]),
Microglia = as.numeric(props["Microglia"]),
Oligodendrocyte = as.numeric(props["Oligodendrocyte"]),
n_cells = length(cell_indices),
stringsAsFactors = FALSE
)
})
spot_annotations <- do.call(rbind, spot_cell_types)
spot_annotations[is.na(spot_annotations)] <- 0
head(spot_annotations)
#> spot Neuron Astrocyte Microglia Oligodendrocyte n_cells
#> 1 Spot1 1.0 0 0.0 0 1
#> 2 Spot2 0.0 0 0.0 0 0
#> 3 Spot3 0.0 1 0.0 0 1
#> 4 Spot4 0.0 0 0.0 0 0
#> 5 Spot5 0.0 0 0.0 0 0
#> 6 Spot6 0.5 0 0.5 0 24.5.1 Reading the composition result
Each spot now has a cell-type composition rather than a single label. The plots below show, per spot, the proportion of each Xenium cell type it contains. This is the quantity a deconvolution method would otherwise have to estimate; here it is counted directly, because the two assays share a coordinate system and the cells inside each spot are known.
# Plot cell type proportions per spot
par(mfrow = c(2, 2))
for (cell_type in c("Neuron", "Astrocyte", "Microglia", "Oligodendrocyte")) {
plot(visium_pts$x, visium_pts$y,
pch = 21, cex = 3,
bg = rgb(spot_annotations[[cell_type]], 0, 0, alpha = 0.7),
main = paste(cell_type, "proportion"),
xlab = "x (µm)", ylab = "y (µm)")
}
4.6 Spatial aggregation: Summarize Xenium by regions
Aggregate Xenium gene expression by tissue regions (e.g., cortical layers).
First, define regions, then use aggregateByRegion().
# Define two regions: "cortex" (y > 500) and "subcortex" (y <= 500)
# In practice, these would be anatomical annotations
region_polygons <- list(
cortex = st_polygon(list(cbind(
c(0, 1000, 1000, 0, 0),
c(500, 500, 1000, 1000, 500)
))),
subcortex = st_polygon(list(cbind(
c(0, 1000, 1000, 0, 0),
c(0, 0, 500, 500, 0)
)))
)
regions_sf <- st_sf(
region_id = c("cortex", "subcortex"),
geometry = st_sfc(region_polygons)
)
# Convert to DataFrame for MASE
anatomical_regions_df <- DataFrame(
geometry = st_geometry(regions_sf),
instance_id = c("cortex", "subcortex")
)
# Add to MASE, as a shapes layer named "anatomical_regions" (not "regions": that
# name is reserved for spatialMap's own `region` column, which names the points
# layer each row's instance_id comes from, "xenium_coords" here)
spatialShapes(mase) <- ShapesLayerList(
visium_spots = visium_shapes,
anatomical_regions = anatomical_regions_df
)# Annotate Xenium centroids with anatomical regions, then aggregate expression.
# Run on the full MASE so the anatomical_regions layer is retained
# (subsetByAssay drops shapes not referenced in spatialMap).
mase <- annotateWithRegions(mase,
points = "xenium_coords", shapes = "anatomical_regions")
xenium_by_region <- aggregateByRegion(mase, by = "anatomical_regions",
assays = "xenium", FUN = "sum")
xenium_by_region[["xenium"]][1:5, ]
#> cortex subcortex
#> Gene1 1256 1285
#> Gene2 1215 1220
#> Gene3 1311 1248
#> Gene4 1194 1276
#> Gene5 1171 12864.7 Comparing expression between the two platforms
With the shared genes identified, the question is whether the two platforms agree on the genes they both measure.
# Extract counts for shared genes
xenium_shared <- assay(experiments(mase)$xenium, "counts")[shared_genes, ]
visium_shared <- assay(experiments(mase)$visium, "counts")[shared_genes, ]
# Calculate mean expression per gene
xenium_means <- rowMeans(xenium_shared)
visium_means <- rowMeans(visium_shared)
# Plot correlation
plot(log1p(xenium_means), log1p(visium_means),
pch = 16, col = rgb(0, 0, 0, 0.5),
xlab = "log1p(mean Xenium counts)",
ylab = "log1p(mean Visium counts)",
main = "Gene expression correlation")
abline(0, 1, col = "red", lty = 2)
# Correlation
cor_val <- cor(xenium_means, visium_means, method = "spearman")
text(0.5, max(log1p(visium_means)) * 0.9,
paste("Spearman rho =", round(cor_val, 3)),
pos = 4)
4.8 Integration insights
Key findings from this analysis:
Spatial resolution: Xenium provides single-cell resolution, while Visium captures tissue-level patterns with multiple cells per spot.
Gene overlap: Shared genes allow direct comparison and validation between technologies.
Cell type transfer: Xenium cell type information can be transferred to Visium spots, revealing the cellular composition that a deconvolution method would otherwise have to estimate.
Regional patterns: Both assays can be aggregated by anatomical regions for comparative analysis.
Complementary data: Xenium’s high resolution complements Visium’s broad transcriptomic coverage.
5 Workflow 2: Coordinate transformation and alignment
When integrating data from different spatial platforms, coordinate systems often differ. This workflow demonstrates how to align Xenium and Visium data using spatial transformations.
5.1 Background
Xenium and Visium may use different coordinate systems: Xenium reports microns (µm) from a reference origin, while Visium reports grid indices (row/column) or pixel coordinates. To overlay them, you need to:
- Convert Visium grid indices to physical coordinates (microns)
- Apply transformations (translation, scaling, rotation) to align modalities
- Store transformed coordinates in MASE
5.2 Converting Visium grid to microns
Visium spots are arranged in a hexagonal grid. Each spot is ~55 µm diameter, with center-to-center spacing of ~100 µm.
# Visium spot coordinates (grid indices from Space Ranger)
visium_grid <- DataFrame(
row_idx = c(0, 0, 1, 1, 2, 2),
col_idx = c(0, 1, 0, 1, 0, 1),
instance_id = paste0("Spot", 1:6)
)
# Convert to microns (approximate)
# Spot spacing: ~100 µm horizontally, ~87 µm vertically (hexagonal)
visium_microns <- DataFrame(
x = visium_grid$col_idx * 100, # Column → x
y = visium_grid$row_idx * 87, # Row → y (hexagonal offset)
instance_id = visium_grid$instance_id
)
visium_microns
#> DataFrame with 6 rows and 3 columns
#> x y instance_id
#> <numeric> <numeric> <character>
#> 1 0 0 Spot1
#> 2 100 0 Spot2
#> 3 0 87 Spot3
#> 4 100 87 Spot4
#> 5 0 174 Spot5
#> 6 100 174 Spot65.3 Applying affine transformations
MASE allows storing the original and transformed coordinates side by side as named
layers of one object that stay attached to the assays through spatialMap; that
step comes right after this one, in Storing transformed coordinates. Computing
the transform itself is plain R: an affine transform is a matrix multiply, and
MASE does not impose a particular package or framework for computing it.
An affine transform composes translation, rotation and scaling into a single matrix:
# Define affine matrix: translate by (50, 30), scale by 1.2, no rotation
# Format: 2x3 matrix [a b c; d e f] where:
# x' = a*x + b*y + c
# y' = d*x + e*y + f
affine_matrix <- matrix(c(
1.2, 0.0, 50, # x' = 1.2*x + 0*y + 50 (scale + translate)
0.0, 1.2, 30 # y' = 0*x + 1.2*y + 30 (scale + translate)
), nrow = 2, byrow = TRUE)
# Apply transformation
apply_affine <- function(pts, mat) {
x_new <- mat[1, 1] * pts$x + mat[1, 2] * pts$y + mat[1, 3]
y_new <- mat[2, 1] * pts$x + mat[2, 2] * pts$y + mat[2, 3]
DataFrame(x = x_new, y = y_new, instance_id = pts$instance_id)
}
# Simulate Xenium coordinates
xenium_original <- DataFrame(
x = c(100, 150, 200, 250, 300),
y = c(100, 150, 200, 150, 100),
instance_id = paste0("Cell", 1:5)
)
# Transform Xenium to align with Visium
xenium_aligned <- apply_affine(xenium_original, affine_matrix)
xenium_aligned
#> DataFrame with 5 rows and 3 columns
#> x y instance_id
#> <numeric> <numeric> <character>
#> 1 170 150 Cell1
#> 2 230 210 Cell2
#> 3 290 270 Cell3
#> 4 350 210 Cell4
#> 5 410 150 Cell55.4 Landmark-based alignment
When you have known landmarks (tissue features visible in both modalities), use them to compute an optimal transformation:
# Example: three landmarks in both Xenium and Visium
# In practice, these would be manually identified features (blood vessels, etc.)
xenium_landmarks <- matrix(c(
100, 100,
200, 150,
300, 200
), ncol = 2, byrow = TRUE)
visium_landmarks <- matrix(c(
50, 50,
150, 100,
250, 150
), ncol = 2, byrow = TRUE)The workflow around those landmarks is the same whichever package computes the
transformation. Identify at least three corresponding points in both modalities, then
fit a transformation to them: Procrustes analysis for a rigid fit, least squares for a
full affine one, or a dedicated registration package such as RNiftyReg, imager or
EBImage. Apply the result to every point, check it visually, and store the
transformed coordinates back into the object as their own layer.
5.5 Storing transformed coordinates
Store both original and aligned coordinates as separate layers for flexibility:
# Create MASE with both coordinate systems
# Original Xenium coordinates
pts_original <- PointsLayerList(xenium_original = xenium_original)
# Aligned Xenium coordinates
pts_aligned <- PointsLayerList(xenium_aligned = xenium_aligned)
# Visium coordinates
pts_visium <- PointsLayerList(visium_coords = visium_microns)Keeping the original and aligned coordinates as separate layers of one object is deliberate: the two can be compared directly, the transformation can be refined and re-applied without re-reading anything, and downstream code chooses a coordinate system by naming a layer rather than by tracking which object holds which version.
5.6 Checking the alignment visually
Plotting the two modalities before and after the transformation is the quickest way to tell whether the alignment is sound.
par(mfrow = c(1, 2))
# Before alignment
plot(xenium_original$x, xenium_original$y,
pch = 16, col = "red", cex = 1.5,
main = "Before alignment",
xlab = "x (µm)", ylab = "y (µm)",
xlim = c(0, 400), ylim = c(0, 350))
points(visium_microns$x, visium_microns$y,
pch = 21, bg = "blue", cex = 3)
legend("topright", legend = c("Xenium (cells)", "Visium (spots)"),
col = c("red", "blue"), pch = c(16, 21), pt.cex = c(1.5, 2))
# After alignment
plot(xenium_aligned$x, xenium_aligned$y,
pch = 16, col = "red", cex = 1.5,
main = "After alignment",
xlab = "x (µm)", ylab = "y (µm)",
xlim = c(0, 400), ylim = c(0, 350))
points(visium_microns$x, visium_microns$y,
pch = 21, bg = "blue", cex = 3)
legend("topright", legend = c("Xenium (aligned)", "Visium (spots)"),
col = c("red", "blue"), pch = c(16, 21), pt.cex = c(1.5, 2))
6 Workflow 3: Multi-sample analysis
Once several patients or tissue sections are involved, the thing that usually goes wrong is bookkeeping: which columns of which assay belong to which patient, and which geometry belongs with them. Handling that by hand across a list of objects is where orphaned coordinates and mismatched subsets come from.
In MASE the answer is that specimens are first-class. Each patient is a row of
colData, and every assay column and every spatial element is tied back to one
of those rows through sampleMap and spatialMap.
specimens_multi <- DataFrame(
patient_id = c(rep("P001", 150), rep("P002", 120), rep("P003", 180)),
tissue = "mouse_brain",
condition = c(rep("control", 150), rep("treated", 120), rep("treated", 180)),
row.names = paste0("Sample", 1:450)
)
# sampleMap is built exactly as in the single-sample example above, with one
# row per assay column across all three patients.Because those maps exist, the questions you actually want to ask are ordinary subsetting operations, and each returns a complete, self-consistent object rather than a matrix that has lost track of its coordinates:
mase[, colData(mase)$patient_id == "P001", ] # one patient, every assay
mase[, colData(mase)$condition == "control", ] # one arm of the study
subsetByAssay(mase, "xenium") # one modality, every patientComparing conditions then means running the same per-patient analysis over those subsets, with the confidence that each carries its own geometry.
7 Workflow 4: Iterative analysis pattern
The three workflows above each made a specific point. This one is the shape most sessions actually take, and the reason to show it is the ordering: every step narrows the object, and doing them in this order means each subsequent step runs on less data while still being correct.
Start by reading the vendor output, then drop the specimens that failed quality control. Doing QC first is a habit worth keeping, because everything after it is more expensive per observation:
mase <- readXeniumMASE(data_dir = "path/to/xenium_output", sample_id = "sample1")
mase_qc <- mase[, colData(mase)$qc_pass, ]Restricting to a region of interest is the next cut. It is a spatial operation rather than an assay one, but because the links are maintained it narrows the assays too, so the annotation step that follows sees only the cells you care about:
library(sf)
roi <- st_polygon(list(matrix(
c(100, 200, 500, 200, 500, 600, 100, 600, 100, 200),
ncol = 2, byrow = TRUE)))
mase_roi <- subsetByPolygon(mase_qc, roi)Then annotate and aggregate. These two belong together and in this order:
aggregateByRegion() summarises over an annotation, so it needs
annotateWithRegions() to have recorded one on spatialMap first. Running them
on the narrowed object rather than the original is the payoff for the two cuts
above:
mase_annotated <- annotateWithRegions(mase_roi,
points = "centroids", shapes = "cells")
cell_expr <- aggregateByRegion(mase_annotated, by = "cells", FUN = "sum")cell_expr is an ordinary list of matrices, so differential expression,
clustering and the rest proceed with whatever tools you normally use. The object
is the thing that kept the coordinates, the specimens and the counts aligned up
to that point; it does not need to be involved afterwards.
8 Summary
Across these workflows the same property did the work each time. Two platforms
could be compared because both were mapped onto one coordinate system; a region
of interest could be cut without losing track of which cells it contained; three
patients could be held in one object and separated again on demand. In each case
what made it straightforward was that the links between assay columns, specimens
and spatial elements are recorded once, in sampleMap and spatialMap, and are
maintained by the operations rather than by the analyst.
That is also the limit of what the container claims. It does not implement
deconvolution, registration or differential expression; Workflow 2 computed its
coordinate transform in plain R for that reason, with MASE only storing the
result. What it offers is that the inputs to those methods, and their outputs,
stay consistent with each other and with the tissue they came from, and that a
single-assay slice can be handed to SpatialExperiment or
SpatialFeatureExperiment whenever a tool expects one.
9 See also
- Introduction to MultiAssaySpatialExperiment: a worked example from construction through aggregation
- Working with MultiAssaySpatialExperiment: construction, subsetting, annotation and aggregation
- Design of MultiAssaySpatialExperiment: slots, mapping tables and the relational schema
- MultiAssaySpatialExperiment cheatsheet: one-page API reference
10 Session info
sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> 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] stats4 stats graphics grDevices utils datasets methods
#> [8] base
#>
#> other attached packages:
#> [1] ggplot2_4.0.3 sf_1.1-2
#> [3] SingleCellExperiment_1.35.2 MultiAssaySpatialExperiment_0.99.12
#> [5] MultiAssayExperiment_1.39.1 SummarizedExperiment_1.43.0
#> [7] Biobase_2.73.2 GenomicRanges_1.65.4
#> [9] Seqinfo_1.3.2 IRanges_2.47.5
#> [11] S4Vectors_0.51.9 BiocGenerics_0.59.12
#> [13] generics_0.1.4 MatrixGenerics_1.25.0
#> [15] matrixStats_1.5.0 BiocStyle_2.41.0
#>
#> loaded via a namespace (and not attached):
#> [1] gtable_0.3.6 rjson_0.2.23 xfun_0.60
#> [4] bslib_0.12.0 lattice_0.23-1 vctrs_0.7.3
#> [7] tools_4.6.1 tibble_3.3.1 proxy_0.4-29
#> [10] pkgconfig_2.0.3 Matrix_1.7-6 KernSmooth_2.23-27
#> [13] RColorBrewer_1.1-3 S7_0.2.2 lifecycle_1.0.5
#> [16] compiler_4.6.1 farver_2.1.2 tinytex_0.60
#> [19] codetools_0.2-20 htmltools_0.5.9 class_7.3-24
#> [22] sass_0.4.10 yaml_2.3.12 pillar_1.11.1
#> [25] jquerylib_0.1.4 classInt_0.4-11 DelayedArray_0.39.6
#> [28] cachem_1.1.0 magick_2.9.1 abind_1.4-8
#> [31] tidyselect_1.2.1 digest_0.6.39 dplyr_1.2.1
#> [34] bookdown_0.48 fastmap_1.2.0 grid_4.6.1
#> [37] cli_3.6.6 SparseArray_1.13.2 magrittr_2.0.5
#> [40] S4Arrays_1.13.0 dichromat_2.0-1 e1071_1.7-17
#> [43] withr_3.0.3 scales_1.4.0 rmarkdown_2.32
#> [46] XVector_0.53.0 otel_0.2.0 SpatialExperiment_1.23.0
#> [49] evaluate_1.0.5 knitr_1.51 rlang_1.3.0
#> [52] Rcpp_1.1.2 glue_1.8.1 DBI_1.3.0
#> [55] BiocManager_1.30.27 jsonlite_2.0.0 R6_2.6.1
#> [58] units_1.0-1