Working with MultiAssaySpatialExperiment

Introduction

This vignette covers the core operations for working with MultiAssaySpatialExperiment (MASE) objects:

  1. Construction: Building MASE objects from matrices, existing Bioconductor objects, or vendor output (see MultiAssaySpatialExperiment Use Cases)
  2. Subsetting: Filtering by specimens, assays, columns, and spatial regions
  3. Spatial operations: Annotation and aggregation by spatial regions
  4. Labels ↔︎ shapes: Converting between raster masks and vector polygons

Prerequisite: If you are new to MASE, read the Introduction to MultiAssaySpatialExperiment vignette first to understand the MASE structure and components.

library(MultiAssaySpatialExperiment)
library(SummarizedExperiment)
library(SingleCellExperiment)
library(S4Vectors)
library(sf)

Building MASE objects

There are several approaches to constructing a MASE object:

  1. Minimal MASE: Start with assays + metadata, add spatial elements later
  2. From scratch: Build all components manually
  3. From existing objects: Convert SpatialExperiment or wrap SingleCellExperiment objects in an ExperimentList
  4. From vendor output: Use platform readers (readXeniumMASE(), etc.) — covered in the MultiAssaySpatialExperiment Use Cases vignette
  5. With prepMASE: Wrap prepMultiAssay() for name harmonization and spatial FK cleanup (see below)

prepMASE and buildSpatialMap

For manual construction, buildSpatialMap() assembles a spatialMap from sampleMap (analogous to MAE’s listToMap()), and prepMASE() wraps prepMultiAssay() while harmonizing spatial slots:

spmap <- buildSpatialMap(sample_map, region = "cells", element_type = "shapes")
prepared <- prepMASE(
  ExperimentList(rna = counts, protein = protein_counts),
  colData = specimens,
  sample_map,
  shapes = ShapesLayerList(cells = cell_boundaries),
  spatialMap = spmap
)
mase <- do.call(MultiAssaySpatialExperiment, prepared)
drops(prepared$metadata$drops)  # rows removed during harmonization, if any

Optional: lazy Parquet I/O via BiocDuckDB

MultiAssaySpatialExperiment does not import BiocDuckDB. For large on-disk datasets, the BiocDuckDB package registers methods on MASE generics (readParquetForMASE, lazy spatialPoints / spatialShapes, etc.) so you can persist and query MASE objects without loading full layers into memory. See the BiocDuckDB documentation for Parquet read/write workflows.

Minimal MASE

The simplest MASE has one experiment, specimen metadata (colData), and a sampleMap. No spatial elements are required initially. In the example below, each cell is treated as its own specimen for simplicity; in practice, colData rows usually represent biological replicates or tissue sections.

# Create a simple expression matrix
counts <- matrix(rpois(50, lambda = 10), nrow = 10, ncol = 5,
                dimnames = list(paste0("Gene", 1:10),
                               paste0("Cell", 1:5)))

# Specimen metadata
specimens <- DataFrame(
  patient_id = c("P1", "P1", "P1", "P2", "P2"),
  tissue = c("cortex", "cortex", "cortex", "medulla", "medulla"),
  row.names = paste0("Cell", 1:5)
)

# Sample map: link assay columns to specimens
sample_map <- DataFrame(
  assay = factor("rna", "rna"),
  primary = paste0("Cell", 1:5),
  colname = paste0("Cell", 1:5)
)

# Construct MASE
mase_minimal <- MultiAssaySpatialExperiment(
  experiments = ExperimentList(rna = counts),
  colData = specimens,
  sampleMap = sample_map
)

mase_minimal
#> A MultiAssaySpatialExperiment object of 1 listed
#>  experiment with a user-defined name and respective class.
#>  Containing an ExperimentList class object of length 1:
#>  [1] rna: matrix with 10 rows and 5 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: 0 elements
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: NULL

Adding spatial coordinates

To add point coordinates (e.g., cell centroids), create a PointsLayerList and a spatialMap linking assay columns to spatial instances.

