diff --git a/.github/workflows/prepare_test_data.yaml b/.github/workflows/prepare_test_data.yaml index 34d15cca..ee14a1b8 100644 --- a/.github/workflows/prepare_test_data.yaml +++ b/.github/workflows/prepare_test_data.yaml @@ -66,6 +66,13 @@ jobs: # OMAP10 for format v0.x.x curl -o OMAP10_small.zip "https://zenodo.org/api/records/18196366/files-archive" + # ------- + # nf-core/mcmicro 2.0.0 `test` profile output (WSI, single sample, mccellpose segmentation). + # Derived from nf-core/test-datasets:mcmicro (lunaphore_example, MIT-licensed; source data + # Wuennemann et al. 2024, syn51471704, CC BY 4.0). Reproduce with: + # nextflow run nf-core/mcmicro -r 2.0.0 -profile test,docker --outdir out + curl -L -o mcmicro-nfcore-wsi-test.zip "https://github.com/mc2-center/spatialdata-io/releases/download/test-data/mcmicro-nfcore-wsi-test.zip" + - name: Unzip files run: | cd ./data diff --git a/src/spatialdata_io/__main__.py b/src/spatialdata_io/__main__.py index f0a5f9f7..4412fb77 100644 --- a/src/spatialdata_io/__main__.py +++ b/src/spatialdata_io/__main__.py @@ -252,6 +252,18 @@ def iss_wrapper( "--input", "-i", type=click.Path(exists=True), help="Path to the mcmicro project directory.", required=True ) @click.option("--output", "-o", type=click.Path(), help="Path to the output.zarr file.", required=True) +@click.option( + "--pipeline", + type=click.Choice(["labsyspharm", "nfcore"]), + default=None, + help="Which mcmicro pipeline produced the output. [default: auto-detect]", +) +@click.option( + "--markers-file", + type=click.Path(), + default=None, + help="Path to the markers CSV. [default: auto-detect]", +) @click.option( "--imread-kwargs", type=str, @@ -273,6 +285,8 @@ def iss_wrapper( def mcmicro_wrapper( input: str, output: str, + pipeline: str | None = None, + markers_file: str | None = None, imread_kwargs: str = "{}", image_models_kwargs: str = "{}", labels_models_kwargs: str = "{}", @@ -282,6 +296,8 @@ def mcmicro_wrapper( sdata = mcmicro( input, + pipeline=pipeline, + markers_file=markers_file, imread_kwargs=_parse_json_param(imread_kwargs, "imread_kwargs"), image_models_kwargs=_parse_json_param(image_models_kwargs, "image_models_kwargs"), labels_models_kwargs=_parse_json_param(labels_models_kwargs, "labels_models_kwargs"), diff --git a/src/spatialdata_io/_constants/_constants.py b/src/spatialdata_io/_constants/_constants.py index 46f983b1..8337a692 100644 --- a/src/spatialdata_io/_constants/_constants.py +++ b/src/spatialdata_io/_constants/_constants.py @@ -275,11 +275,19 @@ class StereoseqKeys(ModeEnum): RESOLUTION = "resolution" +@unique +class McmicroPipeline(ModeEnum): + """The two mcmicro pipeline flavors whose output layouts the reader understands.""" + + LABSYSPHARM = "labsyspharm" # original labsyspharm/mcmicro Nextflow pipeline + NFCORE = "nfcore" # nf-core/mcmicro pipeline + + @unique class McmicroKeys(ModeEnum): """Keys for *Mcmicro* formatted dataset.""" - # files and directories + # files and directories (labsyspharm/mcmicro layout) QUANTIFICATION_DIR = "quantification" MARKERS_FILE = "markers.csv" IMAGES_DIR_WSI = "registration" @@ -292,7 +300,21 @@ class McmicroKeys(ModeEnum): COREOGRAPH_CENTROIDS = "qc/coreograph/centroidsY-X.txt" COREOGRAPH_TMA_MAP = "qc/coreograph/TMA_MAP.tif" + # files and directories (nf-core/mcmicro layout) + NFCORE_IMAGES_DIR_WSI = "registration/ashlar" + NFCORE_IMAGES_DIR_TMA = "tma_dearray" + NFCORE_QUANTIFICATION_DIR = "quantification/mcquant" + NFCORE_COREOGRAPH_CENTROIDS = "tma_dearray/centroidsY-X.txt" + NFCORE_PIPELINE_INFO = "pipeline_info" + NFCORE_PARAMS_GLOB = "params*.json" + NFCORE_BACKSUB_DIR = "backsub" + NFCORE_BACKSUB_MARKERS_GLOB = "*_backsub.csv" + NFCORE_MARKERSHEET_GLOBS = "summary/*markersheet*.tsv;prelude/*markersheet*.tsv" + NFCORE_TMA_MAP_STEM = "TMA_MAP" + # metadata + MARKER_NAME = "marker_name" + CHANNEL_NUMBER = "channel_number" COORDS_X = "X_centroid" COORDS_Y = "Y_centroid" INSTANCE_KEY = "CellID" diff --git a/src/spatialdata_io/readers/mcmicro.py b/src/spatialdata_io/readers/mcmicro.py index 68171da1..0a0a9702 100644 --- a/src/spatialdata_io/readers/mcmicro.py +++ b/src/spatialdata_io/readers/mcmicro.py @@ -1,6 +1,8 @@ from __future__ import annotations +import json import re +import warnings from pathlib import Path from types import MappingProxyType from typing import TYPE_CHECKING, Any @@ -16,8 +18,8 @@ from spatialdata.transformations import Identity, Translation, set_transformation from yaml.loader import SafeLoader -from spatialdata_io._constants._constants import McmicroKeys -from spatialdata_io.readers._utils._utils import _set_reader_metadata +from spatialdata_io._constants._constants import McmicroKeys, McmicroPipeline +from spatialdata_io.readers._utils._utils import _set_reader_metadata, parse_channels if TYPE_CHECKING: from collections.abc import Mapping @@ -48,20 +50,36 @@ def _get_transformation( def mcmicro( path: str | Path, + pipeline: str | McmicroPipeline | None = None, + markers_file: str | Path | None = None, imread_kwargs: Mapping[str, Any] = MappingProxyType({}), image_models_kwargs: Mapping[str, Any] = MappingProxyType({}), labels_models_kwargs: Mapping[str, Any] = MappingProxyType({}), ) -> SpatialData: """Read a *Mcmicro* output into a SpatialData object. + Supports both the original `labsyspharm/mcmicro `_ Nextflow pipeline + and the newer `nf-core/mcmicro `_ pipeline. The two produce + different output directory layouts; the reader auto-detects which one it is looking at, + and either whole-slide-image (WSI) or tissue-microarray (TMA) mode. Multiple samples and + multiple segmentation methods per run are supported. + .. seealso:: - - `Mcmicro pipeline `_. + - `Mcmicro pipeline `_. + - `nf-core/mcmicro pipeline `_. Parameters ---------- path - Path to the dataset. + Path to the dataset (the mcmicro output directory). + pipeline + Which pipeline produced the output, ``"labsyspharm"`` or ``"nfcore"``. When ``None`` + (default) the pipeline is auto-detected from the directory contents. + markers_file + Optional path to the markers CSV (with a ``marker_name`` column). When ``None`` the + reader looks in the conventional locations for the detected pipeline, falling back to + the channel names embedded in the registration OME-TIFF. imread_kwargs Keyword arguments to pass to the image reader. image_models_kwargs @@ -74,65 +92,70 @@ def mcmicro( :class:`spatialdata.SpatialData` """ path = Path(path) - params = _load_params(path) - tma: bool = params["workflow"]["tma"] - tma_centroids: pd.DataFrame | None = None + pipeline = _detect_pipeline(path) if pipeline is None else McmicroPipeline(pipeline) - markers = pd.read_csv(path / McmicroKeys.MARKERS_FILE) - markers.index = markers.marker_name - assert markers.channel_number.is_monotonic_increasing - marker_names = markers.marker_name.tolist() + tma = _is_tma(path, pipeline) - images = {} + registration_dir = _registration_dir(path, pipeline) + if not registration_dir.exists(): + raise ValueError(f"{path} does not contain a `{registration_dir.name}` registration directory") + + registration_samples = sorted(registration_dir.glob("*" + McmicroKeys.IMAGE_SUFFIX)) + if not registration_samples: + raise ValueError(f"No `{McmicroKeys.IMAGE_SUFFIX}` images found in {registration_dir}") + + # sample ids are the registration image stems; used to attribute masks / tables to a sample. + sample_ids = [_image_stem(s) for s in registration_samples] + + # markers DataFrame (labsyspharm) and channel names for the registration images + markers_df, image_channel_names = _load_markers( + path, pipeline, markers_file, registration_samples[0], imread_kwargs + ) + + tma_centroids: pd.DataFrame | None = None if tma: - centroids_file = path / McmicroKeys.COREOGRAPH_CENTROIDS + centroids_file = _centroids_file(path, pipeline) tma_centroids = pd.read_csv(centroids_file, header=None, names=["y", "x"], index_col=False, sep=" ") tma_centroids.index = tma_centroids.index + 1 - image_dir = path / McmicroKeys.IMAGES_DIR_WSI - if not image_dir.exists(): - raise ValueError(f"{path} does not contain {McmicroKeys.IMAGES_DIR_WSI} directory") - - samples = list(image_dir.glob("*" + McmicroKeys.IMAGE_SUFFIX)) - if len(samples) > 1: - raise ValueError("Only one sample per dataset is supported.") - - # if tma is true, the image in `samples` is the global image with all the tma cores; let's use it and then - # reassign `samples` to each individual core + images = {} if tma: - assert len(samples) == 1 - data = imread(samples[0], **imread_kwargs) - # , scale_factors=[2, 2] + # the registration image is the global mosaic containing all cores; keep it as `tma_map` + # and then parse each individual core from the dearray directory. + data = imread(registration_samples[0], **imread_kwargs) data = Image2DModel.parse(data, transformations=_get_transformation(), rgb=None, **image_models_kwargs) images["tma_map"] = data - image_dir = path / McmicroKeys.IMAGES_DIR_TMA - samples = list(image_dir.glob("*" + McmicroKeys.IMAGE_SUFFIX)) - image_dir_masks = image_dir / "masks" - samples_masks = list(image_dir_masks.glob("*")) + dearray_dir = _dearray_dir(path, pipeline) + # per-core images have integer stems (1.tif, 2.tif, ...); skip TMA_MAP.tif and any extras. + all_dearray = sorted(set(dearray_dir.glob("*" + McmicroKeys.IMAGE_SUFFIX)) | set(dearray_dir.glob("*.tif"))) + samples = sorted((s for s in all_dearray if _image_stem(s).isdigit()), key=lambda s: int(_image_stem(s))) + samples_masks = sorted((dearray_dir / "masks").glob("*")) + else: + samples = registration_samples + samples_masks = [] + + core_ids = [_image_stem(s) for s in samples] if tma else [] for sample in samples: - core_id: str = sample.with_name(sample.stem).with_suffix("").stem - if tma: - image_id = f"core_{core_id}" - else: - image_id = core_id + core_id = _image_stem(sample) + image_id = f"core_{core_id}" if tma else core_id data = imread(sample, **imread_kwargs) - data = Image2DModel.parse(data, c_coords=marker_names, rgb=None, **image_models_kwargs) + c_coords = image_channel_names if len(image_channel_names) == data.shape[0] else None + data = Image2DModel.parse(data, c_coords=c_coords, rgb=None, **image_models_kwargs) transformations = _get_transformation( tma=int(core_id) if tma else None, tma_centroids=tma_centroids, raster_data=data ) set_transformation(data, transformation=transformations, set_all=True) images[f"{image_id}_image"] = data - # in exemplar-001 the raw images are aligned with the illumination images, not with the registration image + # in exemplar-001 the raw images are aligned with the illumination images, not with the + # registration image. These directories only exist for the labsyspharm layout. raw_dir = path / McmicroKeys.RAW_DIR if raw_dir.exists(): - raw_images = list(raw_dir.glob("*")) - for raw_image in raw_images: - raw_name = raw_image.with_name(raw_image.stem).with_suffix("").stem - + for raw_image in raw_dir.glob("*"): + raw_name = _image_stem(raw_image) data = imread(raw_image, **imread_kwargs) images[raw_name] = Image2DModel.parse( data, transformations={raw_name: Identity()}, rgb=None, **image_models_kwargs @@ -140,58 +163,98 @@ def mcmicro( illumination_dir = path / McmicroKeys.ILLUMINATION_DIR if illumination_dir.exists(): - illumination_images = list(illumination_dir.glob("*")) - for illumination_image in illumination_images: - illumination_name = illumination_image.with_name(illumination_image.stem).with_suffix("").stem + for illumination_image in illumination_dir.glob("*"): + illumination_name = _image_stem(illumination_image) raw_name = illumination_name.removesuffix(McmicroKeys.ILLUMINATION_SUFFIX_DFP) raw_name = raw_name.removesuffix(McmicroKeys.ILLUMINATION_SUFFIX_FFP) - data = imread(illumination_image, **imread_kwargs) images[illumination_name] = Image2DModel.parse( data, transformations={raw_name: Identity()}, rgb=None, **image_models_kwargs ) - samples_labels = list((path / McmicroKeys.LABELS_DIR).glob("*/*" + McmicroKeys.IMAGE_SUFFIX)) - - labels = {} - for labels_path in samples_labels: - if not tma: - # TODO: when support python >= 3.9 chance to str.removesuffix(McmicroKeys.IMAGE_SUFFIX) - segmentation_stem = labels_path.with_name(labels_path.stem).with_suffix("").stem - - data = imread(labels_path, **imread_kwargs).squeeze() - data = Labels2DModel.parse(data, transformations=_get_transformation(), **labels_models_kwargs) - labels[f"{image_id}_{segmentation_stem}"] = data - else: - segmentation_stem = labels_path.with_name(labels_path.stem).with_suffix("").stem - core_id_search = re.search(r"\d+$", labels_path.parent.stem) - if core_id_search is None: - raise ValueError(f"Cannot infer core_id from {labels_path.parent}") - else: - core_id = core_id_search.group() - assert core_id is not None + labels = _get_labels( + path, + pipeline, + tma=tma, + tma_centroids=tma_centroids, + sample_ids=sample_ids, + core_ids=core_ids, + wsi_image_id=sample_ids[0] if not tma else None, + samples_masks=samples_masks, + imread_kwargs=imread_kwargs, + labels_models_kwargs=labels_models_kwargs, + ) - data = imread(labels_path, **imread_kwargs).squeeze() - data = Labels2DModel.parse(data, **labels_models_kwargs) - transformations = _get_transformation(tma=int(core_id), tma_centroids=tma_centroids, raster_data=data) - set_transformation(data, transformation=transformations, set_all=True) - labels[f"core_{core_id}_{segmentation_stem}"] = data + tables_dict = _get_tables(path, pipeline, markers_df, tma, sample_ids, core_ids, list(labels)) - if tma: - for mask_path in samples_masks: - mask_stem = mask_path.stem - core_id = mask_stem.split("_")[0] + sdata = SpatialData(images=images, labels=labels, tables=tables_dict) + return _set_reader_metadata(sdata, "mcmicro") - data = imread(mask_path, **imread_kwargs).squeeze() - data = Labels2DModel.parse(data, **labels_models_kwargs) - transformations = _get_transformation(tma=int(core_id), tma_centroids=tma_centroids, raster_data=data) - set_transformation(data, transformation=transformations, set_all=True) - labels[f"core_{McmicroKeys.IMAGES_DIR_TMA}_{mask_stem}"] = data - tables_dict = _get_tables(path, markers, tma) +def _detect_pipeline(path: Path) -> McmicroPipeline: + """Infer which mcmicro pipeline produced ``path``.""" + if (path / McmicroKeys.PARAMS_FILE).exists(): + return McmicroPipeline.LABSYSPHARM + if (path / McmicroKeys.NFCORE_PIPELINE_INFO).exists() or (path / McmicroKeys.NFCORE_IMAGES_DIR_WSI).exists(): + return McmicroPipeline.NFCORE + raise ValueError( + f"Could not detect the mcmicro pipeline for {path}. Pass `pipeline='labsyspharm'` or " + f"`pipeline='nfcore'` explicitly." + ) + + +def _is_tma(path: Path, pipeline: McmicroPipeline) -> bool: + if pipeline == McmicroPipeline.LABSYSPHARM: + return bool(_load_params(path)["workflow"]["tma"]) + params = _load_nfcore_params(path) + if params is not None and McmicroKeys.NFCORE_IMAGES_DIR_TMA.value in params: + return bool(params[McmicroKeys.NFCORE_IMAGES_DIR_TMA.value]) + return (path / McmicroKeys.NFCORE_IMAGES_DIR_TMA).exists() + + +def _registration_dir(path: Path, pipeline: McmicroPipeline) -> Path: + if pipeline == McmicroPipeline.LABSYSPHARM: + return path / McmicroKeys.IMAGES_DIR_WSI + return path / McmicroKeys.NFCORE_IMAGES_DIR_WSI + + +def _dearray_dir(path: Path, pipeline: McmicroPipeline) -> Path: + if pipeline == McmicroPipeline.LABSYSPHARM: + return path / McmicroKeys.IMAGES_DIR_TMA + return path / McmicroKeys.NFCORE_IMAGES_DIR_TMA - sdata = SpatialData(images=images, labels=labels, tables=tables_dict) - return _set_reader_metadata(sdata, "mcmicro") + +def _centroids_file(path: Path, pipeline: McmicroPipeline) -> Path: + if pipeline == McmicroPipeline.LABSYSPHARM: + return path / McmicroKeys.COREOGRAPH_CENTROIDS + return path / McmicroKeys.NFCORE_COREOGRAPH_CENTROIDS + + +def _image_stem(sample: Path) -> str: + """Return the file stem with a (possibly double) suffix like ``.ome.tif`` stripped.""" + return sample.with_name(sample.stem).with_suffix("").stem + + +def _match_sample(name: str, sample_ids: list[str]) -> str: + """Return the sample id that best matches ``name`` (longest one that is a prefix/substring).""" + matches = [s for s in sample_ids if name.startswith(s) or s in name] + if matches: + return max(matches, key=len) + return name + + +def _match_core(name: str, core_ids: list[str]) -> str: + """Return the core id whose numeric token appears in ``name``. + + Core mask/table filenames embed the core number as a token, but also contain other digits + (e.g. ``exemplar-002_1_mask``), so match against the known core ids rather than the first digit. + """ + tokens = re.split(r"[_\-.]", name) + for t in tokens: + if t in core_ids: + return t + match = re.search(r"\d+", name) + return match.group() if match else name def _load_params(path: Path) -> Any: @@ -201,70 +264,327 @@ def _load_params(path: Path) -> Any: return params +def _load_nfcore_params(path: Path) -> dict[str, Any] | None: + """Load the most recent ``pipeline_info/params*.json`` if present.""" + params_files = sorted((path / McmicroKeys.NFCORE_PIPELINE_INFO).glob(McmicroKeys.NFCORE_PARAMS_GLOB.value)) + if not params_files: + return None + with open(params_files[-1]) as fp: + return json.load(fp) + + +def _read_marker_sheet(sheet_path: Path) -> pd.DataFrame | None: + """Read a wide marker sheet, or return ``None`` if the file is not one. + + Some nf-core outputs contain a long-format MultiQC validation report named ``*markersheet*`` + that is not a usable marker table; such files (lacking a ``marker_name`` column) are skipped. + """ + sep = "\t" if sheet_path.suffix == ".tsv" else "," + markers = pd.read_csv(sheet_path, sep=sep) + if McmicroKeys.MARKER_NAME.value not in markers.columns: + return None + if McmicroKeys.CHANNEL_NUMBER.value not in markers.columns: + markers[McmicroKeys.CHANNEL_NUMBER.value] = range(1, len(markers) + 1) + markers.index = markers[McmicroKeys.MARKER_NAME.value] + return markers + + +def _load_markers( + path: Path, + pipeline: McmicroPipeline, + markers_file: str | Path | None, + fallback_image: Path, + imread_kwargs: Mapping[str, Any], +) -> tuple[pd.DataFrame | None, list[str]]: + """Load the marker table and the channel names to attach to the registration images. + + Returns a ``(markers_df, channel_names)`` tuple. ``markers_df`` (used to build the labsyspharm + tables) may be ``None`` for nf-core, where table variables are derived from the quantification + CSV columns instead. ``channel_names`` are attached to the images when their length matches the + number of image channels. + """ + if pipeline == McmicroPipeline.LABSYSPHARM: + markers_path = Path(markers_file) if markers_file is not None else path / McmicroKeys.MARKERS_FILE + if not markers_path.exists(): + # historically required; surface a clear error rather than silently falling back. + raise FileNotFoundError(f"Expected a markers file at {markers_path}") + markers = _read_marker_sheet(markers_path) + if markers is None: + raise ValueError(f"Markers file {markers_path} has no `{McmicroKeys.MARKER_NAME.value}` column") + assert markers[McmicroKeys.CHANNEL_NUMBER.value].is_monotonic_increasing + return markers, markers[McmicroKeys.MARKER_NAME.value].tolist() + + # nf-core: the registration image (full ashlar mosaic) and the quantification tables can have + # different channel/marker sets (background subtraction drops channels). Collect every candidate + # marker sheet, then pick the channel names for the image by matching the image's channel count. + n_channels = int(imread(fallback_image, **imread_kwargs).shape[0]) + sheet_paths: list[Path] = [] + if markers_file is not None: + sheet_paths.append(Path(markers_file)) + else: + sheet_paths.append(path / McmicroKeys.MARKERS_FILE) + sheet_paths.extend( + sorted((path / McmicroKeys.NFCORE_BACKSUB_DIR).glob(McmicroKeys.NFCORE_BACKSUB_MARKERS_GLOB.value)) + ) + for glob_pattern in McmicroKeys.NFCORE_MARKERSHEET_GLOBS.value.split(";"): + sheet_paths.extend(sorted(path.glob(glob_pattern))) + sheets = [s for s in (_read_marker_sheet(p) for p in sheet_paths if p.exists()) if s is not None] + + # candidate name lists for the image, in priority order + candidate_names = [_safe_parse_channels(fallback_image)] + candidate_names += [s[McmicroKeys.MARKER_NAME.value].tolist() for s in sheets] + image_names = next((n for n in candidate_names if len(n) == n_channels), None) + if image_names is None: + image_names = next((n for n in candidate_names if n), None) + if image_names is None: + image_names = [f"channel_{i}" for i in range(n_channels)] + warnings.warn( + f"No markers file found; using integer channel names for {fallback_image.name}.", + UserWarning, + stacklevel=2, + ) + + # keep a marker sheet (any is fine: table variables are matched to it by name) for table var + markers_df = next((s for s in sheets if len(s) == n_channels), sheets[0] if sheets else None) + return markers_df, image_names + + +def _safe_parse_channels(image: Path) -> list[str]: + try: + return parse_channels(image) + except Exception: + return [] + + +def _labsyspharm_module(parent_stem: str, ref_id: str | None) -> str: + """Extract the segmentation module from a `-` directory name.""" + if ref_id is not None and parent_stem.endswith(f"-{ref_id}"): + return parent_stem[: -(len(ref_id) + 1)] + return parent_stem + + +def _get_labels( + path: Path, + pipeline: McmicroPipeline, + tma: bool, + tma_centroids: pd.DataFrame | None, + sample_ids: list[str], + core_ids: list[str], + wsi_image_id: str | None, + samples_masks: list[Path], + imread_kwargs: Mapping[str, Any], + labels_models_kwargs: Mapping[str, Any], +) -> dict[str, Any]: + labels: dict[str, Any] = {} + + if pipeline == McmicroPipeline.LABSYSPHARM: + samples_labels = list((path / McmicroKeys.LABELS_DIR).glob("*/*" + McmicroKeys.IMAGE_SUFFIX)) + for labels_path in samples_labels: + segmentation_stem = _image_stem(labels_path) + if not tma: + # the segmentation subdir is `-`; include the module so that + # multiple segmenters (e.g. unmicst + ilastik) don't collide on the same key. + module = _labsyspharm_module(labels_path.parent.stem, wsi_image_id) + data = imread(labels_path, **imread_kwargs).squeeze() + data = Labels2DModel.parse(data, transformations=_get_transformation(), **labels_models_kwargs) + labels[f"{wsi_image_id}_{module}_{segmentation_stem}"] = data + else: + core_id_search = re.search(r"\d+$", labels_path.parent.stem) + if core_id_search is None: + raise ValueError(f"Cannot infer core_id from {labels_path.parent}") + core_id = core_id_search.group() + module = _labsyspharm_module(labels_path.parent.stem, core_id) + data = imread(labels_path, **imread_kwargs).squeeze() + data = Labels2DModel.parse(data, **labels_models_kwargs) + transformations = _get_transformation(tma=int(core_id), tma_centroids=tma_centroids, raster_data=data) + set_transformation(data, transformation=transformations, set_all=True) + labels[f"core_{core_id}_{module}_{segmentation_stem}"] = data + else: # nf-core: segmentation// + for segmenter_dir in sorted(p for p in (path / McmicroKeys.LABELS_DIR).glob("*") if p.is_dir()): + segmenter = segmenter_dir.stem + mask_files = sorted(set(segmenter_dir.glob("*.tif")) | set(segmenter_dir.glob("*.tiff"))) + for mask_path in mask_files: + data = imread(mask_path, **imread_kwargs).squeeze() + if not tma: + sample_id = _match_sample(_image_stem(mask_path), sample_ids) + data = Labels2DModel.parse(data, transformations=_get_transformation(), **labels_models_kwargs) + labels[f"{sample_id}_{segmenter}"] = data + else: + core_id = _match_core(_image_stem(mask_path), core_ids) + data = Labels2DModel.parse(data, **labels_models_kwargs) + transformations = _get_transformation( + tma=int(core_id), tma_centroids=tma_centroids, raster_data=data + ) + set_transformation(data, transformation=transformations, set_all=True) + labels[f"core_{core_id}_{segmenter}"] = data + + # per-core dearray masks (shared by both layouts in TMA mode) + if tma: + for mask_path in samples_masks: + mask_stem = mask_path.stem + core_id = mask_stem.split("_")[0] + data = imread(mask_path, **imread_kwargs).squeeze() + data = Labels2DModel.parse(data, **labels_models_kwargs) + transformations = _get_transformation(tma=int(core_id), tma_centroids=tma_centroids, raster_data=data) + set_transformation(data, transformation=transformations, set_all=True) + labels[f"core_dearray_{mask_stem}"] = data + + return labels + + def _get_tables( path: Path, - marker_df: pd.DataFrame, + pipeline: McmicroPipeline, + marker_df: pd.DataFrame | None, tma: bool, + sample_ids: list[str], + core_ids: list[str], + label_keys: list[str], ) -> dict[str, AnnData]: - var = marker_df.marker_name.tolist() coords = [McmicroKeys.COORDS_X.value, McmicroKeys.COORDS_Y.value] - - table_paths = list((path / McmicroKeys.QUANTIFICATION_DIR).glob("*.csv")) - regions = [] - adatas = None - tables_dict = {} - for table_path in table_paths: - if not tma: - table_name = table_path.stem - adata, region = _create_anndata(csv_path=table_path, markers=marker_df, var=var, coords=coords, tma=tma) - table = TableModel.parse( - adata, region=region, region_key="region", instance_key=McmicroKeys.INSTANCE_KEY.value + instance_key = McmicroKeys.INSTANCE_KEY.value + tables_dict: dict[str, AnnData] = {} + + if pipeline == McmicroPipeline.LABSYSPHARM: + # WSI: one table per CSV; TMA: all cores concatenated into a single table (legacy behavior). + regions: list[str] = [] + adatas = None + for table_path in (path / McmicroKeys.QUANTIFICATION_DIR).glob("*.csv"): + name, region, var, var_df = _table_spec( + table_path, None, pipeline, tma, marker_df, coords, sample_ids, core_ids, label_keys ) - tables_dict[table_name] = table - else: - adata, region = _create_anndata(csv_path=table_path, markers=marker_df, var=var, coords=coords, tma=tma) - regions.append(region) - # TODO: check validity of output - if not adatas: - adatas = adata + adata = _build_anndata(table_path, var_df, var, coords, region) + if not tma: + tables_dict[name] = TableModel.parse( + adata, region=region, region_key="region", instance_key=instance_key + ) else: - adatas = ad.concat([adatas, adata], index_unique="_") - table = TableModel.parse( - adatas, region=regions, region_key="region", instance_key=McmicroKeys.INSTANCE_KEY.value + regions.append(region) + adatas = adata if adatas is None else ad.concat([adatas, adata], index_unique="_") + tables_dict["segmentation_table"] = TableModel.parse( + adatas, region=regions, region_key="region", instance_key=instance_key + ) + return tables_dict + + # nf-core: quantification/mcquant//.csv + entries = [(p, p.parent.stem) for p in (path / McmicroKeys.NFCORE_QUANTIFICATION_DIR).glob("*/*.csv")] + if not tma: + for table_path, segmenter in entries: + name, region, var, var_df = _table_spec( + table_path, segmenter, pipeline, tma, marker_df, coords, sample_ids, core_ids, label_keys ) - tables_dict["segmentation_table"] = table - + tables_dict[name] = TableModel.parse( + _build_anndata(table_path, var_df, var, coords, region), + region=region, + region_key="region", + instance_key=instance_key, + ) + return tables_dict + + # nf-core TMA: one table per segmenter, concatenating that segmenter's cores. + by_segmenter: dict[str, list[Path]] = {} + for table_path, segmenter in entries: + by_segmenter.setdefault(segmenter, []).append(table_path) + for segmenter, csvs in by_segmenter.items(): + regions = [] + adatas = None + for table_path in sorted(csvs): + _, region, var, var_df = _table_spec( + table_path, segmenter, pipeline, tma, marker_df, coords, sample_ids, core_ids, label_keys + ) + adata = _build_anndata(table_path, var_df, var, coords, region) + regions.append(region) + adatas = adata if adatas is None else ad.concat([adatas, adata], index_unique="_") + tables_dict[f"{segmenter}_table"] = TableModel.parse( + adatas, region=regions, region_key="region", instance_key=instance_key + ) return tables_dict -def _create_anndata( +def _table_spec( + csv_path: Path, + segmenter: str | None, + pipeline: McmicroPipeline, + tma: bool, + marker_df: pd.DataFrame | None, + coords: list[str], + sample_ids: list[str], + core_ids: list[str], + label_keys: list[str], +) -> tuple[str, str, list[str], pd.DataFrame]: + """Return ``(table_name, region_value, var, var_df)`` for a quantification CSV.""" + if pipeline == McmicroPipeline.LABSYSPHARM: + assert marker_df is not None + var = marker_df[McmicroKeys.MARKER_NAME.value].tolist() + # filename form: --_, e.g. exemplar-001--unmicst_cell + labels_basename = csv_path.stem.split("_")[-1] + sample_id_search = re.search(r"^(.*?)--", csv_path.stem) + sample_id = sample_id_search.groups()[0] if sample_id_search else None + if not sample_id: + raise ValueError( + f"Csv filename should be in form --_, got {csv_path.stem} " + ) + rest = csv_path.stem[len(sample_id) + 2 :] # drop "--" + module = rest[: -(len(labels_basename) + 1)] if rest.endswith("_" + labels_basename) else rest + prefix = f"core_{sample_id}" if tma else sample_id + region_value = f"{prefix}_{module}_{labels_basename}" + return csv_path.stem, region_value, var, marker_df + + # nf-core: variables are the CSV columns between CellID and the coordinate columns. + header = pd.read_csv(csv_path, nrows=0).columns.tolist() + var = _nfcore_marker_columns(header, coords) + var_df = ( + marker_df.reindex(var) + if marker_df is not None + else pd.DataFrame(index=pd.Index(var, name=McmicroKeys.MARKER_NAME.value)) + ) + # in TMA mode the CSV corresponds to a core; in WSI mode it corresponds to a sample + unit = _match_core(csv_path.stem, core_ids) if tma else _match_sample(csv_path.stem, sample_ids) + # table names share a namespace with images/labels, so suffix to avoid colliding with the label + table_name = f"{unit}_{segmenter}_table" + region_value = _match_label_key(unit, segmenter, tma, label_keys) + return table_name, region_value, var, var_df + + +def _nfcore_marker_columns(header: list[str], coords: list[str]) -> list[str]: + """Marker columns are those after ``CellID`` and before the first coordinate column.""" + instance = McmicroKeys.INSTANCE_KEY.value + start = header.index(instance) + 1 if instance in header else 0 + end = min((header.index(c) for c in coords if c in header), default=len(header)) + return header[start:end] + + +def _match_label_key(sample_id: str, segmenter: str | None, tma: bool, label_keys: list[str]) -> str: + """Link a quantification table to the label element produced by the same segmenter. + + Segmentation and quantification directories do not always share the same segmenter name + (e.g. ``deepcell_mesmer`` vs ``mesmer``), so we match by prefix and then by substring. + """ + prefix = f"core_{sample_id}" if tma else sample_id + exact = f"{prefix}_{segmenter}" + if exact in label_keys: + return exact + candidates = [k for k in label_keys if k.startswith(f"{prefix}_")] + for k in candidates: + seg = k[len(prefix) + 1 :] + if segmenter and (segmenter in seg or seg in segmenter): + return k + return exact + + +def _build_anndata( csv_path: Path, markers: pd.DataFrame, - var: list[str], # mypy has a bug with ellips, should be list[str, ...] + var: list[str], coords: list[str], - tma: bool, -) -> tuple[AnnData, str]: - labels_basename = csv_path.stem.split("_")[-1] - pattern = r"^(.*?)--" - sample_id_search = re.search(pattern, csv_path.stem) - sample_id = sample_id_search.groups()[0] if sample_id_search else None - if not sample_id: - raise ValueError( - f"Csv filename should be in form --_, got {csv_path.stem} " - ) + region_value: str, +) -> AnnData: table = pd.read_csv(csv_path) - - if not tma: - region_value = sample_id + "_" + labels_basename - else: - region_value = "core_" + sample_id + "_" + labels_basename - table[McmicroKeys.INSTANCE_KEY] = table[McmicroKeys.INSTANCE_KEY] adata = AnnData( - table[var].to_numpy(), + table[var].to_numpy().astype(np.float32), obs=table.drop(columns=var + coords), var=markers, obsm={"spatial": table[coords].to_numpy()}, - dtype=float, ) adata.obs["region"] = pd.Categorical([region_value] * len(adata)) - return adata, region_value + return adata diff --git a/tests/test_mcmicro.py b/tests/test_mcmicro.py new file mode 100644 index 00000000..2f1358c8 --- /dev/null +++ b/tests/test_mcmicro.py @@ -0,0 +1,321 @@ +"""Tests for the mcmicro reader across the labsyspharm and nf-core pipeline layouts. + +The reader accepts two different mcmicro output directory layouts (the original +``labsyspharm/mcmicro`` Nextflow pipeline and the newer ``nf-core/mcmicro`` pipeline), in +either whole-slide-image (WSI) or tissue-microarray (TMA) mode. These tests build minimal +synthetic output trees on disk (following the ``test_macsima.py`` pattern) and assert that the +reader auto-detects the pipeline and produces the expected image/label/table keys. +""" + +from pathlib import Path +from tempfile import TemporaryDirectory + +import numpy as np +import pandas as pd +import pytest +import tifffile +import yaml +from click.testing import CliRunner +from spatialdata import get_extent, read_zarr + +from spatialdata_io.__main__ import mcmicro_wrapper +from spatialdata_io._constants._constants import McmicroPipeline +from spatialdata_io.readers.mcmicro import _detect_pipeline, mcmicro + +N_CHANNELS = 2 +MARKER_NAMES = ["DAPI", "CD3"] + + +def _write_image(path: Path, n_channels: int = N_CHANNELS, size: int = 8) -> None: + path.parent.mkdir(parents=True, exist_ok=True) + rng = np.random.default_rng(0) + data = rng.integers(0, 255, size=(n_channels, size, size), dtype=np.uint16) + tifffile.imwrite(path, data) + + +def _write_mask(path: Path, size: int = 8) -> None: + path.parent.mkdir(parents=True, exist_ok=True) + data = np.zeros((size, size), dtype=np.uint32) + data[:4, :4] = 1 + data[4:, 4:] = 2 + tifffile.imwrite(path, data) + + +def _write_markers(path: Path) -> None: + path.parent.mkdir(parents=True, exist_ok=True) + pd.DataFrame({"channel_number": range(1, N_CHANNELS + 1), "marker_name": MARKER_NAMES}).to_csv(path, index=False) + + +def _write_quant(path: Path, n_cells: int = 2) -> None: + path.parent.mkdir(parents=True, exist_ok=True) + df = pd.DataFrame( + { + "CellID": range(1, n_cells + 1), + **{name: np.arange(n_cells, dtype=float) for name in MARKER_NAMES}, + "X_centroid": np.arange(n_cells, dtype=float), + "Y_centroid": np.arange(n_cells, dtype=float), + "Area": np.arange(n_cells, dtype=float), + } + ) + df.to_csv(path, index=False) + + +# --------------------------------------------------------------------------------------------- # +# tree builders +# --------------------------------------------------------------------------------------------- # +def _build_labsyspharm_wsi(root: Path, sample: str = "test") -> None: + (root / "qc").mkdir(parents=True, exist_ok=True) + with open(root / "qc" / "params.yml", "w") as fp: + yaml.safe_dump({"workflow": {"tma": False}}, fp) + _write_markers(root / "markers.csv") + _write_image(root / "registration" / f"{sample}.ome.tif") + _write_mask(root / "segmentation" / f"unmicst-{sample}" / "cell.ome.tif") + _write_quant(root / "quantification" / f"{sample}--unmicst_cell.csv") + + +def _build_nfcore_wsi(root: Path, samples: list[str], markers: bool = True) -> None: + (root / "pipeline_info").mkdir(parents=True, exist_ok=True) + if markers: + _write_markers(root / "markers.csv") + for sample in samples: + _write_image(root / "registration" / "ashlar" / f"{sample}.ome.tif") + _write_mask(root / "segmentation" / "cellpose" / f"{sample}_mask.ome.tif") + _write_quant(root / "quantification" / "mcquant" / "cellpose" / f"{sample}.csv") + + +def _write_centroids(path: Path, n_cores: int) -> None: + path.parent.mkdir(parents=True, exist_ok=True) + lines = [f"{10 * (i + 1)} {20 * (i + 1)}" for i in range(n_cores)] + path.write_text("\n".join(lines) + "\n") + + +def _build_labsyspharm_tma(root: Path, cores: list[int]) -> None: + (root / "qc").mkdir(parents=True, exist_ok=True) + with open(root / "qc" / "params.yml", "w") as fp: + yaml.safe_dump({"workflow": {"tma": True}}, fp) + _write_markers(root / "markers.csv") + _write_image(root / "registration" / "tma.ome.tif") + _write_centroids(root / "qc" / "coreograph" / "centroidsY-X.txt", len(cores)) + for core in cores: + _write_image(root / "dearray" / f"{core}.ome.tif") + _write_mask(root / "dearray" / "masks" / f"{core}_mask.tif") + _write_mask(root / "segmentation" / f"unmicst-{core}" / "cell.ome.tif") + _write_quant(root / "quantification" / f"{core}--unmicst_cell.csv") + + +def _build_nfcore_tma(root: Path, cores: list[int], sample: str = "exemplar-002") -> None: + """Mirror the real nf-core TMA layout. + + Numbered cores + a TMA_MAP, two segmenters with tool-specific mask names (cellpose + ``{core}_cp_masks.tif``, mesmer ``{sample}_{core}_mask.tif``), and per-core-per-segmenter + quantification with a prelude markersheet. + """ + (root / "pipeline_info").mkdir(parents=True, exist_ok=True) + (root / "prelude").mkdir(parents=True, exist_ok=True) + pd.DataFrame({"channel_number": range(1, N_CHANNELS + 1), "marker_name": MARKER_NAMES}).to_csv( + root / "prelude" / "markers_markersheet_mqc.tsv", sep="\t", index=False + ) + _write_image(root / "registration" / "ashlar" / f"{sample}.ome.tif") + _write_centroids(root / "tma_dearray" / "centroidsY-X.txt", len(cores)) + _write_mask(root / "tma_dearray" / "TMA_MAP.tif") # must be ignored as a core image + for core in cores: + _write_image(root / "tma_dearray" / f"{core}.tif") + _write_mask(root / "tma_dearray" / "masks" / f"{core}_mask.tif") + _write_mask(root / "segmentation" / "cellpose" / f"{core}_cp_masks.tif") + _write_mask(root / "segmentation" / "deepcell_mesmer" / f"{sample}_{core}_mask.tif") + _write_quant(root / "quantification" / "mcquant" / "cellpose" / f"{core}_{core}_cp_masks.csv") + _write_quant(root / "quantification" / "mcquant" / "mesmer" / f"{core}_{sample}_{core}_mask.csv") + + +# --------------------------------------------------------------------------------------------- # +# detection +# --------------------------------------------------------------------------------------------- # +def test_detect_labsyspharm(tmp_path: Path) -> None: + _build_labsyspharm_wsi(tmp_path) + assert _detect_pipeline(tmp_path) == McmicroPipeline.LABSYSPHARM + + +def test_detect_nfcore(tmp_path: Path) -> None: + _build_nfcore_wsi(tmp_path, ["s1"]) + assert _detect_pipeline(tmp_path) == McmicroPipeline.NFCORE + + +def test_detect_raises_when_ambiguous(tmp_path: Path) -> None: + with pytest.raises(ValueError, match="Could not detect"): + _detect_pipeline(tmp_path) + + +# --------------------------------------------------------------------------------------------- # +# WSI +# --------------------------------------------------------------------------------------------- # +def test_labsyspharm_wsi(tmp_path: Path) -> None: + _build_labsyspharm_wsi(tmp_path, "test") + sdata = mcmicro(tmp_path) + assert set(sdata.images) == {"test_image"} + assert "test_unmicst_cell" in sdata.labels + assert "test--unmicst_cell" in sdata.tables + assert sdata.tables["test--unmicst_cell"].obs["region"].cat.categories.tolist() == ["test_unmicst_cell"] + assert list(sdata.images["test_image"].c.values) == MARKER_NAMES + + +def test_labsyspharm_wsi_multi_segmenter(tmp_path: Path) -> None: + """Two segmenters writing the same `cell.ome.tif` must not collide on the label key.""" + (tmp_path / "qc").mkdir(parents=True) + with open(tmp_path / "qc" / "params.yml", "w") as fp: + yaml.safe_dump({"workflow": {"tma": False}}, fp) + _write_markers(tmp_path / "markers.csv") + _write_image(tmp_path / "registration" / "test.ome.tif") + for seg in ["unmicst", "ilastik"]: + _write_mask(tmp_path / "segmentation" / f"{seg}-test" / "cell.ome.tif") + _write_quant(tmp_path / "quantification" / f"test--{seg}_cell.csv") + sdata = mcmicro(tmp_path) + assert {"test_unmicst_cell", "test_ilastik_cell"} <= set(sdata.labels) + assert {"test--unmicst_cell", "test--ilastik_cell"} == set(sdata.tables) + assert sdata.tables["test--ilastik_cell"].obs["region"].cat.categories.tolist() == ["test_ilastik_cell"] + + +def test_nfcore_wsi_single_sample(tmp_path: Path) -> None: + _build_nfcore_wsi(tmp_path, ["s1"]) + sdata = mcmicro(tmp_path) + assert set(sdata.images) == {"s1_image"} + assert "s1_cellpose" in sdata.labels + assert "s1_cellpose_table" in sdata.tables + assert sdata.tables["s1_cellpose_table"].obs["region"].cat.categories.tolist() == ["s1_cellpose"] + # WSI uses an identity transform: image and its label share the same global extent, anchored at 0 + img_ext = get_extent(sdata["s1_image"], coordinate_system="global") + lab_ext = get_extent(sdata["s1_cellpose"], coordinate_system="global") + assert img_ext["x"][0] == 0 and img_ext["y"][0] == 0 + assert img_ext == lab_ext + + +def test_nfcore_wsi_multi_sample(tmp_path: Path) -> None: + _build_nfcore_wsi(tmp_path, ["s1", "s2", "s3"]) + sdata = mcmicro(tmp_path) + assert set(sdata.images) == {"s1_image", "s2_image", "s3_image"} + assert {"s1_cellpose", "s2_cellpose", "s3_cellpose"} <= set(sdata.labels) + assert {"s1_cellpose_table", "s2_cellpose_table", "s3_cellpose_table"} == set(sdata.tables) + + +def test_nfcore_wsi_multi_segmenter(tmp_path: Path) -> None: + """Mirror the real S3 multiseg layout: tool-specific mask names and mismatched dir names.""" + sample = "exemplar-001-10cycles" + (tmp_path / "pipeline_info").mkdir(parents=True) + _write_markers(tmp_path / "markers.csv") + _write_image(tmp_path / "registration" / "ashlar" / f"{sample}.ome.tif") + # segmentation dirs and mask filenames vary by tool + _write_mask(tmp_path / "segmentation" / "cellpose" / f"{sample}_backsub.ome_cp_masks.tif") + _write_mask(tmp_path / "segmentation" / "mccellpose" / f"{sample}_mask.ome.tif") + _write_mask(tmp_path / "segmentation" / "deepcell_mesmer" / f"{sample}_mask.tif") + # quantification dir names differ from segmentation dirs (mesmer vs deepcell_mesmer) + _write_quant(tmp_path / "quantification" / "mcquant" / "cellpose" / f"{sample}_backsub_{sample}_backsub.csv") + _write_quant(tmp_path / "quantification" / "mcquant" / "mccellpose" / f"{sample}.csv") + _write_quant(tmp_path / "quantification" / "mcquant" / "mesmer" / f"{sample}.csv") + + sdata = mcmicro(tmp_path) + assert set(sdata.images) == {f"{sample}_image"} + assert {f"{sample}_cellpose", f"{sample}_mccellpose", f"{sample}_deepcell_mesmer"} <= set(sdata.labels) + assert {f"{sample}_cellpose_table", f"{sample}_mccellpose_table", f"{sample}_mesmer_table"} == set(sdata.tables) + # the mesmer table must link (by substring) to the deepcell_mesmer label element + mesmer_region = sdata.tables[f"{sample}_mesmer_table"].obs["region"].cat.categories.tolist() + assert mesmer_region == [f"{sample}_deepcell_mesmer"] + + +def test_nfcore_wsi_markers_fallback(tmp_path: Path) -> None: + _build_nfcore_wsi(tmp_path, ["s1"], markers=False) + # drop quantification so the derived channel names (which need not match the marker columns + # written by mcquant) don't collide with the synthetic CSV columns + import shutil + + shutil.rmtree(tmp_path / "quantification") + with pytest.warns(UserWarning, match="No markers file found"): + sdata = mcmicro(tmp_path) + # OME metadata carries no channel names for our synthetic tiffs -> integer channel names + assert len(sdata.images["s1_image"].c) == N_CHANNELS + assert list(sdata.images["s1_image"].c.values) == ["channel_0", "channel_1"] + + +def test_explicit_pipeline_override(tmp_path: Path) -> None: + # nf-core tree without pipeline_info/ still auto-detects via registration/ashlar, + # and an explicit override is honored regardless. + _build_nfcore_wsi(tmp_path, ["s1"]) + (tmp_path / "pipeline_info").rmdir() + assert _detect_pipeline(tmp_path) == McmicroPipeline.NFCORE + sdata = mcmicro(tmp_path, pipeline="nfcore") + assert set(sdata.images) == {"s1_image"} + + +# --------------------------------------------------------------------------------------------- # +# TMA +# --------------------------------------------------------------------------------------------- # +def test_labsyspharm_tma(tmp_path: Path) -> None: + _build_labsyspharm_tma(tmp_path, [1, 2]) + sdata = mcmicro(tmp_path) + assert "tma_map" in sdata.images + assert {"core_1_image", "core_2_image"} <= set(sdata.images) + assert {"core_1_unmicst_cell", "core_2_unmicst_cell"} <= set(sdata.labels) + assert "segmentation_table" in sdata.tables + + +def test_nfcore_tma(tmp_path: Path) -> None: + sample = "exemplar-002" + _build_nfcore_tma(tmp_path, [1, 2, 3, 4], sample=sample) + sdata = mcmicro(tmp_path) + # global mosaic + one image per core (TMA_MAP.tif must be excluded) + assert "tma_map" in sdata.images + assert {f"core_{c}_image" for c in (1, 2, 3, 4)} | {"tma_map"} == set(sdata.images) + # per-core, per-segmenter labels (core id parsed correctly despite "exemplar-002" in the name) + assert {f"core_{c}_cellpose" for c in (1, 2, 3, 4)} <= set(sdata.labels) + assert {f"core_{c}_deepcell_mesmer" for c in (1, 2, 3, 4)} <= set(sdata.labels) + # one table per segmenter, each concatenating its 4 cores + assert {"cellpose_table", "mesmer_table"} == set(sdata.tables) + mesmer_regions = sorted(sdata.tables["mesmer_table"].obs["region"].cat.categories.tolist()) + assert mesmer_regions == [f"core_{c}_deepcell_mesmer" for c in (1, 2, 3, 4)] + # each core is translated to its centroid, so cores occupy distinct, non-origin global extents + ext = {c: get_extent(sdata[f"core_{c}_image"], coordinate_system="global") for c in (1, 2, 3, 4)} + origins = {(ext[c]["x"][0], ext[c]["y"][0]) for c in (1, 2, 3, 4)} + assert len(origins) == 4 # no two cores overlap at the same origin + assert all(ext[c]["x"][0] > 0 and ext[c]["y"][0] > 0 for c in (1, 2, 3, 4)) + # a core image and its label share the same translated extent (aligned) + assert get_extent(sdata["core_1_cellpose"], coordinate_system="global") == ext[1] + + +# --------------------------------------------------------------------------------------------- # +# real public dataset (downloaded by .github/workflows/prepare_test_data.yaml) +# --------------------------------------------------------------------------------------------- # +# Real nf-core/mcmicro 2.0.0 `test` profile output. See prepare_test_data.yaml for provenance. +_REAL_NFCORE = Path("./data/mcmicro-nfcore-wsi-test") + + +@pytest.mark.skipif(not _REAL_NFCORE.is_dir(), reason="mcmicro-nfcore-wsi-test data not downloaded") +def test_nfcore_real_data() -> None: + # the `test` profile output has no usable marker sheet, so channels fall back to integer names + with pytest.warns(UserWarning, match="No markers file found"): + sdata = mcmicro(_REAL_NFCORE) + assert _detect_pipeline(_REAL_NFCORE) == McmicroPipeline.NFCORE + assert "TEST1_image" in sdata.images + assert "TEST1_mccellpose" in sdata.labels + table = sdata.tables["TEST1_mccellpose_table"] + assert table.obs["region"].cat.categories.tolist() == ["TEST1_mccellpose"] + assert table.n_obs > 0 + + +@pytest.mark.skipif(not _REAL_NFCORE.is_dir(), reason="mcmicro-nfcore-wsi-test data not downloaded") +def test_cli_nfcore_real_data(runner: CliRunner) -> None: + with TemporaryDirectory() as outdir: + output_zarr = Path(outdir) / "data.zarr" + result = runner.invoke(mcmicro_wrapper, ["--input", str(_REAL_NFCORE), "--output", str(output_zarr)]) + assert result.exit_code == 0, result.output + _ = read_zarr(output_zarr) + + +# --------------------------------------------------------------------------------------------- # +# CLI round-trip +# --------------------------------------------------------------------------------------------- # +@pytest.mark.parametrize("builder", [_build_labsyspharm_wsi, lambda root: _build_nfcore_wsi(root, ["s1"])]) +def test_cli_mcmicro(runner: CliRunner, tmp_path: Path, builder) -> None: # type: ignore[no-untyped-def] + builder(tmp_path) + with TemporaryDirectory() as outdir: + output_zarr = Path(outdir) / "data.zarr" + result = runner.invoke(mcmicro_wrapper, ["--input", str(tmp_path), "--output", str(output_zarr)]) + assert result.exit_code == 0, result.output + _ = read_zarr(output_zarr)