Cohort primer¶

In [1]:
slug = "human-colon"
release = "2026-10-28"
use_gpu = True
out_dir = "bundles"
smoke = 0
primers_ref = "feat-primers"
cohorts_path = None
catalog_path = None
representatives = 1
dev_note = ""
In [2]:
# Parameters
slug = "human-pancreatic-cancer"
release = "2026-11-01"
out_dir = "bundles"
primers_ref = "4d7af340849c00077f9d21c5b54af9e689a4a070"
dev_note = ""

Human pancreatic cancer¶

Cohort human-pancreatic-cancer, catalog release 2026-11-01.

This primer describes a cohort; it is not a benchmark. It does not rank papers, technologies or methods. The batch correction (Harmony on the section) is shown for orientation only: well-mixed embeddings are not proof of a correct biology. This cohort's technologies may also differ by paper and lab. Where one paper used one technology, a technology effect can't be separated from a paper effect.

The analysis below runs the primers package at a fixed reference. To install the same code:

pip install "sddb-primers @ git+https://github.com/theislab/spatialdata-db-client@4d7af340849c00077f9d21c5b54af9e689a4a070#subdirectory=primers"

Setup¶

The imports, and a helper that links each primers function to its source at the reference that ran.

In [5]:
import inspect
import json
import platform
import time

import anndata as ad
import numpy as np
import pandas as pd
import scanpy as sc
import spatialdata_plot  # noqa: F401  (registers the .pl accessor)
import zarr
from primers.cohort import hd_table, load_cohort, smoke_spec
from primers.embed import CAVEAT, MAX_FLOOR_LOSS, MIN_SHARED_GENES, N_HVG, N_NEIGHBORS, N_PCS, SEED, embed_views
from primers.export import export_bundle, mixing_by_view, with_qc_status
from primers.markers import load_markers, panel_id
from primers.qc import MIN_GENES_FOR_QC, apply_qc, load_thresholds, qc_floors
from spatialdata.transformations import get_transformation

import sddb
from sddb.remote import storage_options

start = time.time()
out_dir = Path(out_dir)
if use_gpu:
    import cupy as cp
    import rapids_singlecell as rsc
    import rmm
    from rmm.allocators.cupy import rmm_cupy_allocator

    rmm.reinitialize(managed_memory=True)
    cp.cuda.set_allocator(rmm_cupy_allocator)


def show_call(func, what):
    """Show what a primers function does and link its source at ``primers_ref``."""
    _, line = inspect.getsourcelines(func)
    module = func.__module__.split(".")[-1]
    url = f"{REPO}/blob/{primers_ref}/primers/primers/{module}.py#L{line}"
    display(Markdown(f"`{func.__module__}.{func.__name__}`: {what} [Source at `{primers_ref}`]({url})."))


def show_table(df, digits=3):
    """Display ``df`` rounded, with missing values written as n/a, and long text columns in full."""
    shown = df.round(digits).astype(object)
    with pd.option_context("display.max_colwidth", None):
        display(shown.where(df.notna(), "n/a"))


catalog = sddb.Catalog.from_file(catalog_path) if catalog_path else sddb.Catalog(version=release)

Query and cohort¶

The cohort is a stated catalog query, plus the uids it includes and minus the uids it excludes. The spec below is the one in cohorts.yaml.

In [6]:
if smoke:
    show_call(
        smoke_spec,
        "Keeps only the smallest kept member of each of a few papers, taken in turn from each family, for a "
        "quick end-to-end check.",
    )
    n_spec_excluded = len(spec.get("exclude_uids") or [])
    spec = smoke_spec(spec, catalog, smoke)
    print(
        f"Smoke run on {smoke} papers: {len(spec['exclude_uids']) - n_spec_excluded} more uids are excluded "
        f"({len(spec['exclude_uids'])} in all)."
    )
for key in ("slug", "query", "include_uids", "exclude_uids", "families", "marker_set", "rationale"):
    value = spec.get(key)
    if smoke and key == "exclude_uids":
        value = f"{len(value)} uids"
    print(f"{key}: {value}")
slug: human-pancreatic-cancer
query: {'organism': 'human', 'tissue': 'pancreas', 'disease': ['pancreatic ductal adenocarcinoma', 'pancreatic adenocarcinoma', 'exocrine pancreatic carcinoma', 'malignant pancreatic neoplasm', 'pancreatic cancer'], 'n_obs__gt': 0}
include_uids: None
exclude_uids: None
families: ['spot', 'imaging']
marker_set: pancreas
rationale: Pancreatic ductal adenocarcinoma on Visium (one large 107-section study plus smaller ones) and Xenium; one disease context.
In [7]:
show_call(
    load_cohort,
    "Runs the query on the catalog, adds the included uids, removes the excluded ones, and loads each kept "
    "member's table keyed by gene symbol, joined within its technology family.",
)
adatas, cohort = load_cohort(spec, catalog=catalog)

primers.cohort.load_cohort: Runs the query on the catalog, adds the included uids, removes the excluded ones, and loads each kept member's table keyed by gene symbol, joined within its technology family. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

In [8]:
columns = ["uid", "title", "technology", "family", "paper", "paper_source", "license", "n_obs", "status"]
show_table(cohort[columns].rename(columns={"license": "licence"}))
show_table(cohort["status"].value_counts().rename_axis("status").reset_index(name="members"))
print(f"{int((cohort['status'] == 'kept').sum())} kept, {int((cohort['status'] != 'kept').sum())} dropped.")
uid title technology family paper paper_source licence n_obs status
0 hODmVoyNFahbWd1X0000 GSM8641107 Visium spot 10.1016/j.celrep.2025.116827 doi ncbi-geo-policy 3571 kept
1 hD3XfYwTQkHOJJ3D0000 GSM8641108 Visium spot 10.1016/j.celrep.2025.116827 doi ncbi-geo-policy 3996 kept
2 L2p8xQ6uDeIrtam50000 GSM8641074 Visium spot 10.1016/j.celrep.2025.116827 doi ncbi-geo-policy 2691 kept
3 lnP5VqQCyDjNfHHl0000 human-pancreas-preview-data-xenium-human-multi-tissue-and-cancer-panel-1-standard Xenium imaging_cell 10x Genomics demo datasets 10x_demo CC-BY-4.0 103901 kept
4 GdkN2quGouCTy2AT0000 GSM9051173 Visium spot 10.1016/j.canlet.2026.218348 doi ncbi-geo-policy 4611 kept
... ... ... ... ... ... ... ... ... ...
174 7VBtB072n0akCvCa0000 GSM8641056 Visium spot 10.1016/j.celrep.2025.116827 doi ncbi-geo-policy 2247 kept
175 uX6VZ6cou3rvFSnM0000 GSM8641039 Visium spot 10.1016/j.celrep.2025.116827 doi ncbi-geo-policy 3207 kept
176 m5ngpxS7i9eHWZdX0000 GSM8443452 Visium spot geo-gse274103 geo_series ncbi-geo-policy 4983 kept
177 1SOfUUKmmUMNFwj90000 GSM8443451 Visium spot geo-gse274103 geo_series ncbi-geo-policy 4134 kept
178 kcsfSI6N1ESWcR8P0000 GSM8452875 Visium spot geo-gse274557 geo_series ncbi-geo-policy 3177 kept