# Suppose we have two assays: RNA and protein
rna_counts <- matrix(rpois(50, 10), nrow = 10, ncol = 5,
                    dimnames = list(paste0("Gene", 1:10), paste0("Cell", 1:5)))
prot_counts <- matrix(rpois(15, 5), nrow = 3, ncol = 5,
                     dimnames = list(paste0("Protein", 1:3), paste0("Cell", 1:5)))

# Both assays share the same specimens
specimens2 <- DataFrame(
  patient_id = c("P1", "P1", "P1", "P2", "P2"),
  row.names = paste0("Cell", 1:5)
)

# Sample map for two assays
sample_map2 <- DataFrame(
  assay = factor(rep(c("rna", "protein"), each = 5), c("rna", "protein")),
  primary = rep(paste0("Cell", 1:5), 2),
  colname = rep(paste0("Cell", 1:5), 2)
)

# Point coordinates (shared by both assays)
coords2 <- DataFrame(
  x = runif(5, 0, 100),
  y = runif(5, 0, 100),
  instance_id = paste0("Cell", 1:5)
)

# Spatial map: link assay columns to spatial points
spatial_map2 <- DataFrame(
  assay = factor(rep(c("rna", "protein"), each = 5), c("rna", "protein")),
  colname = rep(paste0("Cell", 1:5), 2),
  element_type = "points",
  region = factor(rep("coords", 10), "coords"),
  instance_id = rep(paste0("Cell", 1:5), 2)
)

# Construct MASE with spatialMap
mase_with_map <- MultiAssaySpatialExperiment(
  experiments = ExperimentList(rna = rna_counts, protein = prot_counts),
  colData = specimens2,
  sampleMap = sample_map2,
  points = PointsLayerList(coords = coords2),
  spatialMap = spatial_map2
)

mase_with_map
#> A MultiAssaySpatialExperiment object of 2 listed
#>  experiments with user-defined names and respective classes.
#>  Containing an ExperimentList class object of length 2:
#>  [1] rna: matrix with 10 rows and 5 columns
#>  [2] protein: matrix with 3 rows and 5 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: 1 element
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: present

Adding spatial shapes

Shapes (polygons) typically represent cell boundaries, tissue regions, or annotation areas. Use the sf package to create geometries.

# Create simple circular boundaries around each point
centroids <- data.frame(
  x = runif(5, 0, 100),
  y = runif(5, 0, 100),
  instance_id = paste0("Cell", 1:5)
)

# Convert to sf points, then buffer to circles
centroids_sf <- st_as_sf(centroids, coords = c("x", "y"))
circles <- st_buffer(centroids_sf, dist = 5)

# Convert back to DataFrame with geometry column
boundaries <- DataFrame(
  geometry = st_geometry(circles),
  instance_id = paste0("Cell", 1:5)
)

# Add to MASE
mase_with_shapes <- mase_with_map
spatialShapes(mase_with_shapes) <- ShapesLayerList(boundaries = boundaries)

mase_with_shapes
#> A MultiAssaySpatialExperiment object of 2 listed
#>  experiments with user-defined names and respective classes.
#>  Containing an ExperimentList class object of length 2:
#>  [1] rna: matrix with 10 rows and 5 columns
#>  [2] protein: matrix with 3 rows and 5 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: 1 element
#>   spatialShapes: 1 element
#>   imgData: NULL
#>   spatialMap: present

Building from multiple SingleCellExperiment objects

Combine multiple SingleCellExperiment objects into a MASE:

# Create two SingleCellExperiment objects (simulating different technologies)
sce1 <- SingleCellExperiment(
  assays = list(counts = matrix(rpois(40, 10), nrow = 10, ncol = 4,
                               dimnames = list(paste0("Gene", 1:10),
                                              paste0("Cell", 1:4)))),
  colData = DataFrame(tech = "rna", cell_id = paste0("Cell", 1:4))
)

