Add support for spatial tiling of vector data to SpatialData

[cornhundred]

We have been thinking about spatial tiling approaches to make SpatialData objects compatible with a visualization scheme that loads tiled data at high zoom levels to enable visual exploration of large datasets (e.g., >1B transcripts), as we do with our Celldega front end (https://broadinstitute.github.io/celldega/). We have also had some initial discussions along these lines in a previous hackathon (2025-03 hackathon owkin > TileJSON for SpatialData/OMENGFF @ 💬), have read about recent work from the scverse (https://osf.io/preprints/biohackrxiv/s6bph_v1), and have reached out on a Zulip chat.

Based on our testing we have a minimal set of proposed changes to the SpatialData specification that appears to work well with our visualization spatial tiling approach. We would like to get feedback from the SpatialData community on whether this approach might make sense generally for SpatialData.

Proposed SpatialData storage/schema changes and Celldega compatibility demo

The SpatialData logical model is unchanged: the store still contains standard images, points, shapes, and an AnnData table and can be opened with spatialdata.read_zarr. The proposal changes or more narrowly specifies parts of their on-disk representation to support efficient visualization.

Area	Current SpatialData convention	Proposed profile demonstrated here
Spatial indexing	Points and shapes are stored in Parquet, but there is no standard mapping from a spatial region to particular files or row groups.	Assign points and shapes to a shared regular grid, store each logical tile in the corresponding Parquet row group, and record the grid-to-file mapping in spatial_tiling metadata (zarr.json).
Transcript columns	Point coordinates use standard x and y columns; no projection-oriented physical column order is prescribed.	Keep the standard point representation while placing commonly projected columns (x, y, and feature_name) first (for more efficient parquet_wasm column projection/filtering - note this is demonstrated with this unmerged pull request fix kylebarron/parquet-wasm#811) update this pull request has been merged.
Cell geometries	Shapes may use WKB or the already-supported GeoArrow geometry encoding; GeoArrow is not required by default.	Store cell boundaries as native geoarrow.polygon columns so Arrow-compatible clients can read them without decoding WKB.
Expression matrix	AnnData permits CSR, CSC, dense matrices, and additional layers, but SpatialData does not prescribe a gene-major representation.	Keep the conventional CSR X matrix and add the equivalent CSC matrix as layers["X_csc"] for efficient access to individual genes.

The column ordering, GeoArrow encoding, and CSC layer use capabilities already available in the underlying formats. The new convention introduced by this POC is the shared spatial tile and Parquet row-group mapping.

Please see

• repo for a proof-of-concept example visualization of a modified SpatialData using an experimental branch of Celldega: https://github.com/cornhundred/POC_SpatialData_spatial-tile_viz_with_Celldega
• example modified SpatialData https://huggingface.co/datasets/cornhundred/SpatialData_with_spatial_tiles/tree/main

Please see this repo for a Claude-coded proof-of-concept set of modifications to SpatialData-io in order to produce this new specification: scverse/spatialdata-io#423 [note that this functionality would need to be moved into SpatialData in order to provide general support for spatial transcriptomics technologies]

@jaspreetishar @LucaMarconato @keller-mark

[LucaMarconato]

Thanks again @cornhundred for starting the discussion. I'll start with a comment on the 4th point.

Expression matrix

I agree. Like you, I also think that storing the csc matrix as an additional layer would be beneficial. Another argument in favor of adding the csc representation comes from workflows for Visium HD, Stereoseq, or other similar data. See for example our Visium HD notebook (you can search for .X.tocsc()): we show how using a csc matrix for the expression information in Visium HD can speed up visualization. In the notebook we have to make such conversion manually, while having this info available by default would simplify usage.

Another place where such request was advanced is in this Zulip thread (point 3): #data-structures > anndata storage conventions @ 💬. Philipp seems to be supportive as well. I'd be good to coordinate with Philipp and Ilan and even upstream this convention to anndata.

[keller-mark]

Can you elaborate further on the spatial tiling and the required metadata? Does the parquet standard allow each row group to be a different size? What is the advantage or disadvantage of using row groups vs. multi-part parquet files (should we use one or the other or both to achieve the spatial tiling, or is there any impact/consequence on the multi-part parquet usage that SpatialData currently does)?

[cornhundred]

Hi @keller-mark, yes we are happy to elaborate. In our current implementation, the basic contract is:

• one regular spatial tile = one physical Parquet row group

Parquet row groups can contain different numbers of rows.

We also use multiple row groups per Parquet file, so row groups and multipart Parquet are solving different problems:

• Row groups are the spatial I/O unit: each transcript belongs to exactly one spatial tile/row group.
• Multipart Parquet controls the number of row groups/footer size per file: we pack a configurable number of row groups into each file and then start a new file. In our browser implementation this is useful for keeping Parquet footers manageable for parquet-wasm.

For example, the proof-of-concept modified SpatialData dataset (https://huggingface.co/datasets/cornhundred/SpatialData_with_spatial_tiles/tree/main/skin_adapt_v2.zarr) has a 171 × 82 grid = 14,022 spatial tiles/row groups, stored across 36 Parquet files with at most 400 row groups per file. The layout is recorded in the root zarr.json under attributes.spatial_tiling, including the grid definition, file list, max_row_groups_per_file, and total_row_groups.

Because the regular grid and row-group numbering are deterministic, the client can determine the storage location directly:

(x, y) -> (tile_x, tile_y) -> tile_id -> (parquet_file, local_row_group)

where:

file_index = floor(tile_id / max_row_groups_per_file)

local_row_group = tile_id % max_row_groups_per_file

This is how the current Celldega implementation works in PR #337: the viewport is converted directly to tile coordinates in visibleTiles.js, tile coordinates are mapped to row-group IDs in row_group_tile_reader.js, and in multipart mode those row-group IDs are mapped to the corresponding Parquet file and local row-group index.

We still rely on the Parquet footer/metadata to locate the byte ranges for those row groups, but we do not need to inspect row-group spatial statistics or bounds to determine which row groups intersect the viewport.

PyArrow currently caps a requested single row group at 64 * 1024 * 1024 = 67,108,864 rows. In practice this is far above anything we expect for an interactive spatial tile.

So our current preference is to use both row groups and multipart Parquet: row groups encode the spatial tiling, while multipart files keep the number of row groups/footer size per file manageable. This should fit naturally with SpatialData's existing multipart Parquet usage, with the added requirement that row-group boundaries preserve the spatial partitioning.

We are also open to Morton/Z-ordering if it can preserve the same general rendering properties. The main thing we care about is keeping the browser-side path close to:

viewport -> physical row groups -> Arrow buffers -> GPU

Ideally, a viewport query should resolve quickly to a unique set of physical row groups, each row group should be read at most once, and little or no CPU-side filtering/restructuring should be needed before rendering.

[keller-mark]

Thank you very much! I am very much in support of using the row groups to implement spatial tiling, without the need for the morton/z-ordering.