179 rows × 9 columns

status members
0 kept 178
1 dropped: no paper (no study_id, series or collection) 1
178 kept, 1 dropped.

The batch variable is the paper. A member's paper is its DOI. Without a DOI it is its GEO series, else its first collection, else its study id. A DOI that one member of a series or collection carries is given to the other members of that group. All 10x Genomics demo datasets count as one paper, because they come from one company and not from independent labs. Harmony runs on the section (dataset_uid), not on the paper. 10x often publishes one section in several Space Ranger releases, and only the newest release is kept, and the others are dropped as re-releases. Serial sections of one 10x block count as one section in the mixing statistics.

Quality control¶

Each technology group has its own filters. The table gives the values and the reason for each group.

In [9]:
show_call(load_thresholds, "Reads the QC threshold groups from thresholds.yaml.")
thresholds = load_thresholds()
show_table(
    pd.DataFrame(
        [
            {
                "group": name,
                "technologies": ", ".join(g["technologies"]),
                "min_counts": g["min_counts"],
                "min_genes": g["min_genes"],
                "max_pct_mt": g.get("max_pct_mt"),
                "flag: min median counts": g["flag"]["min_median_counts"],
                "flag: min share kept": g["flag"]["min_share_kept"],
                "reason": g["reason"],
            }
            for name, g in thresholds.items()
        ]
    ).set_index("group")
)
print(f"MIN_GENES_FOR_QC = {MIN_GENES_FOR_QC}: QC, and every view, needs at least this many genes.")

primers.qc.load_thresholds: Reads the QC threshold groups from thresholds.yaml. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

technologies min_counts min_genes max_pct_mt flag: min median counts flag: min share kept reason
group
spot Visium, SpatialTranscriptomics 500 250 30.0 1000 0.5 A Visium or ST spot holds about 1 to 10 cells, so a spot with under 500 UMIs or 250 genes is mostly background or tissue edge; solid tissues often reach 20 to 30% mitochondrial reads, so the cap is 30%.
binned_spot VisiumHD 200 50 30.0 200 0.5 On the 5 Visium HD sections of oliveira-2024-colorectal, 16 um bins at 200 counts kept 74% of bins with the same marker signal as 8 um bins at Oliveira et al.'s 100 (Nat Genet 57:1512-1523), which kept 60%, and found T cells in all 5 sections against 3 (web tasks/2026-11-01-t094-hd-bin-size.md).
imaging_cell Xenium, Merscope, MERFISH 10 3 n/a 30 0.5 Targeted panels give tens to hundreds of decoded transcripts per segmented cell, and a cell with under 10 transcripts or 3 genes cannot be given a type; the panels have no mitochondrial genes.
MIN_GENES_FOR_QC = 50: QC, and every view, needs at least this many genes.
In [10]:
show_call(
    apply_qc,
    "Leaves out control probes, computes counts, genes and the mitochondrial share per observation, keeps "
    "the observations that pass their technology's filters, and describes each dataset.",
)
qcs = []
for fam in list(adatas):
    adatas[fam], q = apply_qc(adatas[fam], thresholds)
    qcs.append(q.assign(family=fam))
qc = pd.concat(qcs, ignore_index=True)
show_table(qc)

primers.qc.apply_qc: Leaves out control probes, computes counts, genes and the mitochondrial share per observation, keeps the observations that pass their technology's filters, and describes each dataset. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

dataset_uid technology n_before n_after median_counts median_genes pct_mt_median flagged flag_reason status family
0 0877Bxvk1AFUpoD90000 Visium 2261 2153 2411.0 1336.0 0.0 False kept spot
1 0OZKzzcWGDkeLfdg0000 Visium 3249 3238 6461.0 2789.0 0.0 False kept spot
2 0fNGoWEF8xCnfaN30000 Visium 3091 3091 13164.0 4966.0 0.0 False kept spot
3 0fpFB3kKbfLq5fNx0000 Visium 2141 2104 10415.0 4319.0 0.0 False kept spot
4 0gBwEBel4oHwTfbm0000 Visium 4261 4261 13239.0 4889.0 0.0 False kept spot
... ... ... ... ... ... ... ... ... ... ... ...
173 TCU9OYSnDceV2OiJ0000 Xenium 169765 162147 128.0 105.0 n/a False kept imaging_cell
174 UFIoalxHktLSI5fo0000 Xenium 238408 224908 144.0 119.0 n/a False kept imaging_cell
175 lnP5VqQCyDjNfHHl0000 Xenium 103901 103862 109.0 39.0 n/a False kept imaging_cell
176 mzxQ72GZIk9tyOpx0000 Xenium 95675 64548 58.0 49.0 n/a False kept imaging_cell
177 n9ZxbmGR8CCbUhgN0000 Xenium 133272 132229 215.0 167.0 n/a False kept imaging_cell

178 rows × 11 columns

In [11]:
MAX_VIOLIN_DATASETS = 40
for fam, adata in adatas.items():
    obs = adata.obs[["dataset_uid", "total_counts", "n_genes"]].copy()
    sizes = obs["dataset_uid"].value_counts()
    if len(sizes) > MAX_VIOLIN_DATASETS:
        print(f"{fam}: {len(sizes)} datasets, so the {MAX_VIOLIN_DATASETS} largest are shown.")
        obs = obs[obs["dataset_uid"].isin(sizes.index[:MAX_VIOLIN_DATASETS])]
    obs["dataset_uid"] = obs["dataset_uid"].astype(str).astype("category")
    display(Markdown(f"Counts and genes per observation after QC, family `{fam}`."))
    sc.pl.violin(ad.AnnData(obs=obs), ["total_counts", "n_genes"], groupby="dataset_uid", rotation=90, stripplot=False)