sce2 <- SingleCellExperiment(
  assays = list(counts = matrix(rpois(12, 5), nrow = 3, ncol = 4,
                               dimnames = list(paste0("Protein", 1:3),
                                              paste0("Cell", 1:4)))),
  colData = DataFrame(tech = "protein", cell_id = paste0("Cell", 1:4))
)

# Specimen metadata (shared)
specimens_sce <- DataFrame(
  patient = rep("P1", 4),
  tissue = rep("cortex", 4),
  row.names = paste0("Cell", 1:4)
)

# Sample maps
samplemap_sce <- DataFrame(
  assay = factor(rep(c("rna", "protein"), each = 4), c("rna", "protein")),
  primary = rep(paste0("Cell", 1:4), 2),
  colname = rep(paste0("Cell", 1:4), 2)
)

# Construct MASE
mase_from_sce <- MultiAssaySpatialExperiment(
  experiments = ExperimentList(rna = sce1, protein = sce2),
  colData = specimens_sce,
  sampleMap = samplemap_sce
)

mase_from_sce
#> A MultiAssaySpatialExperiment object of 2 listed
#>  experiments with user-defined names and respective classes.
#>  Containing an ExperimentList class object of length 2:
#>  [1] rna: SingleCellExperiment with 10 rows and 4 columns
#>  [2] protein: SingleCellExperiment with 3 rows and 4 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: 0 elements
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: NULL

Validation checklist

Before constructing a MASE, ensure:

  1. ExperimentList:
    • All experiments have unique column names
  2. colData:
    • Row names are specimen identifiers (primary IDs)
    • No duplicate row names
  3. sampleMap:
    • assay column matches names in ExperimentList
    • primary column matches row names in colData
    • colname column matches column names in experiments
    • No duplicate (assay, colname) pairs
  4. points/shapes (if provided):
    • instance_id column exists
    • For points: x, y columns exist
    • For shapes: geometry column exists (sfc)
  5. spatialMap (if provided):
    • assay + colname exist in sampleMap
    • element_type is “points” or “shapes”
    • region matches names in points or shapes
    • instance_id exists in the corresponding spatial element

Check validity manually:

validObject(mase_with_map)
#> [1] TRUE

Subsetting operations

MASE extends MultiAssayExperiment with spatial slots (points, shapes, images, labels, imgData, and spatialMap). When you subset a MASE, these slots are automatically updated to stay consistent with the assays and specimen metadata.

The examples below use a separate MASE object from the construction section above, so you can run them independently.

Overview of subsetting operations

Operation What it filters Spatial propagation
subsetByColData(x, y) Specimens (primary identifiers) spatialMap, imgData
subsetByRow(x, y, ...) Rows of experiments Spatial layers unchanged
subsetByColumn(x, y) Assay columns (names, indices, or logical vectors per assay) spatialMap, imgData, points, shapes, images, labels
subsetByAssay(x, y) Assays Points, shapes, images, labels, and spatialMap rows tied to retained assays
x[i, j, k, ..., drop] Row, column, assay indices Composes the above
subsetByBoundingBox(x, xmin, xmax, ymin, ymax, ...) Points/shapes by rectangle Assays via spatialMap
subsetByPolygon(x, polygon, ...) Points/shapes by polygon Assays via spatialMap

Note: subsetByAssay() keeps only spatial layers referenced in spatialMap for the retained assays. Auxiliary shape layers added for annotation (but not linked in spatialMap) are dropped — run annotateWithRegions() on the full MASE when you need those layers.

Basic subsetting

# Create a simple MASE for subsetting examples
mat <- matrix(rnorm(25), nrow = 5, ncol = 5,
  dimnames = list(paste0("G", 1:5), paste0("S", 1:5)))
ex <- SummarizedExperiment(
  assays = list(counts = mat),
  colData = DataFrame(region = rep(c("core", "margin"), length.out = 5))
)
pts <- DataFrame(
  x = c(1, 2, 3, 4, 5),
  y = c(1, 2, 3, 4, 5),
  instance_id = paste0("S", 1:5))

