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 = ""

# Parameters
slug = "mouse-liver"
release = "2026-11-01"
out_dir = "bundles"
primers_ref = "4d7af340849c00077f9d21c5b54af9e689a4a070"
dev_note = ""

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

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)

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: mouse-liver
query: {'organism': 'mouse', 'tissue': 'liver'}
include_uids: None
exclude_uids: None
families: ['spot', 'imaging']
marker_set: liver
rationale: Mouse liver zonation in health, steatohepatitis, cancer and infection on Visium (about 20 studies), MERSCOPE and Xenium.

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)

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.")

180 kept, 9 dropped.

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.")

MIN_GENES_FOR_QC = 50: QC, and every view, needs at least this many genes.

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)

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: 155 datasets, so the 40 largest are shown.

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.")

Marker set 'liver', cited from 10.1038/s41467-018-06318-7, 10.1016/j.cell.2021.12.018.

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

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 05:17:09 | [INFO] root_attr: version

2026-11-04 05:17:09 | [INFO] root_attr: multiscales

2026-11-04 05:17:09 | [INFO] root_attr: omero

2026-11-04 05:17:09 | [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 05:17:09 | [INFO] resolution: s0

2026-11-04 05:17:09 | [INFO]  - shape ('c', 'y', 'x') = (3, 574, 600)

2026-11-04 05:17:09 | [INFO]  - chunks =  ['3', '574', '600']

2026-11-04 05:17:09 | [INFO]  - dtype = uint8

2026-11-04 05:17:09 | [INFO] root_attr: version

2026-11-04 05:17:09 | [INFO] root_attr: multiscales

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())

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']}")

spot: 420308 observations, 155 sections; 2000 highly variable genes of the 17613 measured, expressed genes (seurat_v3 on raw counts, section as batch); its 8 gene panels have 19,465 to 55,401 genes each and share 17,632.
imaging_cell:Merscope: 1168533 observations, 10 sections; all 494 measured, expressed genes of the view. Xenium 298 ∩ Merscope 494 = 91 shared genes (under 200), so one view per technology
imaging_cell:Xenium: 444157 observations, 4 sections; all 298 measured, expressed genes of the view; its 2 gene panels have 440 to 5,050 genes each and share 298. Xenium 298 ∩ Merscope 494 = 91 shared genes (under 200), so one view per technology
member krMwUHseny73yHnJ0000 dropped: panel shares 79 genes with the rest of its technology
member Faa4M6Rd9CJDCocV0000 dropped: panel shares 79 genes with the rest of its technology
member Y88qjuxR4R0ixx110000 dropped: panel shares 79 genes with the rest of its technology
member wB5aXqgrrpIg6I3t0000 dropped: panel shares 79 genes with the rest of its technology
member SMyiKkh2eTPkglFk0000 dropped: panel shares 79 genes with the rest of its technology
member HiACteU63IP5EjrN0000 dropped: panel shares 79 genes with the rest of its technology
member 9tbNAMLRngobckm80000 dropped: panel shares 79 genes with the rest of its technology
member AJJ5fmIyVK3GtKBc0000 dropped: panel would reduce 494 shared genes to 100
member tSKwcYj1CWvsOU1Z0000 dropped: panel would reduce 494 shared genes to 100
member SK9NRv677wU5RnqH0000 dropped: panel shares 43 genes with the rest of its technology
member FVBm87J4g9TB8uCW0000 dropped: panel shares 48 genes with the rest of its technology
Merscope panel of AJJ5fmIyVK3GtKBc0000, tSKwcYj1CWvsOU1Z0000: dropped: panel would reduce 494 shared genes to 100 shared genes
Xenium panel of Ee1qQj17f5M2vMme0000: kept: shares 298 of 440 shared genes
Joint view: no. fewer than two families keep a technology
  left out Visium (gate3): loses 100% of its observations below 500 counts on the 89 shared genes (limit 50%)

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."))

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)

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.

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)

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)))

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]}")

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.

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)

papermill primers/notebooks/primer.ipynb primer.ipynb -p slug mouse-liver -p release 2026-11-01 -p primers_ref 4d7af340849c00077f9d21c5b54af9e689a4a070 -p use_gpu True

	uid	title	technology	family	paper	paper_source	licence	n_obs	status