spot: 169 datasets, so the 40 largest are shown.

Counts and genes per observation after QC, family spot.

No description has been provided for this image

Counts and genes per observation after QC, family imaging_cell.

No description has been provided for this image

Marker set¶

Compartments are named only from this marker set, which was fixed before any embedding was computed. A compartment that a panel misses can have a fallback chosen from that panel's genes.

In [12]:
show_call(load_markers, "Reads the marker sets from markers.yaml, each with its citation and fallbacks.")
marker_set = load_markers()[spec["marker_set"]]
print(f"Marker set {spec['marker_set']!r}, cited from {', '.join(marker_set['citation'])}.")
show_table(
    pd.DataFrame({"markers": {c: ", ".join(g) for c, g in marker_set["compartments"].items()}}).rename_axis(
        "compartment"
    )
)
kept_after_qc = set(qc.loc[qc["status"] == "kept", "dataset_uid"].astype(str))
technology_of = cohort.set_index("uid")["technology"].astype(str)
panel_of = {
    uid: panel_id(technology_of[uid], genes)
    for adata in adatas.values()
    for uid, genes in adata.uns.get("member_genes", {}).items()
    if uid in kept_after_qc
}
used = [
    {
        "panel": pid,
        "compartment": comp,
        "genes": "not measured" if fb.get("not_measured") else ", ".join(fb.get("genes") or []),
        "note": fb.get("note", ""),
    }
    for pid, entries in (marker_set.get("fallback") or {}).items()
    if pid in set(panel_of.values())
    for comp, fb in entries.items()
]
if used:
    show_table(pd.DataFrame(used))
else:
    print("No fallback of this marker set applies to the panels of this cohort.")

primers.markers.load_markers: Reads the marker sets from markers.yaml, each with its citation and fallbacks. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

Marker set 'pancreas', cited from 10.1038/s41422-019-0195-y.
markers
compartment
ductal KRT19, KRT7, SOX9, MUC1
acinar PRSS1, CPA1, CELA3A
endocrine INS, GCG, CHGA
fibroblast COL1A1, LUM, DCN
stellate_pericyte RGS5, PDGFRB
endothelial PECAM1, VWF, CDH5
t_cell CD3E, CD3D, CD2
b_plasma MS4A1, CD79A, MZB1
myeloid CD68, CD14, LYZ
panel compartment genes note
0 Xenium-377-9193d582 acinar AQP8 the only acinar gene on this panel
1 Xenium-377-9193d582 fibroblast FBLN1, SFRP2 fibroblast genes on this panel
2 Xenium-5101-5714379c acinar CTRL, SERPINI2, CUZD1 acinar genes on this panel
3 Xenium-5101-5714379c fibroblast POSTN, FAP, CXCL12 fibroblast and CAF genes on this panel

Spatial overview¶