mase <- MultiAssaySpatialExperiment(
  experiments = ExperimentList(assay1 = ex),
  colData = DataFrame(row.names = paste0("S", 1:5)),
  sampleMap = DataFrame(
    assay = factor("assay1"),
    primary = paste0("S", 1:5),
    colname = paste0("S", 1:5)
  ),
  points = PointsLayerList(coords = pts),
  spatialMap = DataFrame(
    assay = factor("assay1"),
    colname = paste0("S", 1:5),
    element_type = "points",
    region = "coords",
    instance_id = paste0("S", 1:5)
  )
)

Subset by specimen (colData)

Filter by primary identifiers (rownames of colData):

mase[, c("S1", "S3", "S5")]
#> A MultiAssaySpatialExperiment object of 1 listed
#>  experiment with a user-defined name and respective class.
#>  Containing an ExperimentList class object of length 1:
#>  [1] assay1: SummarizedExperiment with 5 rows and 3 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: 1 element
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: present

Subset by assay columns

When subsetting by list, each element applies to the corresponding assay:

mase[, list(assay1 = c("S2", "S4"))]
#> harmonizing input:
#>   removing 3 sampleMap rows with 'colname' not in colnames of experiments
#>   removing 3 colData rownames not in sampleMap 'primary'
#> A MultiAssaySpatialExperiment object of 1 listed
#>  experiment with a user-defined name and respective class.
#>  Containing an ExperimentList class object of length 1:
#>  [1] assay1: SummarizedExperiment with 5 rows and 2 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: 1 element
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: present

Filter assay columns by observation metadata (colData of the experiment), using the standard MAE subsetByColumn list interface. This is the idiomatic replacement for ad hoc “query by metadata” patterns:

cdf <- colData(experiments(mase)[["assay1"]])
mase_core <- subsetByColumn(mase, list(assay1 = cdf$region == "core"))
#> harmonizing input:
#>   removing 2 sampleMap rows with 'colname' not in colnames of experiments
#>   removing 2 colData rownames not in sampleMap 'primary'
ncol(experiments(mase_core)[["assay1"]])
#> [1] 3
nrow(spatialPoints(mase_core)[["coords"]])
#> [1] 3

Bracket notation

The [i, j, k, drop] operator composes row, column, and assay subsetting:

mase[1:3, c("S1", "S2"), "assay1", drop = FALSE]
#> A MultiAssaySpatialExperiment object of 1 listed
#>  experiment with a user-defined name and respective class.
#>  Containing an ExperimentList class object of length 1:
#>  [1] assay1: SummarizedExperiment with 3 rows and 2 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: 1 element
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: present

Spatial subsetting

When you have spatialMap linking assay columns to spatial elements, you can subset by bounding box or polygon. Points and shapes within the region are kept; assays are filtered to retained instance_ids via spatialMap.

Subset by bounding box

Filter by a rectangle (xmin, xmax, ymin, ymax):

m2 <- subsetByBoundingBox(mase, xmin = 1.5, xmax = 4.5, ymin = 1.5, ymax = 4.5)
#> harmonizing input:
#>   removing 2 sampleMap rows with 'colname' not in colnames of experiments
#>   removing 2 colData rownames not in sampleMap 'primary'
m2
#> A MultiAssaySpatialExperiment object of 1 listed
#>  experiment with a user-defined name and respective class.
#>  Containing an ExperimentList class object of length 1:
#>  [1] assay1: SummarizedExperiment with 5 rows and 3 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: 1 element
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: present
spatialPoints(m2)[["coords"]]
#> DataFrame with 3 rows and 3 columns
#>           x         y instance_id
#>   <numeric> <numeric> <character>
#> 1         2         2          S2
#> 2         3         3          S3
#> 3         4         4          S4
ncol(experiments(m2)[["assay1"]])
#> [1] 3

Subset by polygon

Filter by an sf polygon:

