fix: external metadata loading and geojson tracking

config/config.yaml

@@ -313,6 +313,10 @@ resources:
     mem_mb: 32768
     runtime: 30
     threads: 2
+  process_geojson:
+    mem_mb: 4096
+    runtime: 30
+    threads: 2
   qc_sweep:
     mem_mb: 49152
     runtime: 180

workflow/Snakefile

@@ -7,6 +7,7 @@ postprocessing → exploration). Rules are split across rules/*.smk and included
 directly below; helper functions live in rules/common.smk.
 """
 
+import glob
 import os
 import sys
 import pandas as pd
@@ -132,6 +133,17 @@ if CORE_SAMPLES is None:
         )
     SAMPLE_IDS = list(SAMPLES.keys())
 
+# Optional per-sample region GeoJSON (ROI-based analysis). Presence is optional per
+# sample; a geojson-like file that matches neither accepted name is a fatal config
+# error (see _find_geojson).
+SAMPLES_WITH_GEOJSON = [s for s in SAMPLE_IDS if _find_geojson(s) is not None]
+for _s in SAMPLE_IDS:
+    if _s not in SAMPLES_WITH_GEOJSON:
+        print(
+            f"[process_geojson] No GeoJSON found for sample '{_s}' — "
+            "no ROI-based analysis will be performed for this sample."
+        )
+
 # Optional per-column colour palettes (sample + design columns) for downstream visualization
 SAMPLE_COLORS = config.get("sample_colors", {}) or {}
 
@@ -256,6 +268,7 @@ elif MODE == "merscope":
 include: "rules/validate_input.smk"
 include: "rules/qc_sweep.smk"
 include: "rules/preprocess_umap.smk"
+include: "rules/process_geojson.smk"
 include: "rules/leiden_analysis.smk"
 include: "rules/generate_annotation_template.smk"
 include: "rules/ingest_ref.smk"
@@ -311,6 +324,7 @@ rule run_preprocessing:
         expand(rules.preprocess_umap.output.adata, sample=SAMPLE_IDS),
         expand(rules.preprocess_umap.output.metadata, sample=SAMPLE_IDS),
         expand(rules.preprocess_umap.output.report, sample=SAMPLE_IDS),
+        expand(rules.process_geojson.output.report, sample=SAMPLES_WITH_GEOJSON),
         *QUPATH_IMAGES,
         *(
             expand(

workflow/rules/annotate_cells.smk

@@ -6,6 +6,7 @@ rule annotate_cells:
     input:
         adata=_annotate_input_adata,
         metadata=rules.preprocess_umap.output.metadata,
+        external_metadata=_external_metadata_input,
         cluster_annotations=config.get("cluster_annotations", "") or [],
         spatial_niche=_annotate_niche_input,
     output:
@@ -28,7 +29,6 @@ rule annotate_cells:
         de_n_genes=ANALYSIS.get("de_n_genes", 10),
         use_precomputed=USE_PRECOMPUTED,
         external_annotation=config.get("external_annotation", {}),
-        precomputed_metadata_dir=config.get("precomputed_metadata_dir", ""),
         annotation_colors=config.get("annotation_colors", {}),
         region_colors=ANALYSIS.get("region_colors", {}),
         dpi=ANALYSIS.get("plot_dpi", 300),

workflow/rules/common.smk

@@ -195,17 +195,64 @@ def _ate_prepare_inputs(wildcards):
                 sample=wildcards.sample
             )
         )
-    # Track the region GeoJSON when present, so editing/renaming it retriggers the
-    # contract build. It is read by filename (regions are optional), not hard-required —
-    # same "declare only when present" pattern as _preprocess_inputs' precomputed_meta.
-    for _suffix in ("_he_background.geojson", "_morphology.geojson", "_he.geojson"):
-        _gj = os.path.join(GEOJ_DIR, f"{wildcards.sample}{_suffix}")
-        if os.path.isfile(_gj):
-            inputs["geojson"] = _gj
-            break
+    # Region annotations are optional. When one is present, track both the source
+    # and its validation report so the contract is rebuilt only after validation.
+    _gj = _find_geojson(wildcards.sample)
+    if _gj:
+        inputs["geojson"] = _gj
+        inputs["geojson_validation"] = _geojson_validation_report(wildcards.sample)
     return inputs
 
 
+def _geojson_suffixes():
+    """Accepted region-annotation filenames for the active head, in priority order."""
+    return {
+        "visiumhd": ("_tissue_hires_image.geojson",),
+        "xenium5k": ("_morphology.geojson",),
+        "atera": ("_he_background.geojson", "_morphology.geojson"),
+        "merscope": ("_morphology.geojson",),
+    }.get(MODE, ())
+
+
+def _find_geojson(sample):
+    """Resolve an optional region GeoJSON using filenames supported by the active head."""
+    if not GEOJ_DIR or not _geojson_suffixes():
+        return None
+
+    accepted = []
+    for _suffix in _geojson_suffixes():
+        accepted.append(os.path.join(GEOJ_DIR, f"{sample}{_suffix}"))
+
+    # A sample-prefixed GeoJSON with any other suffix is almost certainly a typo or
+    # an export for a different coordinate system. Fail before silently ignoring it.
+    discovered = sorted(glob.glob(os.path.join(GEOJ_DIR, f"{sample}_*.geojson")))
+    unsupported = [path for path in discovered if path not in accepted]
+    if unsupported:
+        sys.exit(
+            f"[config error] Unsupported GeoJSON filename(s) for sample '{sample}' "
+            f"in mode '{MODE}': {unsupported}. Expected one of: {accepted}"
+        )
+    return next((path for path in accepted if os.path.isfile(path)), None)
+
+
+def _geojson_validation_report(sample):
+    """Validation marker consumed by a head only when that sample has a GeoJSON."""
+    if _find_geojson(sample) is None:
+        return []
+    return f"{SAMPLES_DIR}/{sample}/validation/geojson_validation.json"
+
+
+def _external_metadata_input(wildcards):
+    """Tracked per-sample external metadata, or no input when the feature is off."""
+    cfg = config.get("external_annotation", {}) or {}
+    if not cfg.get("enabled", False):
+        return []
+    return os.path.join(
+        config.get("precomputed_metadata_dir", "") or "",
+        f"metadata_{wildcards.sample}.tsv",
+    )
+
+
 # ── HEAD (MERSCOPE) input helper ─────────────────────────────────────────────
 def merscope_dir_for(sample):
     """MERSCOPE region directory for a sample, from the merscope.merscope_dir
@@ -310,6 +357,7 @@ def get_all_targets(wildcards):
     targets += expand(rules.preprocess_umap.output.adata, sample=SAMPLE_IDS)
     targets += expand(rules.preprocess_umap.output.metadata, sample=SAMPLE_IDS)
     targets += expand(rules.preprocess_umap.output.report, sample=SAMPLE_IDS)
+    targets += expand(rules.process_geojson.output.report, sample=SAMPLES_WITH_GEOJSON)
     targets += list(QUPATH_IMAGES)  # head QuPath image(s) for the active mode (or [])
     if RUN_LEIDEN_ANALYSIS:
         targets += expand(

workflow/rules/prepare_input_mer.smk

@@ -8,6 +8,8 @@ rule prepare_input_mer:
     input:
         qupath_meta=rules.generate_qupath_mer.output.qupath_meta,  # px<->µm mapping + p0
         background=rules.generate_qupath_mer.output.background,  # grey contract image
+        geojson=lambda wc: _find_geojson(wc.sample) or [],
+        geojson_validation=lambda wc: _geojson_validation_report(wc.sample),
     output:
         h5ad=config["contract"]["unfiltered_h5ad"],
     log:

workflow/rules/prepare_input_vhd.smk

@@ -3,6 +3,8 @@ rule prepare_input_vhd:
     input:
         sr_done=f"{OUTDIR_SR}/{{sample}}/.done",
         hires_png=rules.generate_qupath_vhd.output.qupath_image,
+        geojson=lambda wc: _find_geojson(wc.sample) or [],
+        geojson_validation=lambda wc: _geojson_validation_report(wc.sample),
     output:
         h5ad=config["contract"]["unfiltered_h5ad"],
     log:

workflow/rules/prepare_input_x5k.smk

@@ -7,6 +7,8 @@ rule prepare_input_x5k:
     input:
         done=rules.convert_zarr_x5k.output.done,  # zarr completion marker
         qupath_meta=rules.generate_qupath_x5k.output.qupath_meta,  # geojson px→µm scale
+        geojson=lambda wc: _find_geojson(wc.sample) or [],
+        geojson_validation=lambda wc: _geojson_validation_report(wc.sample),
     output:
         h5ad=config["contract"]["unfiltered_h5ad"],
     log:

workflow/rules/preprocess_umap.smk

@@ -35,7 +35,6 @@ rule preprocess_umap:
         mito_prefix=config.get("contract", {}).get("mito_prefix", ["MT-", "mt-"]),
         random_seed=RANDOM_SEED,
         use_precomputed=USE_PRECOMPUTED,
-        precomputed_metadata_dir=config.get("precomputed_metadata_dir", ""),
         external_enabled=EXTERNAL_ENABLED,
         external_column=(config.get("external_annotation", {}) or {}).get("column", ""),
         keep_unannotated=(config.get("external_annotation", {}) or {}).get(

workflow/rules/process_geojson.smk

@@ -0,0 +1,21 @@
+rule process_geojson:
+    """Validate an optional region GeoJSON without modifying the source file."""
+    input:
+        geojson=lambda wc: _find_geojson(wc.sample),
+    output:
+        report=f"{SAMPLES_DIR}/{{sample}}/validation/geojson_validation.json",
+    log:
+        out=f"{LOGDIR}/process_geojson/{{sample}}.out",
+        err=f"{LOGDIR}/process_geojson/{{sample}}.err",
+    benchmark:
+        f"{LOGDIR}/benchmarks/process_geojson/{{sample}}.tsv"
+    conda:
+        "../envs/visiumhd.yaml"
+    threads: get_resource("process_geojson", "threads")
+    resources:
+        mem_mb=mem_mb_attempt("process_geojson"),
+        runtime=get_resource("process_geojson", "runtime"),
+    params:
+        sample_id=lambda wc: wc.sample,
+    script:
+        "../scripts/process_geojson.py"

workflow/rules/validate_input.smk

@@ -31,6 +31,7 @@ rule validate_input:
     """Assert the unfiltered h5ad satisfies the core contract; record image presence."""
     input:
         h5ad=lambda wc: _CONTRACT["unfiltered_h5ad"].format(sample=wc.sample),
+        external_metadata=_external_metadata_input,
     output:
         report=f"{_OUT}/{{sample}}/validation/input_validation.json",
     log:
@@ -51,6 +52,5 @@ rule validate_input:
         mito_prefix=_CONTRACT.get("mito_prefix", ["MT-", "mt-"]),
         external_enabled=EXTERNAL_ENABLED,
         external_column=(config.get("external_annotation", {}) or {}).get("column", ""),
-        external_meta_dir=config.get("precomputed_metadata_dir", ""),
     script:
         "../scripts/validate_input.py"

workflow/schemas/config.schema.yaml

@@ -351,7 +351,9 @@ properties:
 
   precomputed_metadata_dir:
     type: string
-    description: Directory of precomputed metadata_{sample}.tsv (with use_precomputed_clusters).
+    description: >-
+      Directory of read-only metadata_{sample}.tsv files. Cell IDs may be stored in
+      the first column or in a cell_id column.
 
   external_annotation:
     type: ["object", "null"]

workflow/scripts/annotate_cells.py

@@ -32,6 +32,7 @@
     _here = os.getcwd()
 sys.path.insert(0, _here)
 from composition_barplots import composition_pair, find_niche_column
+from external_metadata import normalized_labels, optional_input_path, read_cell_metadata
 
 
 # ── Logging ──────────────────────────────────────────────────────────────────
@@ -175,7 +176,6 @@ def generate_annotation_plots(adata, annot_key, label, plots_dir, sample_id,
 DE_N_GENES          = int(snakemake.params.de_n_genes)
 USE_PRECOMPUTED     = bool(snakemake.params.use_precomputed)
 EXT_ANNOT_CFG       = snakemake.params.external_annotation
-PRECOMPUTED_DIR     = str(snakemake.params.precomputed_metadata_dir)
 ANNOTATION_COLORS   = snakemake.params.annotation_colors
 REGION_COLORS       = snakemake.params.region_colors
 DPI          = int(getattr(snakemake.params, "dpi", 300))
@@ -184,6 +184,9 @@ def generate_annotation_plots(adata, annot_key, label, plots_dir, sample_id,
 adata_path     = str(snakemake.input.adata)
 metadata_path  = str(snakemake.input.metadata)
 annot_tsv_path = str(snakemake.input.cluster_annotations)
+external_metadata_path = optional_input_path(
+    getattr(snakemake.input, "external_metadata", None)
+)
 out_adata_path = str(snakemake.output.adata_annot)
 plots_dir      = str(snakemake.output.plots_dir)
 
@@ -209,19 +212,16 @@ def generate_annotation_plots(adata, annot_key, label, plots_dir, sample_id,
 
     # Reload precomputed clusters from metadata TSV if configured
     if USE_PRECOMPUTED:
-        _ext_meta = os.path.join(PRECOMPUTED_DIR, f"metadata_{sample_id}.tsv") if PRECOMPUTED_DIR else ""
-        if _ext_meta and os.path.isfile(_ext_meta):
-            _meta_source = _ext_meta
+        if external_metadata_path:
+            _meta_source = external_metadata_path
         elif os.path.isfile(metadata_path):
             _meta_source = metadata_path
         else:
             _meta_source = None
 
         if _meta_source:
             log.info("Reloading clusters from metadata: %s", _meta_source)
-            saved = pd.read_csv(_meta_source, sep="\t", index_col=0, comment="#")
-            saved.index = saved.index.astype(str)   # obs_names are str; numeric cell ids
-                                                    # would be inferred as int and misalign
+            saved = read_cell_metadata(_meta_source)
             leiden_cols = [c for c in saved.columns if c.startswith("leiden_")]
             for col in leiden_cols:
                 if col not in adata.obs.columns:
@@ -302,29 +302,25 @@ def generate_annotation_plots(adata, annot_key, label, plots_dir, sample_id,
         ext_col = EXT_ANNOT_CFG.get("column", "")
         if not ext_col:
             raise ValueError("external_annotation.enabled is true but no 'column' is set in config.")
-        # Resolve metadata source for external annotation
-        _ext_meta = os.path.join(PRECOMPUTED_DIR, f"metadata_{sample_id}.tsv") if PRECOMPUTED_DIR else ""
-        if _ext_meta and os.path.isfile(_ext_meta):
-            _ext_source = _ext_meta
-        elif os.path.isfile(metadata_path):
-            _ext_source = metadata_path
-        else:
-            _ext_source = None
-
-        if _ext_source is None:
+        if not external_metadata_path:
             raise FileNotFoundError(
                 f"external_annotation enabled but no metadata found for sample "
-                f"'{sample_id}' (looked in precomputed_metadata_dir and {metadata_path}).")
+                f"'{sample_id}'.")
 
-        log.info("Loading external annotation from '%s' in %s …", ext_col, _ext_source)
-        saved = pd.read_csv(_ext_source, sep="\t", index_col=0, comment="#")
-        saved.index = saved.index.astype(str)       # see note above
-        if ext_col not in saved.columns:
-            raise ValueError(
-                f"external_annotation column '{ext_col}' not found in {_ext_source} "
-                f"(sample '{sample_id}'). Available columns: {list(saved.columns)}")
+        log.info(
+            "Loading external annotation from '%s' in %s …",
+            ext_col,
+            external_metadata_path,
+        )
+        saved = read_cell_metadata(
+            external_metadata_path, required_columns=(ext_col,)
+        )
+        labels = normalized_labels(saved, ext_col)
         adata.obs["cell_type_external"] = (
-            saved[ext_col].reindex(adata.obs_names).fillna("Unannotated").astype("category"))
+            labels.reindex(adata.obs_names.astype(str))
+            .fillna("Unannotated")
+            .astype("category")
+        )
         ext_enabled = True
         log.info("  External annotation: %d types (%d cells unmatched → Unannotated)",
                  adata.obs["cell_type_external"].nunique(),