0	gdVPTd06WdCwcK1j0000	GSM7536068	Visium	spot	geo-gse236424	geo_series	ncbi-geo-policy	2281	kept
1	MNE10OQVBR2pMb020000	GSM8983432	Visium	spot	10.1038/s41586-025-09297-0	doi	CC-BY-4.0	4992	kept
2	tcfxNQYxXZjl9VZb0000	GSM8983434	Visium	spot	10.1038/s41586-025-09297-0	doi	CC-BY-4.0	4970	kept
3	JRkHlLBW5sD3WmYz0000	GSM8983431	Visium	spot	10.1038/s41586-025-09297-0	doi	CC-BY-4.0	4525	kept
4	oEjIBSePSA2Iczj40000	GSM8983433	Visium	spot	10.1038/s41586-025-09297-0	doi	CC-BY-4.0	4699	kept
...	...	...	...	...	...	...	...	...	...
184	uFyAw54uJKxv6bG60000	GSM9928801	Xenium	imaging_cell	n/a	none	ncbi-geo-policy	283108	dropped: no paper (no study_id, series or collection)
185	tSKwcYj1CWvsOU1Z0000	GSM8442438	Merscope	imaging_cell	geo-gsm8442438	study_id	ncbi-geo-policy	26037	kept
186	SMyiKkh2eTPkglFk0000	GSM9393661	Merscope	imaging_cell	10.1126/sciadv.aee8580	doi	ncbi-geo-policy	571898	kept
187	HiACteU63IP5EjrN0000	GSM9393663	Merscope	imaging_cell	10.1126/sciadv.aee8580	doi	ncbi-geo-policy	188487	kept
188	9tbNAMLRngobckm80000	GSM9393664	Merscope	imaging_cell	10.1126/sciadv.aee8580	doi	ncbi-geo-policy	411534	kept

	markers
compartment
hepatocyte	ALB, APOA1, TTR
hepatocyte_central	CYP2E1, GLUL
hepatocyte_portal	HAL, SDS
cholangiocyte	KRT19, KRT7, SOX9
endothelial	STAB2, CLEC4G, PECAM1
stellate	DCN, RELN, LRAT
kupffer	VSIG4, MARCO, CD5L, TIMD4, CLEC4F
t_nk	CD3E, NKG7
b_cell	MS4A1, CD79A

	panel	compartment	genes	note
0	Merscope-494-ea80599a	hepatocyte	Slc10a1, Baat, Hnf4a	hepatocyte genes on this panel
1	Merscope-494-ea80599a	hepatocyte_portal	Cyp2f2, Pck1	periportal hepatocyte genes on this panel
2	Xenium-477-7d85f056	hepatocyte	Hp	the only hepatocyte gene on this panel
3	Xenium-477-7d85f056	hepatocyte_central	not measured	no pericentral hepatocyte gene on this panel
4	Xenium-477-7d85f056	hepatocyte_portal	Pck1	periportal hepatocyte gene on this panel
5	Xenium-477-7d85f056	endothelial	Kdr	endothelial gene on this panel
6	Xenium-477-7d85f056	stellate	Colec11, Pdgfrb, Des	stellate genes on this panel
7	Xenium-477-7d85f056	b_cell	not measured	no B-cell gene on this panel

	count floor
technology
Visium	500.0
SpatialTranscriptomics	500.0
VisiumHD	200.0
Xenium	10.0
Merscope	10.0
MERFISH	10.0

	observations	top compartment	top mean	runner-up	runner-up mean
cluster
0	41105	hepatocyte	3.314	cholangiocyte	0.104
1	10393	hepatocyte	3.25	cholangiocyte	0.197
10	2539	hepatocyte	3.431	cholangiocyte	0.03
11	2850	hepatocyte	2.886	cholangiocyte	0.097
12	498	hepatocyte	2.831	cholangiocyte	0.343
2	20316	hepatocyte	2.905	cholangiocyte	0.03
3	43013	hepatocyte_central	3.385	cholangiocyte	0.058
4	53244	hepatocyte	2.536	cholangiocyte	0.216
5	36795	hepatocyte	2.618	cholangiocyte	0.059
6	96156	hepatocyte	3.443	cholangiocyte	0.032
7	21494	hepatocyte	3.308	cholangiocyte	0.029
8	50097	hepatocyte_central	3.589	cholangiocyte	0.022
9	41808	hepatocyte_central	4.227	cholangiocyte	0.02

Cohort primer¶

Mouse liver¶

Setup¶

Query and cohort¶

Quality control¶

Marker set¶

Spatial overview¶

Embedding¶

Biology¶

Draft reading, to be reviewed by the Lead¶

Export and reproduce¶

	status	members
0	kept	180
1	dropped: genes from another organism	7
2	dropped: no paper (no study_id, series or collection)	2

	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.

	dataset_uid	technology	n_before	n_after	median_counts	median_genes	pct_mt_median	flagged	flag_reason	status	family