poly <- st_polygon(list(matrix(
  c(1.5, 1.5, 4.5, 1.5, 4.5, 4.5, 1.5, 4.5, 1.5, 1.5),
  ncol = 2, byrow = TRUE)))
m3 <- subsetByPolygon(mase, poly)
#> harmonizing input:
#>   removing 2 sampleMap rows with 'colname' not in colnames of experiments
#>   removing 2 colData rownames not in sampleMap 'primary'
nrow(spatialPoints(m3)[["coords"]])
#> [1] 3
ncol(experiments(m3)[["assay1"]])
#> [1] 3

Multi-region subsets

Select data from multiple non-contiguous regions:

# Define two regions of interest
region1 <- st_polygon(list(matrix(c(1, 1, 2, 1, 2, 2, 1, 2, 1, 1), ncol = 2, byrow = TRUE)))
region2 <- st_polygon(list(matrix(c(4, 4, 5, 4, 5, 5, 4, 5, 4, 4), ncol = 2, byrow = TRUE)))

# Subset by each region
m_r1 <- subsetByPolygon(mase, region1)
#> harmonizing input:
#>   removing 3 sampleMap rows with 'colname' not in colnames of experiments
#>   removing 3 colData rownames not in sampleMap 'primary'
m_r2 <- subsetByPolygon(mase, region2)
#> harmonizing input:
#>   removing 3 sampleMap rows with 'colname' not in colnames of experiments
#>   removing 3 colData rownames not in sampleMap 'primary'

cat("Region 1:", ncol(experiments(m_r1)[[1]]), "observations\n")
#> Region 1: 2 observations
cat("Region 2:", ncol(experiments(m_r2)[[1]]), "observations\n")
#> Region 2: 2 observations

Buffer-based subsets

Find all points within a distance from a reference point:

# Reference point
ref_point <- st_point(c(3, 3))

# Create buffer (e.g., 1.5 units radius)
buffer_region <- st_buffer(ref_point, dist = 1.5)

# Subset by buffer
m_buffer <- subsetByPolygon(mase, buffer_region)
#> harmonizing input:
#>   removing 2 sampleMap rows with 'colname' not in colnames of experiments
#>   removing 2 colData rownames not in sampleMap 'primary'
m_buffer
#> A MultiAssaySpatialExperiment object of 1 listed
#>  experiment with a user-defined name and respective class.
#>  Containing an ExperimentList class object of length 1:
#>  [1] assay1: SummarizedExperiment with 5 rows and 3 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: 1 element
#>   spatialShapes: 0 elements
#>   imgData: NULL
#>   spatialMap: present

Spatial annotation and aggregation

Spatial omics often measures signals at discrete locations (spots, cell centroids) that lie within larger regions (cells, tissue segments). A common workflow is to:

  1. Annotate each measurement point with the region it falls in
  2. Aggregate assay values by region (sum, mean, count)

MASE provides annotateWithRegions for step 1 and aggregateByRegion for step 2. Both use the central spatialMap table to link assay columns to spatial elements.

Example setup

We construct a minimal MASE with four spots (transcriptomic measurements), coordinates, two cell polygons, and a gene expression assay. The goal is to annotate each spot with its containing cell and then aggregate expression per cell.

# Four spots at coordinates
pts_agg <- DataFrame(
  x = c(1.5, 2.5, 2.5, 3.5),
  y = c(1.5, 1.5, 2.5, 2.5),
  instance_id = paste0("S", 1:4))

# Three cell polygons
cell1 <- st_polygon(list(matrix(c(1, 1, 2, 1, 2, 2, 1, 2, 1, 1), ncol = 2, byrow = TRUE)))
cell2 <- st_polygon(list(matrix(c(2, 1, 3, 1, 3, 2, 2, 2, 2, 1), ncol = 2, byrow = TRUE)))
cell3 <- st_polygon(list(matrix(c(2, 2, 3, 2, 3, 3, 2, 3, 2, 2), ncol = 2, byrow = TRUE)))