For each technology, the kept section with the median number of observations is drawn, coloured by one marker gene of the first compartment whose markers (or its panel's fallbacks) the table measures. Only the table, the shapes it annotates and, for spot technologies, the low-resolution image are loaded.

In [13]:
show_call(with_qc_status, "Marks the members that have no cells left after QC as dropped.")
kept_members = with_qc_status(cohort, qc).query("status == 'kept'")
ENSEMBL = r"ENS[A-Z]*G\d{6}"
SYMBOL_COLUMNS = ("symbol", "gene_symbols", "gene_symbol", "gene_name", "feature_name")


def colour_gene(table, panel):
    """(var name, symbol) of the first gene of the first compartment the table measures.

    A compartment's genes are its markers, then the fallback genes for this section's own panel.
    """
    names = pd.Index(table.var_names.astype(str))
    symbols = names
    if names.str.match(ENSEMBL).mean() >= 0.5:
        column = next((c for c in SYMBOL_COLUMNS if c in table.var), None)
        if column is None:
            return None, None
        symbols = pd.Index(table.var[column].astype(str))
    by_symbol = {s.upper(): (n, s) for n, s in zip(names, symbols, strict=True)}
    fallback = (marker_set.get("fallback") or {}).get(panel) or {}
    for comp, genes in marker_set["compartments"].items():
        for gene in [*genes, *((fallback.get(comp) or {}).get("genes") or [])]:
            if gene.upper() in by_symbol:
                return by_symbol[gene.upper()]
    return None, None


def table_region(zarr_url, table_name):
    """The shapes element a table annotates, read from the table's zarr attributes without loading it."""
    group = zarr.open_group(zarr_url, mode="r", storage_options=storage_options(zarr_url) or None)
    attrs = dict(group[f"tables/{table_name}"].attrs)
    if "region" not in attrs and f"tables/{table_name}/uns/spatialdata_attrs" in group:
        attrs = ad.io.read_elem(group[f"tables/{table_name}/uns/spatialdata_attrs"])
    region = attrs.get("region")
    return region if isinstance(region, str) or region is None else (list(region) or [None])[0]


def representative_uids(members):
    """The ``representatives`` kept sections closest to the median number of observations."""
    n_obs = members.set_index("uid")["n_obs"].astype(float)
    return list((n_obs - n_obs.median()).abs().sort_values(kind="stable").index[:representatives])

primers.export.with_qc_status: Marks the members that have no cells left after QC as dropped. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

In [14]:
for tech, members in kept_members.groupby("technology"):
    for uid in representative_uids(members):
        ds = catalog.dataset(uid)
        stored = ds.elements()
        table_name = hd_table(stored) if tech == "VisiumHD" else (getattr(ds, "default_table", None) or "table")
        region = table_region(ds.zarr_url, table_name)
        if region is None or f"shapes/{region}" not in stored:
            display(Markdown(f"{tech}, section `{uid}`: the table annotates no shapes element, so no plot."))
            continue
        # a small tissue image when there is one: lowres, else hires (never the full-resolution image)
        images = [k.split("/", 1)[1] for k in stored if k.startswith("images/")]
        lowres = next((k for k in images if "lowres" in k), None) or next((k for k in images if "hires" in k), None)
        wanted = [f"tables/{table_name}", f"shapes/{region}"] + ([f"images/{lowres}"] if lowres else [])
        sdata = ds.to_sdata(elements=wanted)
        gene, symbol = colour_gene(sdata[table_name], panel_of.get(uid))
        section = members.set_index("uid").loc[uid, "title"]
        title = f"{tech}, {section}: {symbol}" if gene else f"{tech}, {section}"
        display(Markdown(f"{tech}, section `{uid}` ({int(members.set_index('uid').loc[uid, 'n_obs'])} observations)."))
        systems = list(get_transformation(sdata.shapes[region], get_all=True))
        system = "global" if "global" in systems else systems[0]
        # spatialdata-plot 0.4.2 rasterises shapes with datashader in "global" only; elsewhere use matplotlib
        method = None if system == "global" else "matplotlib"
        # fmt: off
        if lowres:
            (
                sdata
                .pl.render_images(lowres)
                .pl.render_shapes(region, color=gene, table_name=table_name, method=method, colorbar_params={"label": symbol})
                .pl.show(
                    coordinate_systems=system,
                    title=title,
                )
            )
        else:
            (
                sdata
                .pl.render_shapes(region, color=gene, table_name=table_name, method=method, colorbar_params={"label": symbol})
                .pl.show(
                    coordinate_systems=system,
                    title=title,
                )
            )
        # fmt: on
2026-11-04 04:15:46 | [INFO] root_attr: version
2026-11-04 04:15:46 | [INFO] root_attr: multiscales
2026-11-04 04:15:46 | [INFO] root_attr: omero
2026-11-04 04:15:46 | [INFO] datasets [{'path': 's0', 'coordinateTransformations': [{'type': 'scale', 'scale': [1.0, 1.0, 1.0]}, {'type': 'translation', 'translation': [0.0, 0.0, 0.0]}]}]
2026-11-04 04:15:46 | [INFO] resolution: s0
2026-11-04 04:15:46 | [INFO]  - shape ('c', 'y', 'x') = (3, 1806, 2000)
2026-11-04 04:15:46 | [INFO]  - chunks =  ['3', '1024 (+ 782)', '1024 (+ 976)']
2026-11-04 04:15:46 | [INFO]  - dtype = uint8
2026-11-04 04:15:46 | [INFO] root_attr: version
2026-11-04 04:15:46 | [INFO] root_attr: multiscales
2026-11-04 04:15:46 | [INFO] root_attr: omero

Visium, section uX6VZ6cou3rvFSnM0000 (3207 observations).

Xenium, section IlEA4ms6c9S1OvBg0005 (140702 observations).

No description has been provided for this image
No description has been provided for this image

Embedding¶

Each technology family is embedded on the genes every member measured. A joint view across families is added only when the gates pass. The parameters of the embedding come first.

In [15]:
show_call(
    qc_floors, "Gives each technology's count floor, the min_counts of its QC group, for gate 3 and the joint view."
)
floors = qc_floors(thresholds)
show_table(
    pd.DataFrame(
        {
            "value": [str(x) for x in (MIN_SHARED_GENES, MAX_FLOOR_LOSS, N_HVG, N_PCS, N_NEIGHBORS, SEED)],
            "meaning": [
                "genes a view needs before a family is split per technology, and the joint view's gate 1",
                "largest share of a technology's observations that may fall below its floor on the shared genes (gate 3)",
                "highly variable genes kept when a view has more genes",
                "principal components",
                "neighbours per observation",
                "random seed",
            ],
        },
        index=["MIN_SHARED_GENES", "MAX_FLOOR_LOSS", "N_HVG", "N_PCS", "N_NEIGHBORS", "SEED"],
    )
)
show_table(pd.Series(floors, name="count floor").rename_axis("technology").to_frame())

primers.qc.qc_floors: Gives each technology's count floor, the min_counts of its QC group, for gate 3 and the joint view. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

value meaning
MIN_SHARED_GENES 200 genes a view needs before a family is split per technology, and the joint view's gate 1
MAX_FLOOR_LOSS 0.5 largest share of a technology's observations that may fall below its floor on the shared genes (gate 3)
N_HVG 2000 highly variable genes kept when a view has more genes
N_PCS 30 principal components
N_NEIGHBORS 30 neighbours per observation
SEED 0 random seed
count floor
technology
Visium 500.0
SpatialTranscriptomics 500.0
VisiumHD 200.0
Xenium 10.0
Merscope 10.0
MERFISH 10.0
In [16]:
show_call(
    embed_views,
    "Builds one view per family (or per technology when a family shares too few genes), plus the joint view "
    "when the gates pass, and embeds each with PCA, UMAP and Harmony on the section.",
)
views = embed_views(adatas, marker_set, floor=floors, use_gpu=use_gpu)
for name, v in views.items():
    p = v.uns["primer"]
    print(f"{name}: {p['n_obs']} observations, {p['n_sections']} sections; {p['genes']['rule']}. {p['note']}".strip())
first = next(iter(views.values())).uns["primer"]
for uid, reason in first["members_dropped_panel"].items():
    print(f"member {uid} {reason}")
for c in {tuple(c["uids"]): c for v in views.values() for c in v.uns["primer"]["panel_costs"]}.values():
    print(f"{c['technology']} panel of {', '.join(c['uids'])}: {c['reason']} shared genes")
for name, reason in first["views_skipped"].items():
    print(f"no view {name}: {reason}")
gates = next((v.uns["primer"]["gates"] for v in views.values() if "gates" in v.uns["primer"]), None)
if gates is not None:
    print(f"Joint view: {'yes' if gates['joint'] else 'no'}. {gates['reason']}".strip())
    for drop in gates["drops"]:
        print(f"  left out {drop['technology']} ({drop['step']}): {drop['reason']}")

primers.embed.embed_views: Builds one view per family (or per technology when a family shares too few genes), plus the joint view when the gates pass, and embeds each with PCA, UMAP and Harmony on the section. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

spot: 524995 observations, 169 sections; 2000 highly variable genes of the 17242 measured, expressed genes (seurat_v3 on raw counts, without batch, because the per-section fit failed); its 5 gene panels have 17,943 to 36,601 genes each and share 17,242.
imaging_cell: 1220850 observations, 9 sections; all 239 measured, expressed genes of the view; its 3 gene panels have 377 to 5,101 genes each and share 239.
Joint view: no. fewer than two families keep a technology
  left out Visium (gate3): loses 92% of its observations below 500 counts on the 229 shared genes (limit 50%)
In [17]:
for name, v in views.items():
    colours = [c for c in ("paper", "technology", "tissue") if c in v.obs and v.obs[c].nunique() > 1]
    display(Markdown(f"View `{name}`, before Harmony."))
    sc.pl.umap(v, color=colours or None, ncols=1)
    if "X_umap_harmony" in v.obsm:
        display(Markdown(f"View `{name}`, after Harmony on the section."))
        sc.pl.embedding(v, basis="X_umap_harmony", color=colours or None, ncols=1)
    else:
        display(Markdown(f"View `{name}` has one section, so it has no Harmony embedding."))

View spot, before Harmony.

No description has been provided for this image

View spot, after Harmony on the section.

No description has been provided for this image

View imaging_cell, before Harmony.

No description has been provided for this image

View imaging_cell, after Harmony on the section.

No description has been provided for this image
In [18]:
show_call(
    mixing_by_view,
    "Measures how evenly the batches mix among each observation's nearest neighbours, per view, before and "
    "after Harmony; the bundle's manifest records the same numbers.",
)
mixing = mixing_by_view(views)
rows, definitions = [], {}
for name, by_batch in mixing.items():
    for label, by_basis in by_batch.items():
        for basis, m in by_basis.items():
            for stat in ("entropy", "other_share"):
                definitions.setdefault(stat, m[stat]["definition"].replace(m["batch_key"], "batch"))
                rows.append(
                    {
                        "view": name,
                        "batch": label,
                        "embedding": "after Harmony" if basis == "X_pca_harmony" else "before Harmony",
                        "statistic": stat,
                        "mean": m[stat]["mean"],
                        "random reference": m[stat]["random"],
                        "limit from proportions": m[stat]["random_limit"],
                        "note": m["note"],
                    }
                )
show_table(pd.DataFrame(rows))
for stat, text in definitions.items():
    print(f"{stat}: {text}")
print(CAVEAT)

primers.export.mixing_by_view: Measures how evenly the batches mix among each observation's nearest neighbours, per view, before and after Harmony; the bundle's manifest records the same numbers. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

view batch embedding statistic mean random reference limit from proportions note
0 spot section before Harmony entropy 0.175 0.964 1.493
1 spot section before Harmony other_share 0.265 0.96 0.994
2 spot section after Harmony entropy 0.639 0.963 1.493
3 spot section after Harmony other_share 0.768 0.96 0.994
4 spot paper before Harmony entropy 0.036 0.474 0.538
5 spot paper before Harmony other_share 0.042 0.52 0.537
6 spot paper after Harmony entropy 0.263 0.473 0.538
7 spot paper after Harmony other_share 0.336 0.519 0.537
8 imaging_cell section before Harmony entropy 0.478 0.905 0.971
9 imaging_cell section before Harmony other_share 0.561 0.848 0.875
10 imaging_cell section after Harmony entropy 0.547 0.905 0.971
11 imaging_cell section after Harmony other_share 0.589 0.848 0.875
12 imaging_cell paper before Harmony entropy 0.117 0.608 0.639
13 imaging_cell paper before Harmony other_share 0.078 0.385 0.397
14 imaging_cell paper after Harmony entropy 0.169 0.606 0.639
15 imaging_cell paper after Harmony other_share 0.136 0.385 0.397
entropy: How evenly the batch labels are spread among each observation's 30 nearest neighbours, divided by log(min(30, number of labels)): 0 when all neighbours share one label, 1 when the labels are as evenly spread as the neighbours allow.
other_share: The share of each observation's 30 nearest neighbours that carry another batch.
Well mixed is not the same as biologically right: a correction can also mix cells that differ for a real reason, such as tissue or disease.

Biology¶

Each view is clustered with Leiden (resolution 0.5, 2 iterations) on its Harmony graph, or on its PCA graph when it is uncorrected. The clusters are then described by the marker genes of the marker set only.

In [19]:
fallback_genes = {}
for entries in (marker_set.get("fallback") or {}).values():
    for comp, fb in (entries or {}).items():
        fallback_genes.setdefault(comp, set()).update(str(g).upper() for g in (fb or {}).get("genes") or [])
compartment_genes = {
    comp: {str(g).upper() for g in genes} | fallback_genes.get(comp, set())
    for comp, genes in marker_set["compartments"].items()
}
readings = {}
for name, v in views.items():
    key = "nn_harmony" if "nn_harmony_distances" in v.obsp else "nn_pca"
    leiden = {"resolution": 0.5, "random_state": SEED, "neighbors_key": key, "key_added": "leiden"}
    if use_gpu:  # rapids-singlecell's n_iterations counts aggregation levels, so its default is kept
        rsc.tl.leiden(v, **leiden)
    else:
        sc.tl.leiden(v, flavor="igraph", n_iterations=2, **leiden)
    expr = v.obsm["marker_expr"]
    m = ad.AnnData(
        X=expr.to_numpy(dtype=np.float32),
        obs=v.obs[["leiden"]].copy(),
        var=pd.DataFrame(index=expr.columns.astype(str)),
    )
    display(Markdown(f"View `{name}`: {m.obs['leiden'].nunique()} clusters on `{key}`, {m.n_vars} marker genes."))
    if m.n_vars == 0 or m.obs["leiden"].nunique() < 2:
        display(Markdown("This view has fewer than two clusters or no marker gene, so there is no dotplot."))
        continue
    by_compartment = {
        comp: [g for g in expr.columns if g.upper() in genes] for comp, genes in compartment_genes.items()
    }
    sc.pl.dotplot(m, {k: g for k, g in by_compartment.items() if g}, groupby="leiden", standard_scale="var")
    scores = (
        pd.DataFrame(
            {
                comp: expr[[g for g in expr.columns if g.upper() in genes]].mean(axis=1)
                for comp, genes in compartment_genes.items()
                if any(g.upper() in genes for g in expr.columns)
            }
        )
        .groupby(m.obs["leiden"].to_numpy())
        .mean()
    )
    ranked = scores.apply(lambda r: r.sort_values(ascending=False), axis=1)
    table = pd.DataFrame(
        {
            "observations": m.obs["leiden"].value_counts().reindex(scores.index),
            "top compartment": scores.idxmax(axis=1),
            "top mean": scores.max(axis=1),
            "runner-up": [r.dropna().index[1] if r.notna().sum() > 1 else "" for _, r in ranked.iterrows()],
            "runner-up mean": [r.dropna().iloc[1] if r.notna().sum() > 1 else np.nan for _, r in ranked.iterrows()],
        }
    ).rename_axis("cluster")
    readings[name] = table
    show_table(table)

View spot: 14 clusters on nn_harmony, 35 marker genes.

No description has been provided for this image
observations top compartment top mean runner-up runner-up mean
cluster
0 104429 fibroblast 1.812 b_plasma 0.071
1 40508 fibroblast 1.772 b_plasma 0.143
10 30947 acinar 3.498 b_plasma 0.032
11 396 fibroblast 1.772 b_plasma 0.027
12 1553 stellate_pericyte 0.54 b_plasma 0.001
13 646 fibroblast 1.049 b_plasma 0.027
2 119502 ductal 1.697 b_plasma 0.044
3 94649 fibroblast 1.775 b_plasma 0.213
4 49650 acinar 1.754 b_plasma 0.114
5 35581 fibroblast 1.666 b_plasma 0.592
6 25931 endocrine 3.198 b_plasma 0.186
7 2787 myeloid 1.023 b_plasma 0.27
8 16666 fibroblast 2.42 b_plasma 0.055
9 1750 myeloid 1.236 b_plasma 0.003

View imaging_cell: 180 clusters on nn_harmony, 15 marker genes.

No description has been provided for this image
observations top compartment top mean runner-up runner-up mean
cluster
0 378689 stellate_pericyte 0.274 b_plasma 0.002
1 153040 t_cell 0.276 b_plasma 0.004
10 15930 stellate_pericyte 0.078 b_plasma 0.003
100 212 ductal 0.0 b_plasma 0.0
101 113 ductal 0.0 b_plasma 0.0
... ... ... ... ... ...
95 167 endothelial 0.002 b_plasma 0.0
96 80 endocrine 0.003 b_plasma 0.0
97 656 ductal 0.0 b_plasma 0.0
98 218 endothelial 2.624 b_plasma 0.0
99 298 endocrine 0.001 b_plasma 0.0

180 rows × 5 columns

Draft reading, to be reviewed by the Lead¶

Each statement below is generated from the reading table of its view and cites the numbers in it. The score of a compartment is the mean, over a cluster's observations, of the mean log-normalised expression of the compartment's marker and fallback genes. Compartments can sit inside others (goblet cells are epithelial), so the top compartment of a cluster does not exclude the runner-up.

In [20]:
for name, table in readings.items():
    lines = [f"View `{name}`:"]
    for cluster, r in table.iterrows():
        text = (
            f"- Cluster {cluster} ({int(r['observations'])} observations) has its highest mean marker expression "
            f"for {r['top compartment']} ({r['top mean']:.2f})"
        )
        if r["runner-up"]:
            text += f", followed by {r['runner-up']} ({r['runner-up mean']:.2f})"
        lines.append(text + ".")
    display(Markdown("\n".join(lines)))

View spot:

  • Cluster 0 (104429 observations) has its highest mean marker expression for fibroblast (1.81), followed by b_plasma (0.07).
  • Cluster 1 (40508 observations) has its highest mean marker expression for fibroblast (1.77), followed by b_plasma (0.14).
  • Cluster 10 (30947 observations) has its highest mean marker expression for acinar (3.50), followed by b_plasma (0.03).
  • Cluster 11 (396 observations) has its highest mean marker expression for fibroblast (1.77), followed by b_plasma (0.03).
  • Cluster 12 (1553 observations) has its highest mean marker expression for stellate_pericyte (0.54), followed by b_plasma (0.00).
  • Cluster 13 (646 observations) has its highest mean marker expression for fibroblast (1.05), followed by b_plasma (0.03).
  • Cluster 2 (119502 observations) has its highest mean marker expression for ductal (1.70), followed by b_plasma (0.04).
  • Cluster 3 (94649 observations) has its highest mean marker expression for fibroblast (1.77), followed by b_plasma (0.21).
  • Cluster 4 (49650 observations) has its highest mean marker expression for acinar (1.75), followed by b_plasma (0.11).
  • Cluster 5 (35581 observations) has its highest mean marker expression for fibroblast (1.67), followed by b_plasma (0.59).
  • Cluster 6 (25931 observations) has its highest mean marker expression for endocrine (3.20), followed by b_plasma (0.19).
  • Cluster 7 (2787 observations) has its highest mean marker expression for myeloid (1.02), followed by b_plasma (0.27).
  • Cluster 8 (16666 observations) has its highest mean marker expression for fibroblast (2.42), followed by b_plasma (0.05).
  • Cluster 9 (1750 observations) has its highest mean marker expression for myeloid (1.24), followed by b_plasma (0.00).

View imaging_cell:

  • Cluster 0 (378689 observations) has its highest mean marker expression for stellate_pericyte (0.27), followed by b_plasma (0.00).
  • Cluster 1 (153040 observations) has its highest mean marker expression for t_cell (0.28), followed by b_plasma (0.00).
  • Cluster 10 (15930 observations) has its highest mean marker expression for stellate_pericyte (0.08), followed by b_plasma (0.00).
  • Cluster 100 (212 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 101 (113 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 102 (100 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 103 (103 observations) has its highest mean marker expression for ductal (2.07), followed by b_plasma (0.00).
  • Cluster 104 (127 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 105 (47 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 106 (146 observations) has its highest mean marker expression for t_cell (0.87), followed by b_plasma (0.00).
  • Cluster 107 (70 observations) has its highest mean marker expression for t_cell (0.87), followed by b_plasma (0.00).
  • Cluster 108 (173 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 109 (102 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 11 (56606 observations) has its highest mean marker expression for endothelial (0.75), followed by b_plasma (0.00).
  • Cluster 110 (362 observations) has its highest mean marker expression for endocrine (0.88), followed by b_plasma (0.00).
  • Cluster 111 (108 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 112 (1045 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 113 (42 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 114 (159 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 115 (47 observations) has its highest mean marker expression for stellate_pericyte (2.00), followed by b_plasma (0.00).
  • Cluster 116 (38 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 117 (134 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 118 (54 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 119 (193 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 12 (3964 observations) has its highest mean marker expression for t_cell (0.10), followed by b_plasma (0.00).
  • Cluster 120 (72 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 121 (311 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 122 (12 observations) has its highest mean marker expression for t_cell (0.66), followed by b_plasma (0.00).
  • Cluster 123 (13 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 124 (55 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 125 (164 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 126 (17 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 127 (130 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 128 (160 observations) has its highest mean marker expression for acinar (0.01), followed by b_plasma (0.00).
  • Cluster 129 (90 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 13 (2221 observations) has its highest mean marker expression for endocrine (0.01), followed by b_plasma (0.00).
  • Cluster 130 (49 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 131 (95 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 132 (10 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 133 (11 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 134 (57 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 135 (198 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 136 (36 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 137 (173 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 138 (34 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 139 (2000 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 14 (766 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 140 (20 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 141 (6 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 142 (6 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 143 (222 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 144 (96 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 145 (328 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 146 (52 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 147 (130 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 148 (38 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 149 (19 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 15 (2476 observations) has its highest mean marker expression for ductal (0.01), followed by b_plasma (0.00).
  • Cluster 150 (314 observations) has its highest mean marker expression for myeloid (1.32), followed by b_plasma (0.00).
  • Cluster 151 (50 observations) has its highest mean marker expression for t_cell (0.67), followed by b_plasma (0.00).
  • Cluster 152 (57 observations) has its highest mean marker expression for myeloid (1.01), followed by b_plasma (0.00).
  • Cluster 153 (12 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 154 (19 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 155 (115 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 156 (63 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 157 (141 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 158 (103 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 159 (12 observations) has its highest mean marker expression for t_cell (0.88), followed by b_plasma (0.00).
  • Cluster 16 (1327 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 160 (3518 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 161 (28 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 162 (232 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 163 (9399 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 164 (15 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 165 (6083 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 166 (66 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 167 (57 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 168 (44 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 169 (2816 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 17 (10657 observations) has its highest mean marker expression for b_plasma (0.35), followed by b_plasma (0.35).
  • Cluster 170 (5774 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 171 (74 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 172 (4181 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 173 (116 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 174 (148 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 175 (44 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 176 (48 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 177 (55 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 178 (69 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 179 (41 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 18 (1338 observations) has its highest mean marker expression for endocrine (0.03), followed by b_plasma (0.00).
  • Cluster 19 (659 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 2 (132271 observations) has its highest mean marker expression for myeloid (0.30), followed by b_plasma (0.00).
  • Cluster 20 (587 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 21 (863 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 22 (1402 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 23 (596 observations) has its highest mean marker expression for myeloid (0.10), followed by b_plasma (0.00).
  • Cluster 24 (782 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 25 (1404 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 26 (424 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 27 (639 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 28 (787 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 29 (296 observations) has its highest mean marker expression for ductal (0.31), followed by b_plasma (0.00).
  • Cluster 3 (90046 observations) has its highest mean marker expression for acinar (0.77), followed by b_plasma (0.01).
  • Cluster 30 (850 observations) has its highest mean marker expression for stellate_pericyte (0.04), followed by b_plasma (0.00).
  • Cluster 31 (2077 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 32 (908 observations) has its highest mean marker expression for t_cell (0.02), followed by b_plasma (0.00).
  • Cluster 33 (624 observations) has its highest mean marker expression for stellate_pericyte (0.43), followed by b_plasma (0.00).
  • Cluster 34 (1176 observations) has its highest mean marker expression for t_cell (0.00), followed by b_plasma (0.00).
  • Cluster 35 (189 observations) has its highest mean marker expression for endothelial (0.02), followed by b_plasma (0.00).
  • Cluster 36 (721 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 37 (1357 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 38 (515 observations) has its highest mean marker expression for ductal (0.01), followed by b_plasma (0.00).
  • Cluster 39 (759 observations) has its highest mean marker expression for myeloid (0.00), followed by b_plasma (0.00).
  • Cluster 4 (37032 observations) has its highest mean marker expression for endocrine (0.97), followed by b_plasma (0.00).
  • Cluster 40 (1321 observations) has its highest mean marker expression for ductal (0.01), followed by b_plasma (0.00).
  • Cluster 41 (677 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 42 (1557 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 43 (600 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 44 (662 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 45 (640 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 46 (645 observations) has its highest mean marker expression for endothelial (0.04), followed by b_plasma (0.00).
  • Cluster 47 (794 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 48 (759 observations) has its highest mean marker expression for t_cell (0.80), followed by b_plasma (0.00).
  • Cluster 49 (148 observations) has its highest mean marker expression for myeloid (0.01), followed by b_plasma (0.00).
  • Cluster 5 (83973 observations) has its highest mean marker expression for acinar (0.19), followed by b_plasma (0.00).
  • Cluster 50 (1383 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 51 (599 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 52 (285 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 53 (520 observations) has its highest mean marker expression for myeloid (0.00), followed by b_plasma (0.00).
  • Cluster 54 (1041 observations) has its highest mean marker expression for stellate_pericyte (2.56), followed by b_plasma (0.00).
  • Cluster 55 (719 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 56 (205 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 57 (712 observations) has its highest mean marker expression for stellate_pericyte (0.01), followed by b_plasma (0.00).
  • Cluster 58 (720 observations) has its highest mean marker expression for acinar (0.00), followed by b_plasma (0.00).
  • Cluster 59 (502 observations) has its highest mean marker expression for myeloid (0.02), followed by b_plasma (0.00).
  • Cluster 6 (81104 observations) has its highest mean marker expression for ductal (0.76), followed by b_plasma (0.00).
  • Cluster 60 (224 observations) has its highest mean marker expression for ductal (0.05), followed by b_plasma (0.00).
  • Cluster 61 (555 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 62 (440 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 63 (58 observations) has its highest mean marker expression for ductal (0.02), followed by b_plasma (0.00).
  • Cluster 64 (821 observations) has its highest mean marker expression for endocrine (0.04), followed by b_plasma (0.00).
  • Cluster 65 (221 observations) has its highest mean marker expression for myeloid (0.02), followed by b_plasma (0.00).
  • Cluster 66 (160 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 67 (2258 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 68 (186 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 69 (185 observations) has its highest mean marker expression for t_cell (0.82), followed by b_plasma (0.00).
  • Cluster 7 (18798 observations) has its highest mean marker expression for t_cell (0.14), followed by b_plasma (0.01).
  • Cluster 70 (291 observations) has its highest mean marker expression for endothelial (2.56), followed by b_plasma (0.00).
  • Cluster 71 (168 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 72 (42 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 73 (65 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 74 (171 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 75 (429 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 76 (306 observations) has its highest mean marker expression for ductal (2.63), followed by b_plasma (0.00).
  • Cluster 77 (257 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 78 (362 observations) has its highest mean marker expression for myeloid (1.31), followed by b_plasma (0.00).
  • Cluster 79 (196 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 8 (19259 observations) has its highest mean marker expression for b_plasma (0.35), followed by b_plasma (0.35).
  • Cluster 80 (341 observations) has its highest mean marker expression for stellate_pericyte (0.01), followed by b_plasma (0.00).
  • Cluster 81 (162 observations) has its highest mean marker expression for myeloid (1.29), followed by b_plasma (0.00).
  • Cluster 82 (180 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 83 (87 observations) has its highest mean marker expression for endothelial (0.01), followed by b_plasma (0.00).
  • Cluster 84 (151 observations) has its highest mean marker expression for acinar (0.00), followed by b_plasma (0.00).
  • Cluster 85 (182 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 86 (163 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 87 (82 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 88 (598 observations) has its highest mean marker expression for stellate_pericyte (0.00), followed by b_plasma (0.00).
  • Cluster 89 (130 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 9 (45815 observations) has its highest mean marker expression for stellate_pericyte (0.29), followed by b_plasma (0.00).
  • Cluster 90 (195 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 91 (183 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 92 (59 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 93 (65 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 94 (369 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 95 (167 observations) has its highest mean marker expression for endothelial (0.00), followed by b_plasma (0.00).
  • Cluster 96 (80 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).
  • Cluster 97 (656 observations) has its highest mean marker expression for ductal (0.00), followed by b_plasma (0.00).
  • Cluster 98 (218 observations) has its highest mean marker expression for endothelial (2.62), followed by b_plasma (0.00).
  • Cluster 99 (298 observations) has its highest mean marker expression for endocrine (0.00), followed by b_plasma (0.00).

Export and reproduce¶

The bundle holds the manifest, the web subsample of each view, the QC and cohort tables and the embeddings.

In [21]:
show_call(
    export_bundle,
    "Writes the manifest, umap.parquet, qc.parquet, cohort.parquet and embedding.h5ad for this cohort. A cohort whose embedding is over 1 GiB gets parts of whole datasets instead (embedding-1.h5ad, embedding-2.h5ad, ...), and anndata.concat over the parts gives the whole cohort.",
)
bundle = export_bundle(
    views,
    cohort,
    qc,
    spec,
    out_dir,
    release,
    runtime_s=time.time() - start,
    mixing=mixing,
    primers_ref=primers_ref,
)
for f in sorted(bundle.iterdir()):
    print(f"{f.relative_to(out_dir)}  {f.stat().st_size} bytes")
manifest = json.loads((bundle / "manifest.json").read_text())
for key in (
    "slug",
    "release",
    "primers_ref",
    "n_datasets",
    "n_papers",
    "technologies",
    "n_cells_before",
    "n_cells_after",
    "n_flagged",
    "path",
    "hardware",
    "runtime_s",
):
    print(f"{key}: {manifest[key]}")

primers.export.export_bundle: Writes the manifest, umap.parquet, qc.parquet, cohort.parquet and embedding.h5ad for this cohort. A cohort whose embedding is over 1 GiB gets parts of whole datasets instead (embedding-1.h5ad, embedding-2.h5ad, ...), and anndata.concat over the parts gives the whole cohort. Source at 4d7af340849c00077f9d21c5b54af9e689a4a070.

ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
ImplicitModificationWarning: Trying to modify attribute `.obs` of view, initializing view as actual.
human-pancreatic-cancer/cohort.parquet  14468 bytes
human-pancreatic-cancer/embedding.h5ad  370794625 bytes
human-pancreatic-cancer/manifest.json  130094 bytes
human-pancreatic-cancer/primer.ipynb  2396447 bytes
human-pancreatic-cancer/qc.parquet  14636 bytes
human-pancreatic-cancer/umap.parquet  2932669 bytes
slug: human-pancreatic-cancer
release: 2026-11-01
primers_ref: 4d7af340849c00077f9d21c5b54af9e689a4a070
n_datasets: 178
n_papers: 14
technologies: ['Visium', 'Xenium']
n_cells_before: 1851600
n_cells_after: 1745845
n_flagged: 5
path: gpu
hardware: NVIDIA B300 SXM6 AC
runtime_s: 1058.8
In [22]:
command = (
    f"papermill primers/notebooks/primer.ipynb primer.ipynb -p slug {slug} -p release {release} "
    f"-p primers_ref {primers_ref} -p use_gpu {use_gpu}" + (f" -p smoke {smoke}" if smoke else "")
)
versions = "\n".join(f"  - {k} {v}" for k, v in manifest["versions"].items())
display(
    Markdown(
        f"- Catalog release: `{release}`\n- primers reference: `{primers_ref}`\n- Versions:\n{versions}\n"
        f"- Hardware: {manifest['hardware']} ({platform.system()})\n- Path: {manifest['path']}\n\n"
        "To run this primer again:"
    )
)
show_command(command)
  • Catalog release: 2026-11-01
  • primers reference: 4d7af340849c00077f9d21c5b54af9e689a4a070
  • Versions:
    • spatialdata-db 0.1.3
    • scanpy 1.12.5
    • anndata 0.13.4
    • spatialdata 0.8.0
    • spatialdata-plot 0.4.2
    • rapids-singlecell-cu13 0.18.1
  • Hardware: NVIDIA B300 SXM6 AC (Linux)
  • Path: gpu

To run this primer again:

papermill primers/notebooks/primer.ipynb primer.ipynb -p slug human-pancreatic-cancer -p release 2026-11-01 -p primers_ref 4d7af340849c00077f9d21c5b54af9e689a4a070 -p use_gpu True