0	0HoM4H8TZuT5Tz7j0000	Visium	4992	4975	13867.5	3254.0	3.851	False		kept	spot
1	0f7FT658lU9q6T7h0000	Visium	1790	1785	2798.0	1208.0	2.712	False		kept	spot
2	10kUOJx7RvBOVFWq0000	Visium	2066	2061	11647.0	3104.5	2.008	False		kept	spot
3	19T45ukApNybn0D90000	Visium	4992	4908	16472.0	4371.0	7.437	False		kept	spot
4	1BO6fIb3KbpF05oL0000	Visium	1988	1985	15316.5	3559.5	3.016	False		kept	spot
...	...	...	...	...	...	...	...	...	...	...	...
175	tU1PnP3OlsnbofIQ0000	Merscope	237950	210877	80.0	46.0	n/a	False		kept	imaging_cell
176	wB5aXqgrrpIg6I3t0000	Merscope	78942	78876	109.0	50.0	n/a	False		kept	imaging_cell
177	yWNtUL28gU5hG5K40000	Merscope	105825	95235	59.0	38.0	n/a	False		kept	imaging_cell
178	ybYk1BaXNYEpMpst0000	Xenium	162483	161614	565.0	117.0	n/a	False		kept	imaging_cell
179	zeJY9d76sQTmICTS0000	Merscope	163186	146299	92.0	53.0	n/a	False		kept	imaging_cell

	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

	view	batch	embedding	statistic	mean	random reference	limit from proportions	note
0	spot	section	before Harmony	entropy	0.142	0.956	1.453
1	spot	section	before Harmony	other_share	0.254	0.959	0.992
2	spot	section	after Harmony	entropy	0.67	0.956	1.453
3	spot	section	after Harmony	other_share	0.821	0.959	0.992
4	spot	paper	before Harmony	entropy	0.035	0.803	0.92
5	spot	paper	before Harmony	other_share	0.069	0.892	0.923
6	spot	paper	after Harmony	entropy	0.42	0.803	0.92
7	spot	paper	after Harmony	other_share	0.575	0.892	0.923
8	imaging_cell:Merscope	section	before Harmony	entropy	0.378	0.868	0.939
9	imaging_cell:Merscope	section	before Harmony	other_share	0.421	0.845	0.874
10	imaging_cell:Merscope	section	after Harmony	entropy	0.657	0.868	0.939
11	imaging_cell:Merscope	section	after Harmony	other_share	0.741	0.845	0.874
12	imaging_cell:Merscope	paper	before Harmony	entropy	n/a	n/a	n/a	one paper, so neither statistic is defined
13	imaging_cell:Merscope	paper	before Harmony	other_share	n/a	n/a	n/a	one paper, so neither statistic is defined
14	imaging_cell:Merscope	paper	after Harmony	entropy	n/a	n/a	n/a	one paper, so neither statistic is defined
15	imaging_cell:Merscope	paper	after Harmony	other_share	n/a	n/a	n/a	one paper, so neither statistic is defined
16	imaging_cell:Xenium	section	before Harmony	entropy	0.333	0.789	0.824
17	imaging_cell:Xenium	section	before Harmony	other_share	0.25	0.647	0.67
18	imaging_cell:Xenium	section	after Harmony	entropy	0.699	0.789	0.824
19	imaging_cell:Xenium	section	after Harmony	other_share	0.575	0.647	0.67
20	imaging_cell:Xenium	paper	before Harmony	entropy	0.003	0.073	0.092
21	imaging_cell:Xenium	paper	before Harmony	other_share	0.001	0.023	0.023
22	imaging_cell:Xenium	paper	after Harmony	entropy	0.026	0.072	0.092
23	imaging_cell:Xenium	paper	after Harmony	other_share	0.009	0.023	0.023

	observations	top compartment	top mean	runner-up	runner-up mean
cluster
0	28098	hepatocyte_portal	0.995	cholangiocyte	0.02
1	185426	hepatocyte_portal	2.001	cholangiocyte	0.083
2	88442	hepatocyte_portal	1.07	cholangiocyte	0.018
3	95991	hepatocyte_portal	1.373	cholangiocyte	0.023
4	68083	hepatocyte_portal	1.089	cholangiocyte	0.019
5	178417	hepatocyte_portal	1.276	cholangiocyte	0.006
6	263096	hepatocyte_portal	2.069	cholangiocyte	0.049
7	13043	hepatocyte_portal	1.059	cholangiocyte	0.011
8	247937	hepatocyte_portal	0.889	cholangiocyte	0.022

	observations	top compartment	top mean	runner-up	runner-up mean
cluster
0	42656	hepatocyte_central	1.628	hepatocyte	0.651
1	27655	hepatocyte_central	2.538	hepatocyte	0.833
2	58370	hepatocyte_central	1.961	hepatocyte	0.722
3	10031	hepatocyte_central	1.901	hepatocyte	0.862
4	14779	hepatocyte_central	1.116	hepatocyte	0.384
5	2571	hepatocyte_central	1.769	hepatocyte	0.772
6	61169	hepatocyte_central	4.189	hepatocyte	1.335
7	97994	hepatocyte_central	3.642	hepatocyte	1.422
8	128932	hepatocyte_portal	1.862	hepatocyte	1.584