shp_df <- DataFrame(
  instance_id = c("cell1", "cell2", "cell3"),
  geometry = st_sfc(cell1, cell2, cell3))

# Gene expression assay (3 genes × 4 spots)
expr <- matrix(c(10, 20, 5,  15,
                 30, 10, 25, 5,
                 5,  15, 20, 30),
  nrow = 3, ncol = 4,
  dimnames = list(paste0("Gene", 1:3), paste0("S", 1:4)))

# Construct MASE
mase_agg <- MultiAssaySpatialExperiment(
  experiments = ExperimentList(rna = expr),
  colData = DataFrame(row.names = paste0("S", 1:4)),
  sampleMap = DataFrame(
    assay = factor("rna"),
    primary = paste0("S", 1:4),
    colname = paste0("S", 1:4)
  ),
  points = PointsLayerList(centroids = pts_agg),
  shapes = ShapesLayerList(cells = shp_df),
  spatialMap = DataFrame(
    assay = factor("rna"),
    colname = paste0("S", 1:4),
    element_type = "points",
    region = "centroids",
    instance_id = paste0("S", 1:4)
  )
)
mase_agg
#> A MultiAssaySpatialExperiment object of 1 listed
#>  experiment with a user-defined name and respective class.
#>  Containing an ExperimentList class object of length 1:
#>  [1] rna: matrix with 3 rows and 4 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: 1 element
#>   spatialShapes: 1 element
#>   imgData: NULL
#>   spatialMap: present

Annotate with regions

annotateWithRegions performs a spatial join: for each point, it finds which shape (if any) matches and records the shape’s instance_id in a new column in spatialMap.

mase_agg_orig <- mase_agg
mase_agg <- annotateWithRegions(mase_agg, points = "centroids", shapes = "cells")
spatialMap(mase_agg)
#> DataFrame with 4 rows and 6 columns
#>      assay     colname element_type      region instance_id       cells
#>   <factor> <character>  <character> <character> <character> <character>
#> 1      rna          S1       points   centroids          S1       cell1
#> 2      rna          S2       points   centroids          S2       cell2
#> 3      rna          S3       points   centroids          S3       cell3
#> 4      rna          S4       points   centroids          S4          NA

The new column cells shows: S1 → cell1, S2 → cell2, S3 → cell3, S4 → NA (no containing cell). To assign S4 to its nearest cell instead, use join = sf::st_nearest_feature.

Aggregate by region

With the annotation in place, aggregateByRegion groups assay columns by region and applies an aggregation function.

Count points per region

FUN = "count" returns the number of points in each shape:

agg_count <- aggregateByRegion(mase_agg, by = "cells", FUN = "count")
agg_count
#> DataFrame with 3 rows and 2 columns
#>        region     count
#>   <character> <integer>
#> 1       cell1         1
#> 2       cell2         1
#> 3       cell3         1

Sum expression per region

FUN = "sum" sums each gene’s expression over spots in the same cell:

agg_sum <- aggregateByRegion(mase_agg, by = "cells", FUN = "sum")
agg_sum[["rna"]]
#>       cell1 cell2 cell3
#> Gene1    10    15    25
#> Gene2    20    30     5
#> Gene3     5    10     5

Mean expression per region

FUN = "mean" averages expression across spots in each cell:

agg_mean <- aggregateByRegion(mase_agg, by = "cells", FUN = "mean")
agg_mean[["rna"]]
#>       cell1 cell2 cell3
#> Gene1    10    15    25
#> Gene2    20    30     5
#> Gene3     5    10     5

Lower-level spatial joins

annotateWithRegions() wraps a point-in-polygon join and writes results to spatialMap. For ad hoc layer-to-layer joins, use spatialJoin() on DataFrame layers with geometry columns (see MultiAssaySpatialExperiment Cheatsheet).

Aggregation strategies

Different strategies serve different purposes:

  • Sum: Preserve total signal (RNA-seq, ATAC-seq counts)
  • Mean: Compare regions with different numbers of points (normalized comparison)
  • Count: Quality control, density analysis, normalization factor
# Compare sum vs mean for Gene1
cells <- colnames(agg_sum[["rna"]])
data.frame(
  cell = cells,
  sum = agg_sum[["rna"]]["Gene1", ],
  mean = agg_mean[["rna"]]["Gene1", ],
  n_spots = agg_count$count[match(cells, agg_count$region)]
)
#>        cell sum mean n_spots
#> cell1 cell1  10   10       1
#> cell2 cell2  15   15       1
#> cell3 cell3  25   25       1

Visualizing aggregated data

After aggregation, visualize results spatially:

# Get aggregated data for Gene1
gene1_expr <- agg_sum[["rna"]]["Gene1", ]

# Get cell polygons as sf
cells_sf <- st_sf(
  instance_id = shp_df$instance_id,
  geometry = shp_df$geometry
)

# Add expression values
cells_sf$Gene1_sum <- gene1_expr[cells_sf$instance_id]

# Plot
par(mfrow = c(1, 2))

# Original points
plot(pts_agg$x, pts_agg$y, 
     pch = 16, cex = 2,
     col = rainbow(4)[rank(expr["Gene1", ])],
     main = "Gene1: Original spots",
     xlab = "x", ylab = "y")
text(pts_agg$x, pts_agg$y, paste0("S", 1:4), pos = 3, cex = 0.7)

# Aggregated cells
plot(st_geometry(cells_sf), 
     col = rainbow(3)[rank(cells_sf$Gene1_sum)],
     main = "Gene1: Aggregated by cell",
     xlab = "x", ylab = "y")
text(st_coordinates(st_centroid(cells_sf)), 
     labels = round(cells_sf$Gene1_sum, 1),
     cex = 0.8)
#> Warning: st_centroid assumes attributes are constant over geometries

Labels ↔︎ shapes interoperability

MASE supports both labels (raster segmentation masks) and shapes (vector polygons). The package stores both slot types; converting between them is a manual workflow using external tools (see below).

  • Rasterization: shapes → labels (convert polygons to a mask)
  • Vectorization: labels → shapes (extract polygon boundaries from a mask)

Rasterization: shapes → labels

Convert cell boundary polygons to a segmentation mask:

# Example: 3 cell polygons
cell1_rast <- st_polygon(list(matrix(c(0, 0, 1, 0, 1, 1, 0, 1, 0, 0), ncol = 2, byrow = TRUE)))
cell2_rast <- st_polygon(list(matrix(c(1, 0, 2, 0, 2, 1, 1, 1, 1, 0), ncol = 2, byrow = TRUE)))
cell3_rast <- st_polygon(list(matrix(c(0, 1, 1, 1, 1, 2, 0, 2, 0, 1), ncol = 2, byrow = TRUE)))

shapes_df <- DataFrame(
  instance_id = c("cell1", "cell2", "cell3"),
  geometry = st_sfc(cell1_rast, cell2_rast, cell3_rast)
)

cat("Rasterization workflow:\n")
#> Rasterization workflow:
cat("1. Define target grid dimensions and extent\n")
#> 1. Define target grid dimensions and extent
cat("2. For each polygon, determine which pixels it covers\n")
#> 2. For each polygon, determine which pixels it covers
cat("3. Assign polygon ID to those pixels\n")
#> 3. Assign polygon ID to those pixels
cat("4. Store as RasterLayerList in MASE\n")
#> 4. Store as RasterLayerList in MASE
cat("\nRecommended packages: stars::st_rasterize, terra::rasterize\n")
#> 
#> Recommended packages: stars::st_rasterize, terra::rasterize

Vectorization: labels → shapes

Extract polygon boundaries from a segmentation mask:

# Example label matrix (3 cells)
label_matrix <- matrix(c(
  1, 1, 2, 2,
  1, 1, 2, 2,
  3, 3, 0, 0,
  3, 3, 0, 0
), nrow = 4, byrow = TRUE)

cat("Vectorization workflow:\n")
#> Vectorization workflow:
cat("1. Load label matrix (from file or RasterLayerList)\n")
#> 1. Load label matrix (from file or RasterLayerList)
cat("2. Extract polygon boundaries for each unique label value\n")
#> 2. Extract polygon boundaries for each unique label value
cat("3. Convert to sf polygons with instance_id\n")
#> 3. Convert to sf polygons with instance_id
cat("4. Store as ShapesLayerList in MASE\n")
#> 4. Store as ShapesLayerList in MASE
cat("\nRecommended packages: stars::st_as_sf, sf::st_polygonize\n")
#> 
#> Recommended packages: stars::st_as_sf, sf::st_polygonize

Use cases

Rasterization (shapes → labels):

  • Generate training data for deep learning models (requires raster masks)
  • Perform fast pixel-based operations (convolution, morphology)
  • Memory-efficient storage for dense segmentations

Vectorization (labels → shapes):

  • Perform spatial operations (area, perimeter, centroid)
  • Spatial subsetting (point-in-polygon, overlap)
  • Visualization with ggplot2 + sf

Next steps

  • MultiAssaySpatialExperiment Use Cases: Realistic workflows with data readers and spatial transformations
  • Introduction to MultiAssaySpatialExperiment: MASE structure, design principles, and terminology
  • MultiAssaySpatialExperiment Cheatsheet: Quick reference for common operations

Session info

sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04 LTS
#> 
#> Matrix products: default
#> BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so;  LAPACK version 3.12.0
#> 
#> locale:
#>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
#>  [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: Etc/UTC
#> 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-1                          
#>  [3] SingleCellExperiment_1.35.2        MultiAssaySpatialExperiment_0.99.2
#>  [5] MultiAssayExperiment_1.39.0        SummarizedExperiment_1.43.0       
#>  [7] Biobase_2.73.1                     GenomicRanges_1.65.1              
#>  [9] Seqinfo_1.3.0                      IRanges_2.47.2                    
#> [11] S4Vectors_0.51.5                   BiocGenerics_0.59.10              
#> [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.11.0             lattice_0.22-9           vctrs_0.7.3             
#>  [7] tools_4.6.1              tibble_3.3.1             proxy_0.4-29            
#> [10] pkgconfig_2.0.3          BiocBaseUtils_1.15.1     Matrix_1.7-5            
#> [13] KernSmooth_2.23-26       RColorBrewer_1.1-3       S7_0.2.2                
#> [16] lifecycle_1.0.5          compiler_4.6.1           farver_2.1.2            
#> [19] codetools_0.2-20         htmltools_0.5.9          sys_3.4.3               
#> [22] buildtools_1.0.0         class_7.3-23             sass_0.4.10             
#> [25] yaml_2.3.12              pillar_1.11.1            jquerylib_0.1.4         
#> [28] classInt_0.4-11          DelayedArray_0.39.3      cachem_1.1.0            
#> [31] magick_2.9.1             abind_1.4-8              tidyselect_1.2.1        
#> [34] digest_0.6.39            dplyr_1.2.1              maketools_1.3.2         
#> [37] fastmap_1.2.0            grid_4.6.1               cli_3.6.6               
#> [40] SparseArray_1.13.2       magrittr_2.0.5           S4Arrays_1.13.0         
#> [43] e1071_1.7-17             withr_3.0.3              scales_1.4.0            
#> [46] rmarkdown_2.31           XVector_0.53.0           otel_0.2.0              
#> [49] SpatialExperiment_1.23.0 evaluate_1.0.5           knitr_1.51              
#> [52] rlang_1.3.0              Rcpp_1.1.2               glue_1.8.1              
#> [55] DBI_1.3.0                BiocManager_1.30.27      jsonlite_2.0.0          
#> [58] R6_2.6.1                 units_1.0-1