From c7b287ac84562544706c67dd372dbb9b1102d386 Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 14:10:28 +0200 Subject: [PATCH 01/11] feat(sanitize): central non-contiguous label-ID sanitation Relabel arbitrary or non-contiguous mask label IDs (e.g. {1,5,17}) to contiguous 1..N as early as possible in the pipeline, so every measurement function downstream can assume clean labels and focus on the math. The caller's array is never mutated and results are reported against the original IDs. - new cp_measure._sanitize: sanitize_masks() built on skimage.segmentation.relabel_sequential (already a dependency), plus a thin @sanitize_labels decorator that finds the label argument by name - featurizer sanitizes once at the top of the per-mask loop and emits original IDs as row labels (was range(1, n+1)) - @sanitize_labels applied to every core/correlation/numba get_* entry so direct callers get the same guarantee (cheap idempotent guard) - raises ValueError on negative / non-integer labels - tests: gapped {1,5,17} == contiguous baseline per function, featurizer end-to-end with original-ID rows, no-mutation, validation, 3D, registry decoration meta-test (partial-aware, multimask excluded) Co-Authored-By: Claude Opus 4.8 (1M context) --- src/cp_measure/_sanitize.py | 107 ++++++++++ src/cp_measure/core/measurecolocalization.py | 6 + src/cp_measure/core/measuregranularity.py | 2 + src/cp_measure/core/measureobjectintensity.py | 3 + .../measureobjectintensitydistribution.py | 3 + src/cp_measure/core/measureobjectsizeshape.py | 4 + src/cp_measure/core/measuretexture.py | 3 + .../core/numba/measureobjectintensity.py | 3 + src/cp_measure/featurizer.py | 27 ++- test/test_sanitize.py | 199 ++++++++++++++++++ 10 files changed, 346 insertions(+), 11 deletions(-) create mode 100644 src/cp_measure/_sanitize.py create mode 100644 test/test_sanitize.py diff --git a/src/cp_measure/_sanitize.py b/src/cp_measure/_sanitize.py new file mode 100644 index 0000000..4e8ba7f --- /dev/null +++ b/src/cp_measure/_sanitize.py @@ -0,0 +1,107 @@ +"""Central input sanitation for non-contiguous mask label IDs. + +cp_measure's measurement functions assume object labels are the contiguous +integers ``1..N`` (the convention documented on +:func:`cp_measure.featurizer.featurize`). Real segmentations do not always +honour that: labels may have gaps (``{1, 5, 17}``) or arbitrary values. This +module provides the single policy that maps arbitrary positive-integer labels to +``1..N`` *before* any math runs, so every downstream function can assume clean +labels and focus on the calculation. + +Two pieces: + +* :func:`sanitize_masks` — the policy. Returns ``(clean, ids)`` where ``clean`` + has labels ``1..N`` and ``ids[i]`` is the *original* label of rank ``i + 1``. + It never mutates the caller's array (the relabel path returns a fresh copy; the + already-clean fast path returns the input unchanged). +* :func:`sanitize_labels` — a thin decorator that applies the policy to whichever + argument of a ``get_*`` function holds the label image, so direct callers get + the same guarantee as the featurizer. + +Relabelling is cheap relative to a single feature (<1 % of a featurized image), +so the featurizer sanitizes once up front and the per-function decorator is a +cheap idempotent guard (a single :func:`numpy.unique`) for direct callers. +""" + +import functools +import inspect +from typing import Callable + +import numpy +from numpy.typing import NDArray +from skimage.segmentation import relabel_sequential + +# Argument names that, across the ``get_*`` functions, hold the label image. +_MASK_PARAMS = ("masks", "labels", "mask") + + +def _is_contiguous(masks: NDArray) -> bool: + """Are the positive labels exactly ``1..N``? + + Uses :func:`numpy.unique` rather than a dense ``bincount`` so the check is + safe and bounded for any dtype, negative values, or large label values. + """ + unique = numpy.unique(masks) + positive = unique[unique > 0] + return positive.size == 0 or (positive[0] == 1 and positive[-1] == positive.size) + + +def sanitize_masks(masks: NDArray) -> tuple[NDArray, NDArray[numpy.int64]]: + """Relabel arbitrary positive labels to contiguous ``1..N``. + + Parameters + ---------- + masks + Integer label array (any number of dimensions). Background is ``0``. + + Returns + ------- + clean + Array with labels ``1..N`` in ascending original-label order. The + already-contiguous input is returned unchanged (no copy); otherwise a + fresh relabelled copy is returned. The input is never mutated. + ids + ``ids[i]`` is the original label whose sanitized value is ``i + 1`` + (ascending). Use it to report results against the caller's IDs. + + Raises + ------ + ValueError + If ``masks`` is not an integer array or contains negative values. + """ + if not numpy.issubdtype(masks.dtype, numpy.integer): + raise ValueError(f"labels must be an integer array, got dtype {masks.dtype!r}") + if masks.size and masks.min() < 0: + raise ValueError("labels must be non-negative") + + if _is_contiguous(masks): + n = int(masks.max(initial=0)) + return masks, numpy.arange(1, n + 1, dtype=numpy.int64) + + clean, _forward, _inverse = relabel_sequential(masks) + ids = numpy.unique(masks) + ids = ids[ids > 0].astype(numpy.int64) + return clean, ids + + +def sanitize_labels(func: Callable) -> Callable: + """Decorate a ``get_*`` function to sanitize its label argument. + + The label argument is detected by name (:data:`_MASK_PARAMS`), so the + decorator is position-independent. Functions without a recognised label + argument (e.g. the two-mask ``multimask`` functions, which take + ``masks1``/``masks2``) are returned unchanged. + """ + sig = inspect.signature(func) + param = next((name for name in _MASK_PARAMS if name in sig.parameters), None) + if param is None: + return func + + @functools.wraps(func) + def wrapper(*args, **kwargs): + bound = sig.bind(*args, **kwargs) + bound.arguments[param], _ids = sanitize_masks(bound.arguments[param]) + return func(*bound.args, **bound.kwargs) + + wrapper._sanitized = True # type: ignore[attr-defined] + return wrapper diff --git a/src/cp_measure/core/measurecolocalization.py b/src/cp_measure/core/measurecolocalization.py index 9c97a5f..221e8d9 100644 --- a/src/cp_measure/core/measurecolocalization.py +++ b/src/cp_measure/core/measurecolocalization.py @@ -86,6 +86,8 @@ import scipy.stats from scipy.linalg import lstsq +from cp_measure._sanitize import sanitize_labels + M_IMAGES = "Across entire image" M_OBJECTS = "Within objects" M_IMAGES_AND_OBJECTS = "Both" @@ -544,6 +546,7 @@ def infer_scale(data: numpy.ndarray) -> int: # --------------------------------------------------------------------------- +@sanitize_labels def get_correlation_pearson( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], @@ -562,6 +565,7 @@ def get_correlation_pearson( return {F_CORRELATION_FORMAT: corrs, F_SLOPE_FORMAT: slopes} +@sanitize_labels def get_correlation_manders_fold( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], @@ -581,6 +585,7 @@ def get_correlation_manders_fold( return {f"{F_MANDERS_FORMAT}_1": m1_list, f"{F_MANDERS_FORMAT}_2": m2_list} +@sanitize_labels def get_correlation_rwc( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], @@ -600,6 +605,7 @@ def get_correlation_rwc( return {f"{F_RWC_FORMAT}_1": r1, f"{F_RWC_FORMAT}_2": r2} +@sanitize_labels def get_correlation_costes( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], diff --git a/src/cp_measure/core/measuregranularity.py b/src/cp_measure/core/measuregranularity.py index d0cffee..0d7fc7e 100644 --- a/src/cp_measure/core/measuregranularity.py +++ b/src/cp_measure/core/measuregranularity.py @@ -54,10 +54,12 @@ import numpy import scipy.ndimage import skimage.morphology +from cp_measure._sanitize import sanitize_labels from cp_measure.utils import _ensure_np_array as fix from numpy.typing import NDArray +@sanitize_labels def get_granularity( mask: NDArray[numpy.integer], pixels: NDArray[numpy.floating], diff --git a/src/cp_measure/core/measureobjectintensity.py b/src/cp_measure/core/measureobjectintensity.py index 3ca99dd..1b425a2 100644 --- a/src/cp_measure/core/measureobjectintensity.py +++ b/src/cp_measure/core/measureobjectintensity.py @@ -3,6 +3,8 @@ import skimage.segmentation from numpy.typing import NDArray +from cp_measure._sanitize import sanitize_labels + __doc__ = """ MeasureObjectIntensity ====================== @@ -136,6 +138,7 @@ def _interp(sorted_arr: numpy.ndarray, n: int, frac: float, legacy: bool) -> flo return float(sorted_arr[lo] * (1.0 - f) + sorted_arr[hi] * f) +@sanitize_labels def get_intensity( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating], diff --git a/src/cp_measure/core/measureobjectintensitydistribution.py b/src/cp_measure/core/measureobjectintensitydistribution.py index f142d2a..9450250 100644 --- a/src/cp_measure/core/measureobjectintensitydistribution.py +++ b/src/cp_measure/core/measureobjectintensitydistribution.py @@ -5,6 +5,7 @@ import numpy.ma import scipy.ndimage import scipy.sparse +from cp_measure._sanitize import sanitize_labels from cp_measure.utils import masks_to_ijv from numpy.typing import NDArray @@ -89,6 +90,7 @@ } +@sanitize_labels def get_radial_distribution( labels: NDArray[numpy.integer], pixels: NDArray[numpy.floating], @@ -305,6 +307,7 @@ def get_radial_distribution( return results +@sanitize_labels def get_radial_zernikes( labels: NDArray[numpy.integer], pixels: NDArray[numpy.floating], diff --git a/src/cp_measure/core/measureobjectsizeshape.py b/src/cp_measure/core/measureobjectsizeshape.py index 0e32b17..be1521b 100644 --- a/src/cp_measure/core/measureobjectsizeshape.py +++ b/src/cp_measure/core/measureobjectsizeshape.py @@ -7,6 +7,7 @@ import numpy import scipy.ndimage import skimage.measure +from cp_measure._sanitize import sanitize_labels from cp_measure.utils import masks_to_ijv __doc__ = """\ @@ -566,6 +567,7 @@ """ +@sanitize_labels def get_sizeshape( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating] | None, @@ -1003,6 +1005,7 @@ def get_sizeshape( return results +@sanitize_labels def get_zernike( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating] | None, @@ -1027,6 +1030,7 @@ def get_zernike( return results +@sanitize_labels def get_feret( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating] | None ) -> dict[str, NDArray[numpy.floating]]: diff --git a/src/cp_measure/core/measuretexture.py b/src/cp_measure/core/measuretexture.py index 659fe58..769a88e 100644 --- a/src/cp_measure/core/measuretexture.py +++ b/src/cp_measure/core/measuretexture.py @@ -149,11 +149,14 @@ from numpy.typing import NDArray import skimage.util +from cp_measure._sanitize import sanitize_labels + F_HARALICK = """AngularSecondMoment Contrast Correlation Variance InverseDifferenceMoment SumAverage SumVariance SumEntropy Entropy DifferenceVariance DifferenceEntropy InfoMeas1 InfoMeas2""".split() +@sanitize_labels def get_texture( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating], diff --git a/src/cp_measure/core/numba/measureobjectintensity.py b/src/cp_measure/core/numba/measureobjectintensity.py index 57cf841..aaff39c 100644 --- a/src/cp_measure/core/numba/measureobjectintensity.py +++ b/src/cp_measure/core/numba/measureobjectintensity.py @@ -14,6 +14,8 @@ import skimage.segmentation from numpy.typing import NDArray +from cp_measure._sanitize import sanitize_labels + from cp_measure.core.measureobjectintensity import ( C_LOCATION, INTEGRATED_INTENSITY, @@ -50,6 +52,7 @@ ) +@sanitize_labels def get_intensity( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating], diff --git a/src/cp_measure/featurizer.py b/src/cp_measure/featurizer.py index a1766a6..dccb3a6 100644 --- a/src/cp_measure/featurizer.py +++ b/src/cp_measure/featurizer.py @@ -20,6 +20,8 @@ import numpy as np +from cp_measure._sanitize import sanitize_masks + # Feature groups that only support 2D spatial data. _2D_ONLY = {"radial_distribution", "radial_zernikes", "zernike", "feret"} @@ -184,9 +186,10 @@ def featurize( masks : numpy.ndarray Integer-labeled masks with shape ``(M, H, W)`` or ``(M, Z, H, W)``. Must have the same ``ndim`` as *image*. - Background is 0; labels must be contiguous integers ``1..N`` - (standard cp_measure convention, see - ``skimage.segmentation.relabel_sequential``). + Background is 0; labels are positive integers. Non-contiguous or + arbitrary IDs (e.g. ``{1, 5, 17}``) are sanitized to ``1..N`` + internally (see :mod:`cp_measure._sanitize`); the input array is not + modified and output rows are reported against the original IDs. config : dict, optional Configuration dictionary produced by :func:`make_featurizer_config`. If ``None``, all features are enabled with default parameters. @@ -238,21 +241,23 @@ def featurize( columns: list[str] | None = None for mask_idx, object_name in enumerate(objects): - mask = masks[mask_idx] - # Assumes contiguous labels 1..max (standard cp_measure contract). - n_labels = int(mask.max()) - if n_labels == 0: + # Sanitize once, as early as possible: relabel arbitrary or + # non-contiguous label IDs to 1..N so every measurement function below + # sees clean labels, and keep the original IDs to report rows against. + # The input array is never mutated. See cp_measure._sanitize. + clean, ids = sanitize_masks(masks[mask_idx]) + if ids.size == 0: continue results: dict[str, np.ndarray] = {} for func, params in shape_feats: - results.update(func(mask, dummy_pixels, **params)) + results.update(func(clean, dummy_pixels, **params)) for ch_idx, ch_name in enumerate(channels): pixels = image[ch_idx] for func, params in channel_feats: - for key, values in func(mask, pixels, **params).items(): + for key, values in func(clean, pixels, **params).items(): results[f"{key}__{ch_name}"] = values n_ch = len(channels) @@ -262,7 +267,7 @@ def featurize( for key, values in func( pixels_1=image[ch_i], pixels_2=image[ch_j], - masks=mask, + masks=clean, **params, ).items(): results[f"{key}__{channels[ch_i]}__{channels[ch_j]}"] = values @@ -285,7 +290,7 @@ def featurize( all_blocks.append(block) all_rows.extend( - (image_id, object_name, label) for label in range(1, n_labels + 1) + (image_id, object_name, int(label)) for label in ids ) if not all_blocks: diff --git a/test/test_sanitize.py b/test/test_sanitize.py new file mode 100644 index 0000000..4610e6c --- /dev/null +++ b/test/test_sanitize.py @@ -0,0 +1,199 @@ +"""Tests for central non-contiguous label sanitation (cp_measure._sanitize).""" + +import numpy as np +import pytest + +from cp_measure._detect import HAS_NUMBA +from cp_measure._sanitize import _is_contiguous, sanitize_labels, sanitize_masks +from cp_measure.bulk import ( + get_core_measurements, + get_core_measurements_3d, + get_correlation_measurements, +) +from cp_measure.featurizer import featurize, make_featurizer_config + + +# -------------------------------------------------------------------------- +# sanitize_masks policy +# -------------------------------------------------------------------------- + + +def _three_objects(labels): + """3 squares in a 64x64 frame, labelled with the given (l1, l2, l3).""" + m = np.zeros((64, 64), dtype=np.int32) + m[5:19, 5:19] = labels[0] + m[5:19, 30:44] = labels[1] + m[35:49, 5:19] = labels[2] + return m + + +def test_gapped_labels_relabelled_to_1_n(): + gapped = _three_objects((1, 5, 17)) + clean, ids = sanitize_masks(gapped) + assert sorted(set(np.unique(clean)) - {0}) == [1, 2, 3] + assert ids.tolist() == [1, 5, 17] + + +def test_input_not_mutated(): + gapped = _three_objects((1, 5, 17)) + before = gapped.copy() + clean, _ids = sanitize_masks(gapped) + assert np.array_equal(gapped, before) + assert clean is not gapped # relabel path returns a fresh copy + + +def test_contiguous_fast_path_returns_input_unchanged(): + contig = _three_objects((1, 2, 3)) + clean, ids = sanitize_masks(contig) + assert clean is contig # no copy when already 1..N + assert ids.tolist() == [1, 2, 3] + + +def test_all_background(): + clean, ids = sanitize_masks(np.zeros((8, 8), dtype=np.int32)) + assert ids.size == 0 + + +def test_3d_gapped(): + m = np.zeros((4, 16, 16), dtype=np.int32) + m[0, 0:4, 0:4] = 3 + m[3, 10:14, 10:14] = 9 + clean, ids = sanitize_masks(m) + assert sorted(set(np.unique(clean)) - {0}) == [1, 2] + assert ids.tolist() == [3, 9] + + +def test_single_pixel(): + m = np.zeros((8, 8), dtype=np.int32) + m[3, 3] = 42 + clean, ids = sanitize_masks(m) + assert ids.tolist() == [42] + assert clean[3, 3] == 1 + + +def test_is_contiguous(): + assert _is_contiguous(_three_objects((1, 2, 3))) + assert not _is_contiguous(_three_objects((1, 5, 17))) + assert _is_contiguous(np.zeros((4, 4), dtype=np.int32)) + + +@pytest.mark.parametrize( + "bad", + [np.array([[0, -1, 2]]), np.array([[0.0, 1.0, 2.0]])], +) +def test_invalid_input_raises(bad): + with pytest.raises(ValueError): + sanitize_masks(bad) + + +# -------------------------------------------------------------------------- +# decorator coverage (every registry function is sanitized) +# -------------------------------------------------------------------------- + + +def _is_sanitized(fn): + # unwrap functools.partial (legacy binding) before checking the flag + return bool( + getattr(fn, "_sanitized", False) + or getattr(getattr(fn, "func", None), "_sanitized", False) + ) + + +def test_all_single_mask_registry_funcs_sanitized(): + registries = [ + get_core_measurements(), + get_core_measurements(legacy=True), # partial-wrapped intensity + get_core_measurements_3d(), + get_correlation_measurements(), + ] + for reg in registries: + for name, fn in reg.items(): + assert _is_sanitized(fn), f"{name} is not @sanitize_labels-wrapped" + + +@pytest.mark.skipif(not HAS_NUMBA, reason="numba not installed") +def test_numba_registry_funcs_sanitized(): + import cp_measure + + cp_measure.set_accelerator("numba") + try: + for name, fn in get_core_measurements().items(): + assert _is_sanitized(fn), f"numba {name} is not sanitized" + finally: + cp_measure.set_accelerator(None) + + +def test_decorator_leaves_unknown_signature_untouched(): + def two_mask(masks1, masks2): # multimask-style: no recognised label arg + return masks1 + + assert sanitize_labels(two_mask) is two_mask + + +# -------------------------------------------------------------------------- +# per-function: gapped result == contiguous baseline (direct-call coverage) +# -------------------------------------------------------------------------- + +CONTIG = _three_objects((1, 2, 3)) +GAPPED = _three_objects((1, 5, 17)) # same geometry, ascending -> same ranks +_RNG = np.random.default_rng(0) +PIX1 = _RNG.random((64, 64)) +PIX2 = _RNG.random((64, 64)) + + +def _allclose_dicts(a, b): + assert list(a.keys()) == list(b.keys()) + for k in a: + np.testing.assert_allclose( + np.asarray(a[k], dtype=float), + np.asarray(b[k], dtype=float), + equal_nan=True, + err_msg=f"mismatch in feature {k!r}", + ) + + +@pytest.mark.parametrize("name", list(get_core_measurements().keys())) +def test_core_func_gapped_matches_contiguous(name): + func = get_core_measurements()[name] + _allclose_dicts(func(GAPPED, PIX1), func(CONTIG, PIX1)) + + +@pytest.mark.parametrize("name", list(get_correlation_measurements().keys())) +def test_correlation_func_gapped_matches_contiguous(name): + func = get_correlation_measurements()[name] + _allclose_dicts( + func(pixels_1=PIX1, pixels_2=PIX2, masks=GAPPED), + func(pixels_1=PIX1, pixels_2=PIX2, masks=CONTIG), + ) + + +# -------------------------------------------------------------------------- +# featurizer end-to-end +# -------------------------------------------------------------------------- + + +def test_featurizer_gapped_reports_original_ids(): + image = np.stack([PIX1, PIX2], axis=0) + config = make_featurizer_config(["DNA", "ER"]) + + gapped_masks = GAPPED[np.newaxis] + contig_masks = CONTIG[np.newaxis] + + data_g, cols_g, rows_g = featurize(image, gapped_masks, config) + data_c, cols_c, rows_c = featurize(image, contig_masks, config) + + assert cols_g == cols_c + # rows carry the ORIGINAL ids, not 1..N + assert rows_g == [(None, "object", 1), (None, "object", 5), (None, "object", 17)] + assert rows_c == [(None, "object", 1), (None, "object", 2), (None, "object", 3)] + # and the actual feature values are identical (same geometry) + np.testing.assert_allclose(data_g, data_c, equal_nan=True) + + +def test_featurizer_does_not_mutate_input(): + image = np.stack([PIX1, PIX2], axis=0) + config = make_featurizer_config(["DNA", "ER"]) + masks = GAPPED[np.newaxis].copy() + before = masks.copy() + featurize(image, masks, config) + assert np.array_equal(masks, before) From 37cde59a8eda0476a3329bd385cc1326b8b13c47 Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 14:34:01 +0200 Subject: [PATCH 02/11] refactor(sanitize): leaner policy, bool support, fix overlap coverage MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Review-driven cleanup of the sanitation layer: - collapse sanitize_masks to a single numpy.unique pass: derive the contiguity test, negative-value guard and original IDs from one sorted array; drop the separate _is_contiguous helper and the redundant max()/min()/second-unique passes - accept boolean masks (treated as one object) instead of raising — restores the pre-existing behaviour of the raw measurement functions - decorate get_correlation_overlap (the one public correlation entry that was missed), so direct callers get correct results on gapped labels - trim the over-long module/function docstrings and the featurizer comment - tests: parametrize the policy/edge cases, exercise the rank remap with non-monotonic labels, assert real multimask functions stay unsanitized, add a bool-mask case (net leaner, 109 passing) Co-Authored-By: Claude Opus 4.8 (1M context) --- src/cp_measure/_sanitize.py | 91 +++------- src/cp_measure/core/measurecolocalization.py | 1 + .../core/numba/measureobjectintensity.py | 1 - src/cp_measure/featurizer.py | 9 +- test/test_sanitize.py | 157 ++++++++---------- 5 files changed, 92 insertions(+), 167 deletions(-) diff --git a/src/cp_measure/_sanitize.py b/src/cp_measure/_sanitize.py index 4e8ba7f..fdb1f01 100644 --- a/src/cp_measure/_sanitize.py +++ b/src/cp_measure/_sanitize.py @@ -1,26 +1,10 @@ -"""Central input sanitation for non-contiguous mask label IDs. +"""Central sanitation of non-contiguous mask label IDs. -cp_measure's measurement functions assume object labels are the contiguous -integers ``1..N`` (the convention documented on -:func:`cp_measure.featurizer.featurize`). Real segmentations do not always -honour that: labels may have gaps (``{1, 5, 17}``) or arbitrary values. This -module provides the single policy that maps arbitrary positive-integer labels to -``1..N`` *before* any math runs, so every downstream function can assume clean -labels and focus on the calculation. - -Two pieces: - -* :func:`sanitize_masks` — the policy. Returns ``(clean, ids)`` where ``clean`` - has labels ``1..N`` and ``ids[i]`` is the *original* label of rank ``i + 1``. - It never mutates the caller's array (the relabel path returns a fresh copy; the - already-clean fast path returns the input unchanged). -* :func:`sanitize_labels` — a thin decorator that applies the policy to whichever - argument of a ``get_*`` function holds the label image, so direct callers get - the same guarantee as the featurizer. - -Relabelling is cheap relative to a single feature (<1 % of a featurized image), -so the featurizer sanitizes once up front and the per-function decorator is a -cheap idempotent guard (a single :func:`numpy.unique`) for direct callers. +cp_measure's measurement functions assume labels are the contiguous integers +``1..N`` (see :func:`cp_measure.featurizer.featurize`). Real segmentations may +use gaps (``{1, 5, 17}``) or arbitrary values; this maps them to ``1..N`` before +any math runs and reports results against the original IDs, without mutating the +caller's array. """ import functools @@ -31,67 +15,34 @@ from numpy.typing import NDArray from skimage.segmentation import relabel_sequential -# Argument names that, across the ``get_*`` functions, hold the label image. +# Argument names that hold the label image across the ``get_*`` functions. _MASK_PARAMS = ("masks", "labels", "mask") -def _is_contiguous(masks: NDArray) -> bool: - """Are the positive labels exactly ``1..N``? - - Uses :func:`numpy.unique` rather than a dense ``bincount`` so the check is - safe and bounded for any dtype, negative values, or large label values. - """ - unique = numpy.unique(masks) - positive = unique[unique > 0] - return positive.size == 0 or (positive[0] == 1 and positive[-1] == positive.size) - - def sanitize_masks(masks: NDArray) -> tuple[NDArray, NDArray[numpy.int64]]: - """Relabel arbitrary positive labels to contiguous ``1..N``. + """Relabel positive labels to contiguous ``1..N``. - Parameters - ---------- - masks - Integer label array (any number of dimensions). Background is ``0``. - - Returns - ------- - clean - Array with labels ``1..N`` in ascending original-label order. The - already-contiguous input is returned unchanged (no copy); otherwise a - fresh relabelled copy is returned. The input is never mutated. - ids - ``ids[i]`` is the original label whose sanitized value is ``i + 1`` - (ascending). Use it to report results against the caller's IDs. - - Raises - ------ - ValueError - If ``masks`` is not an integer array or contains negative values. + Returns ``(clean, ids)`` where ``clean`` has labels ``1..N`` (the input is + returned unchanged when already contiguous) and ``ids[i]`` is the original + label of rank ``i + 1``. Never mutates the input. Raises ``ValueError`` for + non-integer (non-boolean) or negative labels. """ - if not numpy.issubdtype(masks.dtype, numpy.integer): + if masks.dtype != bool and not numpy.issubdtype(masks.dtype, numpy.integer): raise ValueError(f"labels must be an integer array, got dtype {masks.dtype!r}") - if masks.size and masks.min() < 0: + labels = numpy.unique(masks) + if labels.size and labels[0] < 0: raise ValueError("labels must be non-negative") - - if _is_contiguous(masks): - n = int(masks.max(initial=0)) - return masks, numpy.arange(1, n + 1, dtype=numpy.int64) - + ids = labels[labels > 0].astype(numpy.int64) + if ids.size == 0 or (ids[0] == 1 and ids[-1] == ids.size): + return masks, ids # already contiguous (or empty): no copy clean, _forward, _inverse = relabel_sequential(masks) - ids = numpy.unique(masks) - ids = ids[ids > 0].astype(numpy.int64) return clean, ids def sanitize_labels(func: Callable) -> Callable: - """Decorate a ``get_*`` function to sanitize its label argument. - - The label argument is detected by name (:data:`_MASK_PARAMS`), so the - decorator is position-independent. Functions without a recognised label - argument (e.g. the two-mask ``multimask`` functions, which take - ``masks1``/``masks2``) are returned unchanged. - """ + """Decorate a ``get_*`` function to sanitize its label argument (found by + name in :data:`_MASK_PARAMS`); functions with no such argument (e.g. the + two-mask multimask functions) are returned unchanged.""" sig = inspect.signature(func) param = next((name for name in _MASK_PARAMS if name in sig.parameters), None) if param is None: diff --git a/src/cp_measure/core/measurecolocalization.py b/src/cp_measure/core/measurecolocalization.py index 221e8d9..993f5df 100644 --- a/src/cp_measure/core/measurecolocalization.py +++ b/src/cp_measure/core/measurecolocalization.py @@ -627,6 +627,7 @@ def get_correlation_costes( return {f"{F_COSTES_FORMAT}_1": c1, f"{F_COSTES_FORMAT}_2": c2} +@sanitize_labels def get_correlation_overlap( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], diff --git a/src/cp_measure/core/numba/measureobjectintensity.py b/src/cp_measure/core/numba/measureobjectintensity.py index aaff39c..6ae43b8 100644 --- a/src/cp_measure/core/numba/measureobjectintensity.py +++ b/src/cp_measure/core/numba/measureobjectintensity.py @@ -15,7 +15,6 @@ from numpy.typing import NDArray from cp_measure._sanitize import sanitize_labels - from cp_measure.core.measureobjectintensity import ( C_LOCATION, INTEGRATED_INTENSITY, diff --git a/src/cp_measure/featurizer.py b/src/cp_measure/featurizer.py index dccb3a6..2c0fe7f 100644 --- a/src/cp_measure/featurizer.py +++ b/src/cp_measure/featurizer.py @@ -241,10 +241,7 @@ def featurize( columns: list[str] | None = None for mask_idx, object_name in enumerate(objects): - # Sanitize once, as early as possible: relabel arbitrary or - # non-contiguous label IDs to 1..N so every measurement function below - # sees clean labels, and keep the original IDs to report rows against. - # The input array is never mutated. See cp_measure._sanitize. + # Relabel arbitrary IDs to 1..N once up front; keep originals for rows. clean, ids = sanitize_masks(masks[mask_idx]) if ids.size == 0: continue @@ -289,9 +286,7 @@ def featurize( block = np.column_stack([results[c] for c in columns]) all_blocks.append(block) - all_rows.extend( - (image_id, object_name, int(label)) for label in ids - ) + all_rows.extend((image_id, object_name, int(label)) for label in ids) if not all_blocks: raise ValueError("all masks have no labels (all zeros)") diff --git a/test/test_sanitize.py b/test/test_sanitize.py index 4610e6c..885ee21 100644 --- a/test/test_sanitize.py +++ b/test/test_sanitize.py @@ -4,20 +4,17 @@ import pytest from cp_measure._detect import HAS_NUMBA -from cp_measure._sanitize import _is_contiguous, sanitize_labels, sanitize_masks +from cp_measure._sanitize import sanitize_masks from cp_measure.bulk import ( get_core_measurements, get_core_measurements_3d, get_correlation_measurements, + get_multimask_measurements, ) +from cp_measure.core.measurecolocalization import get_correlation_overlap from cp_measure.featurizer import featurize, make_featurizer_config -# -------------------------------------------------------------------------- -# sanitize_masks policy -# -------------------------------------------------------------------------- - - def _three_objects(labels): """3 squares in a 64x64 frame, labelled with the given (l1, l2, l3).""" m = np.zeros((64, 64), dtype=np.int32) @@ -27,68 +24,59 @@ def _three_objects(labels): return m -def test_gapped_labels_relabelled_to_1_n(): - gapped = _three_objects((1, 5, 17)) - clean, ids = sanitize_masks(gapped) - assert sorted(set(np.unique(clean)) - {0}) == [1, 2, 3] - assert ids.tolist() == [1, 5, 17] - - -def test_input_not_mutated(): - gapped = _three_objects((1, 5, 17)) - before = gapped.copy() - clean, _ids = sanitize_masks(gapped) - assert np.array_equal(gapped, before) - assert clean is not gapped # relabel path returns a fresh copy - +# --- sanitize_masks policy ------------------------------------------------- -def test_contiguous_fast_path_returns_input_unchanged(): - contig = _three_objects((1, 2, 3)) - clean, ids = sanitize_masks(contig) - assert clean is contig # no copy when already 1..N - assert ids.tolist() == [1, 2, 3] +_M3D = np.zeros((4, 16, 16), dtype=np.int32) +_M3D[0, 0:4, 0:4] = 3 +_M3D[3, 10:14, 10:14] = 9 +_MPIX = np.zeros((8, 8), dtype=np.int32) +_MPIX[3, 3] = 42 -def test_all_background(): - clean, ids = sanitize_masks(np.zeros((8, 8), dtype=np.int32)) - assert ids.size == 0 - - -def test_3d_gapped(): - m = np.zeros((4, 16, 16), dtype=np.int32) - m[0, 0:4, 0:4] = 3 - m[3, 10:14, 10:14] = 9 - clean, ids = sanitize_masks(m) - assert sorted(set(np.unique(clean)) - {0}) == [1, 2] - assert ids.tolist() == [3, 9] +@pytest.mark.parametrize( + "mask, clean_labels, ids", + [ + (_three_objects((1, 2, 3)), [1, 2, 3], [1, 2, 3]), # contiguous + (_three_objects((1, 17, 5)), [1, 2, 3], [1, 5, 17]), # gapped + (_M3D, [1, 2], [3, 9]), # 3D + (_MPIX, [1], [42]), # single pixel + (np.zeros((8, 8), np.int32), [], []), # all background + ], +) +def test_sanitize_relabels(mask, clean_labels, ids): + clean, got_ids = sanitize_masks(mask) + assert sorted(set(np.unique(clean)) - {0}) == clean_labels + assert got_ids.tolist() == ids -def test_single_pixel(): - m = np.zeros((8, 8), dtype=np.int32) - m[3, 3] = 42 - clean, ids = sanitize_masks(m) - assert ids.tolist() == [42] - assert clean[3, 3] == 1 +@pytest.mark.parametrize( + "mask, same_object", + [ + (_three_objects((1, 2, 3)), True), # contiguous -> returned as-is + (_three_objects((1, 17, 5)), False), # relabelled -> fresh copy + ], +) +def test_no_mutation_and_copy_policy(mask, same_object): + before = mask.copy() + clean, _ids = sanitize_masks(mask) + assert (clean is mask) == same_object + assert np.array_equal(mask, before) -def test_is_contiguous(): - assert _is_contiguous(_three_objects((1, 2, 3))) - assert not _is_contiguous(_three_objects((1, 5, 17))) - assert _is_contiguous(np.zeros((4, 4), dtype=np.int32)) +def test_bool_mask_is_single_object(): + m = np.zeros((8, 8), dtype=bool) + m[2:5, 2:5] = True + _clean, ids = sanitize_masks(m) + assert ids.tolist() == [1] -@pytest.mark.parametrize( - "bad", - [np.array([[0, -1, 2]]), np.array([[0.0, 1.0, 2.0]])], -) +@pytest.mark.parametrize("bad", [np.array([[0, -1, 2]]), np.array([[0.0, 1.0, 2.0]])]) def test_invalid_input_raises(bad): with pytest.raises(ValueError): sanitize_masks(bad) -# -------------------------------------------------------------------------- -# decorator coverage (every registry function is sanitized) -# -------------------------------------------------------------------------- +# --- decorator coverage ---------------------------------------------------- def _is_sanitized(fn): @@ -100,13 +88,12 @@ def _is_sanitized(fn): def test_all_single_mask_registry_funcs_sanitized(): - registries = [ + for reg in ( get_core_measurements(), get_core_measurements(legacy=True), # partial-wrapped intensity get_core_measurements_3d(), get_correlation_measurements(), - ] - for reg in registries: + ): for name, fn in reg.items(): assert _is_sanitized(fn), f"{name} is not @sanitize_labels-wrapped" @@ -123,70 +110,62 @@ def test_numba_registry_funcs_sanitized(): cp_measure.set_accelerator(None) -def test_decorator_leaves_unknown_signature_untouched(): - def two_mask(masks1, masks2): # multimask-style: no recognised label arg - return masks1 - - assert sanitize_labels(two_mask) is two_mask - +def test_multimask_funcs_left_unsanitized(): + # two-mask functions are out of scope and must not be wrapped + assert all(not _is_sanitized(fn) for fn in get_multimask_measurements().values()) -# -------------------------------------------------------------------------- -# per-function: gapped result == contiguous baseline (direct-call coverage) -# -------------------------------------------------------------------------- -CONTIG = _three_objects((1, 2, 3)) -GAPPED = _three_objects((1, 5, 17)) # same geometry, ascending -> same ranks +# --- per-function: gapped result == contiguous baseline -------------------- +# CONTIG and GAPPED share geometry and rank order, so a correct relabel makes +# every feature's output identical. Labels are non-monotonic vs spatial order, +# so this also exercises the rank remap, not just value substitution. +CONTIG = _three_objects((1, 3, 2)) +GAPPED = _three_objects((1, 17, 5)) _RNG = np.random.default_rng(0) PIX1 = _RNG.random((64, 64)) PIX2 = _RNG.random((64, 64)) -def _allclose_dicts(a, b): +def _assert_close_dicts(a, b): assert list(a.keys()) == list(b.keys()) for k in a: np.testing.assert_allclose( - np.asarray(a[k], dtype=float), - np.asarray(b[k], dtype=float), + np.asarray(a[k], float), + np.asarray(b[k], float), equal_nan=True, err_msg=f"mismatch in feature {k!r}", ) -@pytest.mark.parametrize("name", list(get_core_measurements().keys())) +@pytest.mark.parametrize("name", list(get_core_measurements())) def test_core_func_gapped_matches_contiguous(name): func = get_core_measurements()[name] - _allclose_dicts(func(GAPPED, PIX1), func(CONTIG, PIX1)) + _assert_close_dicts(func(GAPPED, PIX1), func(CONTIG, PIX1)) -@pytest.mark.parametrize("name", list(get_correlation_measurements().keys())) -def test_correlation_func_gapped_matches_contiguous(name): - func = get_correlation_measurements()[name] - _allclose_dicts( +@pytest.mark.parametrize( + "func", + list(get_correlation_measurements().values()) + [get_correlation_overlap], + ids=lambda f: f.__name__, +) +def test_correlation_func_gapped_matches_contiguous(func): + _assert_close_dicts( func(pixels_1=PIX1, pixels_2=PIX2, masks=GAPPED), func(pixels_1=PIX1, pixels_2=PIX2, masks=CONTIG), ) -# -------------------------------------------------------------------------- -# featurizer end-to-end -# -------------------------------------------------------------------------- +# --- featurizer end-to-end ------------------------------------------------- def test_featurizer_gapped_reports_original_ids(): image = np.stack([PIX1, PIX2], axis=0) config = make_featurizer_config(["DNA", "ER"]) - - gapped_masks = GAPPED[np.newaxis] - contig_masks = CONTIG[np.newaxis] - - data_g, cols_g, rows_g = featurize(image, gapped_masks, config) - data_c, cols_c, rows_c = featurize(image, contig_masks, config) - + data_g, cols_g, rows_g = featurize(image, GAPPED[np.newaxis], config) + data_c, cols_c, rows_c = featurize(image, CONTIG[np.newaxis], config) assert cols_g == cols_c - # rows carry the ORIGINAL ids, not 1..N assert rows_g == [(None, "object", 1), (None, "object", 5), (None, "object", 17)] assert rows_c == [(None, "object", 1), (None, "object", 2), (None, "object", 3)] - # and the actual feature values are identical (same geometry) np.testing.assert_allclose(data_g, data_c, equal_nan=True) From 8633bab848c8efe08cd891a04a165e47d0b3b327 Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 15:05:19 +0200 Subject: [PATCH 03/11] perf(sanitize): cheaper contiguity check (max + bincount, not unique) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit The decorator ran numpy.unique over the whole image on every call (~10 ms on 1024²), a fixed cost that dominated cheap functions and large/few-object images (manders_fold fell to 0.40×). Replace it with a max() gate plus a bincount presence check (~4x faster), falling back to unique only for pathologically large labels. relabel_sequential still handles the rare gapped path. Measured on 1024²: manders_fold 0.74→0.91×, intensity 0.97×, feret 0.99×; expensive functions (zernike/sizeshape/texture) were already ~1.0×. Negatives still raise; bit-identical results. 109 tests pass. Co-Authored-By: Claude Opus 4.8 (1M context) --- src/cp_measure/_sanitize.py | 18 ++++++++++++------ 1 file changed, 12 insertions(+), 6 deletions(-) diff --git a/src/cp_measure/_sanitize.py b/src/cp_measure/_sanitize.py index fdb1f01..892ce51 100644 --- a/src/cp_measure/_sanitize.py +++ b/src/cp_measure/_sanitize.py @@ -29,12 +29,18 @@ def sanitize_masks(masks: NDArray) -> tuple[NDArray, NDArray[numpy.int64]]: """ if masks.dtype != bool and not numpy.issubdtype(masks.dtype, numpy.integer): raise ValueError(f"labels must be an integer array, got dtype {masks.dtype!r}") - labels = numpy.unique(masks) - if labels.size and labels[0] < 0: - raise ValueError("labels must be non-negative") - ids = labels[labels > 0].astype(numpy.int64) - if ids.size == 0 or (ids[0] == 1 and ids[-1] == ids.size): - return masks, ids # already contiguous (or empty): no copy + mx = int(masks.max(initial=0)) + if mx == 0: + return masks, numpy.empty(0, dtype=numpy.int64) + # Cheap contiguity check (~4x faster than numpy.unique); bincount needs a + # bounded range, so fall back to unique for pathologically large labels. + if mx <= masks.size: + ids = numpy.flatnonzero(numpy.bincount(masks.ravel(), minlength=mx + 1)) + else: + ids = numpy.unique(masks) + ids = ids[ids > 0].astype(numpy.int64) + if ids.size == mx: + return masks, ids # already 1..N: no copy clean, _forward, _inverse = relabel_sequential(masks) return clean, ids From 186dd4b69d607bd6de75550a2598f638ccb11e8e Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 15:16:59 +0200 Subject: [PATCH 04/11] perf(featurizer): sanitize each mask once, call raw feature functions featurize() sanitized each mask up front and then called the @sanitize_labels- decorated functions, so every feature re-ran the contiguity scan on the already-clean mask (~12x per object-type). Unwrap the decorator (via the __wrapped__ the wrapper already sets, partial-aware for the legacy path) so the featurizer pays the sanitation cost exactly once per mask. Direct callers still get the decorator's guarantee. Test asserts sanitize_masks is called once per object-type mask. 110 pass. Co-Authored-By: Claude Opus 4.8 (1M context) --- src/cp_measure/featurizer.py | 13 +++++++++++++ test/test_sanitize.py | 21 +++++++++++++++++++++ 2 files changed, 34 insertions(+) diff --git a/src/cp_measure/featurizer.py b/src/cp_measure/featurizer.py index 2c0fe7f..d6b2dd1 100644 --- a/src/cp_measure/featurizer.py +++ b/src/cp_measure/featurizer.py @@ -15,6 +15,7 @@ from __future__ import annotations +import functools import itertools import warnings @@ -26,6 +27,13 @@ _2D_ONLY = {"radial_distribution", "radial_zernikes", "zernike", "feret"} +def _unwrap(fn): + """Strip the @sanitize_labels wrapper; the featurizer sanitizes once up front.""" + if isinstance(fn, functools.partial): + return functools.partial(_unwrap(fn.func), *fn.args, **fn.keywords) + return getattr(fn, "__wrapped__", fn) + + def make_featurizer_config( channels: list[str] | None = None, *, @@ -227,6 +235,11 @@ def featurize( ) corr_funcs = get_correlation_measurements() + # Each mask is sanitized once below, so call the underlying implementations + # directly instead of re-sanitizing in every feature's decorator. + core_funcs = {name: _unwrap(fn) for name, fn in core_funcs.items()} + corr_funcs = {name: _unwrap(fn) for name, fn in corr_funcs.items()} + if is_3d: _warn_2d_only_in_3d(config) channel_feats = _collect_channel_features(config, core_funcs) diff --git a/test/test_sanitize.py b/test/test_sanitize.py index 885ee21..774ce45 100644 --- a/test/test_sanitize.py +++ b/test/test_sanitize.py @@ -169,6 +169,27 @@ def test_featurizer_gapped_reports_original_ids(): np.testing.assert_allclose(data_g, data_c, equal_nan=True) +def test_featurizer_sanitizes_once_per_mask(monkeypatch): + # featurize() sanitizes up front and calls the raw functions, so the cost is + # paid once per object-type mask — not again in every feature's decorator. + import cp_measure._sanitize as _san + import cp_measure.featurizer as _fz + + calls = [] + real = _san.sanitize_masks + + def counting(m): + calls.append(1) + return real(m) + + monkeypatch.setattr(_san, "sanitize_masks", counting) + monkeypatch.setattr(_fz, "sanitize_masks", counting) + image = np.stack([PIX1, PIX2], axis=0) + config = make_featurizer_config(["DNA", "ER"]) + featurize(image, GAPPED[np.newaxis], config) + assert len(calls) == 1 + + def test_featurizer_does_not_mutate_input(): image = np.stack([PIX1, PIX2], axis=0) config = make_featurizer_config(["DNA", "ER"]) From c9d292222146741ed883989ce81291f97e7bf6f9 Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 15:40:30 +0200 Subject: [PATCH 05/11] test(sanitize): trim to the essentials The per-function equivalence tests parametrized over all ~13 measurement functions were redundant: the decorator is uniform, so the coverage meta-test (every public function wrapped) + one representative direct call (zernike, which crashes on gaps without sanitation) + the featurizer end-to-end test cover the same ground. Also merged the mutation/copy cases, dropped the separate numba meta-test, and removed a dead second monkeypatch. 178 -> 102 lines. Co-Authored-By: Claude Opus 4.8 (1M context) --- test/test_sanitize.py | 183 ++++++++++-------------------------------- 1 file changed, 43 insertions(+), 140 deletions(-) diff --git a/test/test_sanitize.py b/test/test_sanitize.py index 774ce45..cec495d 100644 --- a/test/test_sanitize.py +++ b/test/test_sanitize.py @@ -3,11 +3,9 @@ import numpy as np import pytest -from cp_measure._detect import HAS_NUMBA from cp_measure._sanitize import sanitize_masks from cp_measure.bulk import ( get_core_measurements, - get_core_measurements_3d, get_correlation_measurements, get_multimask_measurements, ) @@ -24,176 +22,81 @@ def _three_objects(labels): return m -# --- sanitize_masks policy ------------------------------------------------- - _M3D = np.zeros((4, 16, 16), dtype=np.int32) -_M3D[0, 0:4, 0:4] = 3 +_M3D[0, :4, :4] = 3 _M3D[3, 10:14, 10:14] = 9 -_MPIX = np.zeros((8, 8), dtype=np.int32) -_MPIX[3, 3] = 42 @pytest.mark.parametrize( - "mask, clean_labels, ids", + "mask, clean, ids", [ - (_three_objects((1, 2, 3)), [1, 2, 3], [1, 2, 3]), # contiguous - (_three_objects((1, 17, 5)), [1, 2, 3], [1, 5, 17]), # gapped + (_three_objects((1, 17, 5)), [1, 2, 3], [1, 5, 17]), # gapped, non-monotonic (_M3D, [1, 2], [3, 9]), # 3D - (_MPIX, [1], [42]), # single pixel (np.zeros((8, 8), np.int32), [], []), # all background + (np.array([[False, True], [True, False]]), [1], [1]), # bool -> one object ], ) -def test_sanitize_relabels(mask, clean_labels, ids): - clean, got_ids = sanitize_masks(mask) - assert sorted(set(np.unique(clean)) - {0}) == clean_labels - assert got_ids.tolist() == ids +def test_sanitize_relabels(mask, clean, ids): + out, out_ids = sanitize_masks(mask) + assert sorted(set(np.unique(out)) - {0}) == clean + assert out_ids.tolist() == ids -@pytest.mark.parametrize( - "mask, same_object", - [ - (_three_objects((1, 2, 3)), True), # contiguous -> returned as-is - (_three_objects((1, 17, 5)), False), # relabelled -> fresh copy - ], -) -def test_no_mutation_and_copy_policy(mask, same_object): - before = mask.copy() - clean, _ids = sanitize_masks(mask) - assert (clean is mask) == same_object - assert np.array_equal(mask, before) +def test_no_mutation_and_copy_policy(): + contig, gapped = _three_objects((1, 2, 3)), _three_objects((1, 17, 5)) + assert sanitize_masks(contig)[0] is contig # already 1..N: no copy + before = gapped.copy() + assert sanitize_masks(gapped)[0] is not gapped # relabelled: fresh copy + assert np.array_equal(gapped, before) # input untouched -def test_bool_mask_is_single_object(): - m = np.zeros((8, 8), dtype=bool) - m[2:5, 2:5] = True - _clean, ids = sanitize_masks(m) - assert ids.tolist() == [1] - - -@pytest.mark.parametrize("bad", [np.array([[0, -1, 2]]), np.array([[0.0, 1.0, 2.0]])]) +@pytest.mark.parametrize("bad", [np.array([[0, -1, 2]]), np.array([[1.0, 2.0]])]) def test_invalid_input_raises(bad): with pytest.raises(ValueError): sanitize_masks(bad) -# --- decorator coverage ---------------------------------------------------- - - def _is_sanitized(fn): - # unwrap functools.partial (legacy binding) before checking the flag return bool( getattr(fn, "_sanitized", False) or getattr(getattr(fn, "func", None), "_sanitized", False) ) -def test_all_single_mask_registry_funcs_sanitized(): - for reg in ( - get_core_measurements(), - get_core_measurements(legacy=True), # partial-wrapped intensity - get_core_measurements_3d(), - get_correlation_measurements(), - ): - for name, fn in reg.items(): - assert _is_sanitized(fn), f"{name} is not @sanitize_labels-wrapped" - - -@pytest.mark.skipif(not HAS_NUMBA, reason="numba not installed") -def test_numba_registry_funcs_sanitized(): - import cp_measure - - cp_measure.set_accelerator("numba") - try: - for name, fn in get_core_measurements().items(): - assert _is_sanitized(fn), f"numba {name} is not sanitized" - finally: - cp_measure.set_accelerator(None) - - -def test_multimask_funcs_left_unsanitized(): - # two-mask functions are out of scope and must not be wrapped - assert all(not _is_sanitized(fn) for fn in get_multimask_measurements().values()) - - -# --- per-function: gapped result == contiguous baseline -------------------- -# CONTIG and GAPPED share geometry and rank order, so a correct relabel makes -# every feature's output identical. Labels are non-monotonic vs spatial order, -# so this also exercises the rank remap, not just value substitution. -CONTIG = _three_objects((1, 3, 2)) -GAPPED = _three_objects((1, 17, 5)) -_RNG = np.random.default_rng(0) -PIX1 = _RNG.random((64, 64)) -PIX2 = _RNG.random((64, 64)) - - -def _assert_close_dicts(a, b): - assert list(a.keys()) == list(b.keys()) - for k in a: - np.testing.assert_allclose( - np.asarray(a[k], float), - np.asarray(b[k], float), - equal_nan=True, - err_msg=f"mismatch in feature {k!r}", - ) - - -@pytest.mark.parametrize("name", list(get_core_measurements())) -def test_core_func_gapped_matches_contiguous(name): - func = get_core_measurements()[name] - _assert_close_dicts(func(GAPPED, PIX1), func(CONTIG, PIX1)) - - -@pytest.mark.parametrize( - "func", - list(get_correlation_measurements().values()) + [get_correlation_overlap], - ids=lambda f: f.__name__, -) -def test_correlation_func_gapped_matches_contiguous(func): - _assert_close_dicts( - func(pixels_1=PIX1, pixels_2=PIX2, masks=GAPPED), - func(pixels_1=PIX1, pixels_2=PIX2, masks=CONTIG), +def test_public_funcs_sanitized_multimask_excluded(): + public = ( + *get_core_measurements().values(), + *get_core_measurements(legacy=True).values(), # partial-wrapped intensity + *get_correlation_measurements().values(), + get_correlation_overlap, # public, not in the registry ) + assert all(_is_sanitized(fn) for fn in public) + # two-mask functions are out of scope and must stay unwrapped + assert not any(_is_sanitized(fn) for fn in get_multimask_measurements().values()) -# --- featurizer end-to-end ------------------------------------------------- +def test_decorated_func_handles_gapped_labels(): + # get_zernike raises on gapped labels without sanitation; the decorator fixes it. + px = np.random.default_rng(0).random((64, 64)) + zernike = get_core_measurements()["zernike"] + gapped = zernike(_three_objects((1, 17, 5)), px) + contig = zernike(_three_objects((1, 3, 2)), px) + for key in gapped: + np.testing.assert_allclose(gapped[key], contig[key], equal_nan=True) -def test_featurizer_gapped_reports_original_ids(): - image = np.stack([PIX1, PIX2], axis=0) - config = make_featurizer_config(["DNA", "ER"]) - data_g, cols_g, rows_g = featurize(image, GAPPED[np.newaxis], config) - data_c, cols_c, rows_c = featurize(image, CONTIG[np.newaxis], config) - assert cols_g == cols_c - assert rows_g == [(None, "object", 1), (None, "object", 5), (None, "object", 17)] - assert rows_c == [(None, "object", 1), (None, "object", 2), (None, "object", 3)] - np.testing.assert_allclose(data_g, data_c, equal_nan=True) - - -def test_featurizer_sanitizes_once_per_mask(monkeypatch): - # featurize() sanitizes up front and calls the raw functions, so the cost is - # paid once per object-type mask — not again in every feature's decorator. - import cp_measure._sanitize as _san - import cp_measure.featurizer as _fz +def test_featurizer_uses_original_ids_and_sanitizes_once(monkeypatch): + import cp_measure.featurizer as fz calls = [] - real = _san.sanitize_masks - - def counting(m): - calls.append(1) - return real(m) - - monkeypatch.setattr(_san, "sanitize_masks", counting) - monkeypatch.setattr(_fz, "sanitize_masks", counting) - image = np.stack([PIX1, PIX2], axis=0) + real = fz.sanitize_masks + monkeypatch.setattr(fz, "sanitize_masks", lambda m: (calls.append(1), real(m))[1]) + img = np.random.default_rng(0).random((2, 64, 64)) config = make_featurizer_config(["DNA", "ER"]) - featurize(image, GAPPED[np.newaxis], config) - assert len(calls) == 1 - -def test_featurizer_does_not_mutate_input(): - image = np.stack([PIX1, PIX2], axis=0) - config = make_featurizer_config(["DNA", "ER"]) - masks = GAPPED[np.newaxis].copy() - before = masks.copy() - featurize(image, masks, config) - assert np.array_equal(masks, before) + data_g, _, rows_g = featurize(img, _three_objects((1, 17, 5))[np.newaxis], config) + assert len(calls) == 1 # cost paid once per object-type mask + data_c, _, _ = featurize(img, _three_objects((1, 3, 2))[np.newaxis], config) + # original IDs in the rows, and identical geometry -> identical values + assert rows_g == [(None, "object", 1), (None, "object", 5), (None, "object", 17)] + np.testing.assert_allclose(data_g, data_c, equal_nan=True) From b21be05ee25e4f80404f71f45762777f5443f3ac Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 15:58:51 +0200 Subject: [PATCH 06/11] refactor: drop per-function label-normalization now centralized With sanitize_masks guaranteeing contiguous 1..N at every entry, the gap-tolerant label bookkeeping inside the feature functions is dead weight. Remove it: - texture/granularity/radial_distribution/radial_zernikes: drop the internal `numpy.unique(masks)` enumeration (a full-image sort, redundant after the decorator's pass). Use len(props) / range_.size / arange(1, max+1) instead. - numba intensity: labels are 1..N, so the label->index map is just `label - 1`. Drop `lut` from flatten_numba (and the edge path) and compute it inline; nobjects = masks.max(). - delete primitives/segment.py: label_to_idx_lut had no remaining consumers. All outputs verified bit-identical on contiguous input (numpy + numba), numba gapped input maps correctly, backend-correctness + full suite (91) pass. Co-Authored-By: Claude Opus 4.8 (1M context) --- src/cp_measure/core/measuregranularity.py | 4 +- .../measureobjectintensitydistribution.py | 8 ++-- src/cp_measure/core/measuretexture.py | 5 +-- .../core/numba/measureobjectintensity.py | 6 +-- src/cp_measure/primitives/__init__.py | 10 ++--- src/cp_measure/primitives/_segment_numba.py | 6 ++- src/cp_measure/primitives/segment.py | 39 ------------------- 7 files changed, 16 insertions(+), 62 deletions(-) delete mode 100644 src/cp_measure/primitives/segment.py diff --git a/src/cp_measure/core/measuregranularity.py b/src/cp_measure/core/measuregranularity.py index 0d7fc7e..07b7c10 100644 --- a/src/cp_measure/core/measuregranularity.py +++ b/src/cp_measure/core/measuregranularity.py @@ -248,8 +248,6 @@ def get_granularity( # Per-object stats use original-scale labels so the cell boundaries are exact. # start_mean uses the raw (non-background-subtracted) original pixels, matching # CellProfiler's ObjectRecord initialisation which also uses im_pixel_data directly. - unique_labels = numpy.unique(orig_mask) - unique_labels = unique_labels[unique_labels > 0] range_ = numpy.arange(1, numpy.max(orig_mask) + 1) current_mean = fix(scipy.ndimage.mean(orig_pixels, orig_mask, range_)) @@ -280,7 +278,7 @@ def get_granularity( # Calculate the means for the objects gss = numpy.zeros((0,)) - if unique_labels.any(): + if range_.size: new_mean = fix(scipy.ndimage.mean(rec_orig, orig_mask, range_)) gss = (current_mean - new_mean) * 100 / start_mean current_mean = new_mean # update running mean for next iteration diff --git a/src/cp_measure/core/measureobjectintensitydistribution.py b/src/cp_measure/core/measureobjectintensitydistribution.py index 9450250..7aad4cd 100644 --- a/src/cp_measure/core/measureobjectintensitydistribution.py +++ b/src/cp_measure/core/measureobjectintensitydistribution.py @@ -142,9 +142,8 @@ def get_radial_distribution( if labels.dtype == bool: labels = labels.astype(numpy.integer) - unique_labels = numpy.unique(labels) - unique_labels = unique_labels[unique_labels > 0] - nobjects = len(unique_labels) + nobjects = int(labels.max()) + unique_labels = numpy.arange(1, nobjects + 1) d_to_edge = centrosome.cpmorphology.distance_to_edge(labels) # Find the point in each object farthest away from the edge. @@ -318,8 +317,7 @@ def get_radial_zernikes( return {} zernike_indexes = centrosome.zernike.get_zernike_indexes(zernike_degree + 1) - unique_labels = numpy.unique(labels) # Will be used later for scipy.ndimage.sum - unique_labels = unique_labels[unique_labels > 0] + unique_labels = numpy.arange(1, int(labels.max()) + 1) # MODIFIED: Delegate index generation to the minimum_enclosing_circle # MODIFIED: We assume non-overlapping labels for now # TODO Support label overlap (i.e., format in ijv) diff --git a/src/cp_measure/core/measuretexture.py b/src/cp_measure/core/measuretexture.py index 769a88e..c5fa590 100644 --- a/src/cp_measure/core/measuretexture.py +++ b/src/cp_measure/core/measuretexture.py @@ -194,9 +194,6 @@ def get_texture( CellProfiler 3 versions it was fixed at 256. The minimum number of levels is 2, the maximum is 256. """ - unique_labels = numpy.unique(masks) - unique_labels = unique_labels[unique_labels > 0] - # Modified to use the number of dimensions in pixels to determine the number of directions n_directions = 13 if pixels.ndim > 2 else 4 @@ -213,7 +210,7 @@ def get_texture( ).astype(numpy.uint8) props = skimage.measure.regionprops(masks, pixels) - features = numpy.empty((n_directions, 13, len(unique_labels))) + features = numpy.empty((n_directions, 13, len(props))) for index, prop in enumerate(props): label_data = prop["intensity_image"] diff --git a/src/cp_measure/core/numba/measureobjectintensity.py b/src/cp_measure/core/numba/measureobjectintensity.py index 6ae43b8..558b51e 100644 --- a/src/cp_measure/core/numba/measureobjectintensity.py +++ b/src/cp_measure/core/numba/measureobjectintensity.py @@ -40,7 +40,6 @@ STD_INTENSITY_EDGE, UPPER_QUARTILE_INTENSITY, ) -from cp_measure.primitives.segment import label_to_idx_lut from cp_measure.primitives._segment_numba import ( flatten_numba, inner_boundary, @@ -74,7 +73,7 @@ def get_intensity( elif pixels.ndim == 3 and masks.ndim == 2: # 3D image, 2D mask masks = masks.reshape(1, *masks.shape) - lut, nobjects = label_to_idx_lut(masks) + nobjects = int(masks.max()) integrated_intensity = numpy.zeros(nobjects) mean_intensity = numpy.zeros(nobjects) @@ -96,7 +95,6 @@ def get_intensity( values, seg0, xc, yc, zc = flatten_numba( numpy.ascontiguousarray(masks), numpy.ascontiguousarray(masked_image), - lut, ) has_objects = values.size > 0 @@ -156,7 +154,7 @@ def get_intensity( else: emask = skimage.segmentation.find_boundaries(masks, mode="inner") > 0 e_values = masked_image[emask].astype(numpy.float64) - e_seg0 = lut[masks[emask]] + e_seg0 = masks[emask] - 1 if e_values.size > 0: ecount, esum, emin, emax = segment_stats(e_values, e_seg0, nobjects) diff --git a/src/cp_measure/primitives/__init__.py b/src/cp_measure/primitives/__init__.py index 3619525..79b534e 100644 --- a/src/cp_measure/primitives/__init__.py +++ b/src/cp_measure/primitives/__init__.py @@ -7,11 +7,11 @@ segment index + per-axis coordinates) which the segment kernels reduce. All spatial structure (2D vs 3D) and any future batch/image axis are encoded in that flat segment index, so a single set of segment kernels covers every case without -a rewrite. The numpy ``label_to_idx_lut`` (in ``segment``) builds the -label→index lookup; the flattening + reductions themselves live in -``_segment_numba``. +a rewrite. Labels are the contiguous ``1..N`` cp_measure guarantees (see +:mod:`cp_measure._sanitize`), so the segment index is ``label - 1``; the +flattening + reductions live in ``_segment_numba``. This is an internal layer with no curated public API: import its building -blocks directly from the concrete modules (``primitives.segment``, -``primitives._segment_numba``) rather than re-exporting them here. +blocks directly from the concrete modules (``primitives._segment_numba``) +rather than re-exporting them here. """ diff --git a/src/cp_measure/primitives/_segment_numba.py b/src/cp_measure/primitives/_segment_numba.py index a786d5f..dcc453a 100644 --- a/src/cp_measure/primitives/_segment_numba.py +++ b/src/cp_measure/primitives/_segment_numba.py @@ -14,12 +14,14 @@ @njit(cache=True) -def flatten_numba(masks, pixels, lut): +def flatten_numba(masks, pixels): """Flatten a labeled (Z, Y, X) image to ``(values, seg0, xc, yc, zc)``. Two grid scans (count, then fill) replace the numpy ``(masks>0)&isfinite`` mask + ``numpy.nonzero`` + fancy-index gathers; coordinates are the loop indices. Background and non-finite pixels are dropped, in C (raster) order. + Labels are the contiguous ``1..N`` cp_measure guarantees (see + :mod:`cp_measure._sanitize`), so the segment index is ``label - 1``. ``masks`` and ``pixels`` must be C-contiguous; ``pixels`` may be any float dtype (kept values are upcast into the float64 ``values`` output). """ @@ -46,7 +48,7 @@ def flatten_numba(masks, pixels, lut): if not np.isfinite(v): continue values[i] = v - seg0[i] = lut[L] + seg0[i] = L - 1 xc[i] = x yc[i] = y zc[i] = z diff --git a/src/cp_measure/primitives/segment.py b/src/cp_measure/primitives/segment.py deleted file mode 100644 index e5b5894..0000000 --- a/src/cp_measure/primitives/segment.py +++ /dev/null @@ -1,39 +0,0 @@ -"""Host-side segment helpers (numpy, backend-agnostic). - -A labeled image is reduced to a flat *segment* form — ``values[M]`` intensities, -``seg0[M]`` 0-based segment indices, and ``xc/yc/zc[M]`` per-pixel coordinates — -which the segment kernels then reduce. The flattening itself lives in the numba -layer (:func:`cp_measure.primitives._segment_numba.flatten_numba`); this module -holds only the numpy label→index lookup that feeds it. -""" - -import numpy -import scipy.ndimage -from numpy.typing import NDArray - - -def label_to_idx_lut( - masks: NDArray[numpy.integer], -) -> tuple[NDArray[numpy.int64], int]: - """Build a ``label -> 0..n-1`` lookup over the sorted positive labels. - - Returns ``(lut, n)`` where ``lut[label]`` is the segment index (and ``-1`` - for absent labels / background) and ``n`` is the object count. Output arrays - are indexed by this segment rank, which equals ``label - 1`` when labels are - the contiguous ``1..n`` that cp_measure expects. - - Present labels are enumerated with ``scipy.ndimage.find_objects`` (one pass) - rather than ``numpy.unique`` (a full-image sort): same ascending label set, - identical LUT, ~3-5x faster on large/3D masks. - """ - if not numpy.issubdtype(masks.dtype, numpy.integer): - masks = masks.astype(numpy.intp, copy=False) - bboxes = scipy.ndimage.find_objects(masks) - labels = numpy.array( - [i + 1 for i, sl in enumerate(bboxes) if sl is not None], dtype=numpy.int64 - ) - n = int(labels.size) - max_label = int(labels[-1]) if n else 0 - lut = numpy.full(max_label + 1, -1, dtype=numpy.int64) - lut[labels] = numpy.arange(n, dtype=numpy.int64) - return lut, n From 01279d416b1f0437f5205323c0e6af833e3aa419 Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 16:06:25 +0200 Subject: [PATCH 07/11] docs(primitives): fix stale segment.py cross-ref in numba kernel docstring flatten_numba (in this module) builds the flat arrays now; segment.py is gone. Co-Authored-By: Claude Opus 4.8 (1M context) --- src/cp_measure/primitives/_segment_numba.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/cp_measure/primitives/_segment_numba.py b/src/cp_measure/primitives/_segment_numba.py index dcc453a..b44249e 100644 --- a/src/cp_measure/primitives/_segment_numba.py +++ b/src/cp_measure/primitives/_segment_numba.py @@ -1,9 +1,9 @@ """Numba segment kernels (single-threaded, cached). These are the numba implementation of the segment-reduce / segment-quantile -primitives. They loop over the flat ``(values, seg0, coords)`` arrays produced -by :mod:`cp_measure.primitives.segment` — no image shape, no batch axis — so one -kernel set covers 2D, 3D, and (future) batched inputs unchanged. +primitives. They loop over the flat ``(values, seg0, coords)`` arrays that +:func:`flatten_numba` builds — no image shape, no batch axis — so one kernel set +covers 2D, 3D, and (future) batched inputs unchanged. All kernels are ``@njit(cache=True)`` and serial: no ``prange``/``nogil``. Parallelism is the job of the (future) batch layer over images, not the kernel. From 7ba1d0f799601ceafba23f64da5b23b7e3847aee Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Thu, 18 Jun 2026 16:47:03 +0200 Subject: [PATCH 08/11] docs(readme): non-contiguous labels now handled internally Sanitation relabels arbitrary IDs to 1..N; drop the relabel_sequential caveat. Co-Authored-By: Claude Opus 4.8 (1M context) --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index 470c25d..d15bfbd 100644 --- a/README.md +++ b/README.md @@ -189,7 +189,7 @@ measurecolocalization.get_correlation_overlap ### Important notes -- **Contiguous labels**: The input labels must be sequential (e.g., `[1,2,3]`, not `[1,3,4]`). You can use `skimage.segmentation.relabel_sequential` to ensure compliance. +- **Labels**: Any positive integer labels work — non-contiguous IDs (e.g. `[1, 3, 4]`) are relabelled to `1..N` internally without modifying your array, and results are reported against your original IDs. - **Fidelity**: If you need to match CellProfiler measurements 1:1, you must convert your image arrays to float values between 0 and 1. For instance, if you have an array of data type uint16, you must divide them all by 65535. This is important for radial distribution measurements. For the four intensity quantile features (`LowerQuartileIntensity`, `MedianIntensity`, `UpperQuartileIntensity`, `MADIntensity`) you additionally need `legacy=True` — see below. - **Speed**: The Granularity measurement is particularly slow (~80% of the compute time). Skip this one it if speed is of utmost importance. - **Legacy percentile convention**: `get_intensity` (numpy and numba backends), `get_core_measurements`, `get_core_measurements_3d`, and `make_featurizer_config` accept `legacy: bool = False`. The default uses `numpy.percentile` 'linear' (`(n-1)*q`) quartiles and the textbook `median(|x - median(x)|)` MAD. Pass `legacy=True` to reproduce the original cp_measure / CellProfiler behavior: `n*q` quartiles and the `(1/ndim)`-quantile MAD (which returns the 33rd percentile in 3D rather than the median). All other intensity features are identical either way. From 5046b5ea8e2f639338d6e9f236b483296e489b65 Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Mon, 22 Jun 2026 16:28:34 +0200 Subject: [PATCH 09/11] refactor(sanitize): sanitize only at entry points, not in raw funcs Address review feedback: raw measurement functions no longer sanitize their label argument. Drop the @sanitize_labels decorator from all core and correlation get_* functions and rename the helper sanitize_labels -> sanitize (the user-facing sanitize(fn) escape hatch), removing the wrapper._sanitized state attribute. Sanitation now lives only at the entry points: the bulk registries (get_core_measurements/_3d, get_correlation_measurements) gain a sanitize=True flag, and featurize fetches the raw functions (sanitize=False) so it keeps sanitizing each mask exactly once instead of decorating-then-unwrapping. This also removes the per-call overhead on fast functions when callers import them directly. Document the contiguous-1..N assumption in every raw function's docstring, update the README label note, and drop the redundant int(label) cast. Tests rewritten accordingly (no more _sanitized flag). Co-Authored-By: Claude Opus 4.8 (1M context) --- README.md | 2 +- src/cp_measure/_sanitize.py | 17 ++++-- src/cp_measure/bulk.py | 33 +++++++++--- src/cp_measure/core/measurecolocalization.py | 36 ++++++++++--- src/cp_measure/core/measuregranularity.py | 9 ++-- src/cp_measure/core/measureobjectintensity.py | 6 ++- .../measureobjectintensitydistribution.py | 14 +++-- src/cp_measure/core/measureobjectsizeshape.py | 27 ++++++---- src/cp_measure/core/measuretexture.py | 9 ++-- .../core/numba/measureobjectintensity.py | 6 ++- src/cp_measure/featurizer.py | 23 +++----- test/test_sanitize.py | 52 +++++++++---------- 12 files changed, 148 insertions(+), 86 deletions(-) diff --git a/README.md b/README.md index d15bfbd..d6a637a 100644 --- a/README.md +++ b/README.md @@ -189,7 +189,7 @@ measurecolocalization.get_correlation_overlap ### Important notes -- **Labels**: Any positive integer labels work — non-contiguous IDs (e.g. `[1, 3, 4]`) are relabelled to `1..N` internally without modifying your array, and results are reported against your original IDs. +- **Labels**: Any positive integer labels work — non-contiguous IDs (e.g. `[1, 3, 4]`) are relabelled to `1..N` without modifying your array, and results are reported against your original IDs. This happens at the entry points: `featurize` always sanitizes, and the bulk registries (`get_core_measurements`, `get_core_measurements_3d`, `get_correlation_measurements`) sanitize by default (pass `sanitize=False` to opt out). The raw measurement functions assume contiguous `1..N` labels; if you import one directly and pass gapped IDs, wrap it with `cp_measure._sanitize.sanitize` first. - **Fidelity**: If you need to match CellProfiler measurements 1:1, you must convert your image arrays to float values between 0 and 1. For instance, if you have an array of data type uint16, you must divide them all by 65535. This is important for radial distribution measurements. For the four intensity quantile features (`LowerQuartileIntensity`, `MedianIntensity`, `UpperQuartileIntensity`, `MADIntensity`) you additionally need `legacy=True` — see below. - **Speed**: The Granularity measurement is particularly slow (~80% of the compute time). Skip this one it if speed is of utmost importance. - **Legacy percentile convention**: `get_intensity` (numpy and numba backends), `get_core_measurements`, `get_core_measurements_3d`, and `make_featurizer_config` accept `legacy: bool = False`. The default uses `numpy.percentile` 'linear' (`(n-1)*q`) quartiles and the textbook `median(|x - median(x)|)` MAD. Pass `legacy=True` to reproduce the original cp_measure / CellProfiler behavior: `n*q` quartiles and the `(1/ndim)`-quantile MAD (which returns the 33rd percentile in 3D rather than the median). All other intensity features are identical either way. diff --git a/src/cp_measure/_sanitize.py b/src/cp_measure/_sanitize.py index 892ce51..f06ec26 100644 --- a/src/cp_measure/_sanitize.py +++ b/src/cp_measure/_sanitize.py @@ -45,10 +45,18 @@ def sanitize_masks(masks: NDArray) -> tuple[NDArray, NDArray[numpy.int64]]: return clean, ids -def sanitize_labels(func: Callable) -> Callable: - """Decorate a ``get_*`` function to sanitize its label argument (found by - name in :data:`_MASK_PARAMS`); functions with no such argument (e.g. the - two-mask multimask functions) are returned unchanged.""" +def sanitize(func: Callable) -> Callable: + """Wrap a ``get_*`` function so its label argument (found by name in + :data:`_MASK_PARAMS`) is relabelled to ``1..N`` before the call; functions + with no such argument (e.g. the two-mask multimask functions) are returned + unchanged. + + Measurement functions are *not* sanitized by default — the bulk entry + points (:func:`cp_measure.bulk.get_core_measurements` and friends) apply + this for you, and :func:`cp_measure.featurizer.featurize` sanitizes once up + front. Apply it yourself only when calling a raw function directly with + gapped or arbitrary label IDs. + """ sig = inspect.signature(func) param = next((name for name in _MASK_PARAMS if name in sig.parameters), None) if param is None: @@ -60,5 +68,4 @@ def wrapper(*args, **kwargs): bound.arguments[param], _ids = sanitize_masks(bound.arguments[param]) return func(*bound.args, **bound.kwargs) - wrapper._sanitized = True # type: ignore[attr-defined] return wrapper diff --git a/src/cp_measure/bulk.py b/src/cp_measure/bulk.py index 1456b83..d5444ce 100644 --- a/src/cp_measure/bulk.py +++ b/src/cp_measure/bulk.py @@ -105,25 +105,44 @@ def _apply_legacy(core: dict[str, Callable], legacy: bool) -> dict[str, Callable } -def get_core_measurements(legacy: bool = False) -> dict[str, Callable]: +def _maybe_sanitize(funcs: dict[str, Callable], sanitize: bool) -> dict[str, Callable]: + """Wrap each function to relabel arbitrary IDs to ``1..N``; no-op when False.""" + if not sanitize: + return funcs + from cp_measure._sanitize import sanitize as _sanitize + + return {name: _sanitize(fn) for name, fn in funcs.items()} + + +def get_core_measurements( + legacy: bool = False, sanitize: bool = True +) -> dict[str, Callable]: """Core per-object measurement functions. ``legacy`` (default False) selects the original CellProfiler ``n*q`` percentile convention for the intensity quartile/MAD features instead of the default ``numpy.percentile`` ``(n-1)*q``; see :func:`cp_measure.core.measureobjectintensity.get_intensity`. + + ``sanitize`` (default True) wraps each function so non-contiguous label IDs + (gaps or arbitrary values) are relabelled to ``1..N`` before the call. Set + False only when the caller has already guaranteed contiguous labels. """ - return _apply_legacy(_dispatch("core"), legacy) + core = _apply_legacy(_dispatch("core"), legacy) + return _maybe_sanitize(core, sanitize) -def get_core_measurements_3d(legacy: bool = False) -> dict[str, Callable]: - """Return only measurements that support 3D input (see ``legacy`` above).""" +def get_core_measurements_3d( + legacy: bool = False, sanitize: bool = True +) -> dict[str, Callable]: + """Return only measurements that support 3D input (see ``legacy``/``sanitize`` above).""" core = {k: _dispatch("core")[k] for k in _3D_FEATURES} - return _apply_legacy(core, legacy) + return _maybe_sanitize(_apply_legacy(core, legacy), sanitize) -def get_correlation_measurements() -> dict[str, Callable]: - return _dispatch("correlation") +def get_correlation_measurements(sanitize: bool = True) -> dict[str, Callable]: + """Correlation/colocalization measurement functions (see ``sanitize`` above).""" + return _maybe_sanitize(_dispatch("correlation"), sanitize) def get_multimask_measurements() -> dict[str, Callable]: diff --git a/src/cp_measure/core/measurecolocalization.py b/src/cp_measure/core/measurecolocalization.py index 993f5df..00b289e 100644 --- a/src/cp_measure/core/measurecolocalization.py +++ b/src/cp_measure/core/measurecolocalization.py @@ -86,7 +86,6 @@ import scipy.stats from scipy.linalg import lstsq -from cp_measure._sanitize import sanitize_labels M_IMAGES = "Across entire image" M_OBJECTS = "Within objects" @@ -546,12 +545,17 @@ def infer_scale(data: numpy.ndarray) -> int: # --------------------------------------------------------------------------- -@sanitize_labels def get_correlation_pearson( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], masks: NDArray[numpy.integer], ) -> dict[str, list[float]]: + """Per-object Pearson correlation and slope between two channels. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ corrs: list[float] = [] slopes: list[float] = [] for fi, si in _iter_label_pixels(pixels_1, pixels_2, masks): @@ -565,13 +569,18 @@ def get_correlation_pearson( return {F_CORRELATION_FORMAT: corrs, F_SLOPE_FORMAT: slopes} -@sanitize_labels def get_correlation_manders_fold( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], masks: NDArray[numpy.integer], thr: int = 15, ) -> dict[str, list[float]]: + """Per-object Manders fold coefficients between two channels. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ m1_list: list[float] = [] m2_list: list[float] = [] for fi, si in _iter_label_pixels(pixels_1, pixels_2, masks): @@ -585,13 +594,18 @@ def get_correlation_manders_fold( return {f"{F_MANDERS_FORMAT}_1": m1_list, f"{F_MANDERS_FORMAT}_2": m2_list} -@sanitize_labels def get_correlation_rwc( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], masks: NDArray[numpy.integer], thr: int = 15, ) -> dict[str, list[float]]: + """Per-object rank-weighted colocalization coefficients between two channels. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ r1: list[float] = [] r2: list[float] = [] for fi, si in _iter_label_pixels(pixels_1, pixels_2, masks): @@ -605,7 +619,6 @@ def get_correlation_rwc( return {f"{F_RWC_FORMAT}_1": r1, f"{F_RWC_FORMAT}_2": r2} -@sanitize_labels def get_correlation_costes( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], @@ -613,6 +626,12 @@ def get_correlation_costes( fast_costes: str = M_FASTER, thr: int = 15, ) -> dict[str, list[float]]: + """Per-object Costes colocalization coefficients between two channels. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ scale = infer_scale(pixels_1) c1: list[float] = [] c2: list[float] = [] @@ -627,13 +646,18 @@ def get_correlation_costes( return {f"{F_COSTES_FORMAT}_1": c1, f"{F_COSTES_FORMAT}_2": c2} -@sanitize_labels def get_correlation_overlap( pixels_1: NDArray[numpy.floating], pixels_2: NDArray[numpy.floating], masks: NDArray[numpy.integer], thr: int = 15, ) -> dict[str, list[float]]: + """Per-object overlap and k1/k2 colocalization coefficients between two channels. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ overlap_list: list[float] = [] k1_list: list[float] = [] k2_list: list[float] = [] diff --git a/src/cp_measure/core/measuregranularity.py b/src/cp_measure/core/measuregranularity.py index 07b7c10..235a9d6 100644 --- a/src/cp_measure/core/measuregranularity.py +++ b/src/cp_measure/core/measuregranularity.py @@ -54,12 +54,10 @@ import numpy import scipy.ndimage import skimage.morphology -from cp_measure._sanitize import sanitize_labels from cp_measure.utils import _ensure_np_array as fix from numpy.typing import NDArray -@sanitize_labels def get_granularity( mask: NDArray[numpy.integer], pixels: NDArray[numpy.floating], @@ -68,7 +66,12 @@ def get_granularity( element_size: int = 10, granular_spectrum_length: int = 16, ) -> dict[str, NDArray[numpy.floating]]: - """ + """Per-object granularity spectrum features. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + 1. (Outcommented) Subsample image 2. Remove background pixels using a greyscale tophat filter 3. Calculate granular spectrum (size distribution) for all masks diff --git a/src/cp_measure/core/measureobjectintensity.py b/src/cp_measure/core/measureobjectintensity.py index 1b425a2..6a3c0db 100644 --- a/src/cp_measure/core/measureobjectintensity.py +++ b/src/cp_measure/core/measureobjectintensity.py @@ -3,7 +3,6 @@ import skimage.segmentation from numpy.typing import NDArray -from cp_measure._sanitize import sanitize_labels __doc__ = """ MeasureObjectIntensity @@ -138,7 +137,6 @@ def _interp(sorted_arr: numpy.ndarray, n: int, frac: float, legacy: bool) -> flo return float(sorted_arr[lo] * (1.0 - f) + sorted_arr[hi] * f) -@sanitize_labels def get_intensity( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating], @@ -147,6 +145,10 @@ def get_intensity( ) -> dict[str, NDArray[numpy.floating]]: """Per-object intensity features. + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Walks each object on its ``scipy.ndimage.find_objects`` bounding box rather than the full image; for each label the per-pixel reductions, quantiles, MAD, max position, and centroid run on the bbox crop only. 2D inputs are promoted diff --git a/src/cp_measure/core/measureobjectintensitydistribution.py b/src/cp_measure/core/measureobjectintensitydistribution.py index 7aad4cd..57464c5 100644 --- a/src/cp_measure/core/measureobjectintensitydistribution.py +++ b/src/cp_measure/core/measureobjectintensitydistribution.py @@ -5,7 +5,6 @@ import numpy.ma import scipy.ndimage import scipy.sparse -from cp_measure._sanitize import sanitize_labels from cp_measure.utils import masks_to_ijv from numpy.typing import NDArray @@ -90,7 +89,6 @@ } -@sanitize_labels def get_radial_distribution( labels: NDArray[numpy.integer], pixels: NDArray[numpy.floating], @@ -101,6 +99,10 @@ def get_radial_distribution( """ Radial features (2D only) + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + zernike_degree : int Maximum zernike moment. @@ -306,13 +308,17 @@ def get_radial_distribution( return results -@sanitize_labels def get_radial_zernikes( labels: NDArray[numpy.integer], pixels: NDArray[numpy.floating], zernike_degree: int = 9, ) -> dict[str, NDArray[numpy.floating]]: - # Radial Zernike features (2D only) + """Per-object radial Zernike features (2D only). + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ if labels.ndim == 3: return {} zernike_indexes = centrosome.zernike.get_zernike_indexes(zernike_degree + 1) diff --git a/src/cp_measure/core/measureobjectsizeshape.py b/src/cp_measure/core/measureobjectsizeshape.py index be1521b..4842773 100644 --- a/src/cp_measure/core/measureobjectsizeshape.py +++ b/src/cp_measure/core/measureobjectsizeshape.py @@ -7,7 +7,6 @@ import numpy import scipy.ndimage import skimage.measure -from cp_measure._sanitize import sanitize_labels from cp_measure.utils import masks_to_ijv __doc__ = """\ @@ -567,7 +566,6 @@ """ -@sanitize_labels def get_sizeshape( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating] | None, @@ -575,7 +573,12 @@ def get_sizeshape( new_features: bool = True, spacing: Optional[tuple[float, ...]] = None, ) -> dict[str, NDArray[numpy.floating]]: - """Compute the measurements for multiple object masks.""" + """Compute the size/shape measurements for multiple object masks. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ # Properties available for both 2d and 3d desired_properties = [ "image", @@ -1005,15 +1008,17 @@ def get_sizeshape( return results -@sanitize_labels def get_zernike( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating] | None, zernike_numbers: int = 9, ) -> dict[str, NDArray[numpy.floating]]: - # - # Zernike features (2D only) - # + """Per-object Zernike shape features (2D only). + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ if masks.ndim == 3: return {} unique_indices = numpy.unique(masks) @@ -1030,11 +1035,15 @@ def get_zernike( return results -@sanitize_labels def get_feret( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating] | None ) -> dict[str, NDArray[numpy.floating]]: - # Feret diameter (2D only) + """Per-object Feret diameter features (2D only). + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + """ if masks.ndim == 3: return {} ijv = masks_to_ijv(masks) diff --git a/src/cp_measure/core/measuretexture.py b/src/cp_measure/core/measuretexture.py index c5fa590..1612b47 100644 --- a/src/cp_measure/core/measuretexture.py +++ b/src/cp_measure/core/measuretexture.py @@ -149,21 +149,24 @@ from numpy.typing import NDArray import skimage.util -from cp_measure._sanitize import sanitize_labels F_HARALICK = """AngularSecondMoment Contrast Correlation Variance InverseDifferenceMoment SumAverage SumVariance SumEntropy Entropy DifferenceVariance DifferenceEntropy InfoMeas1 InfoMeas2""".split() -@sanitize_labels def get_texture( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating], scale: int = 3, gray_levels: int = 256, ) -> dict[str, NDArray[numpy.floating]]: - """ + """Per-object Haralick texture features. + + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Parameters ---------- gray_levels : int, optional (default is 256) diff --git a/src/cp_measure/core/numba/measureobjectintensity.py b/src/cp_measure/core/numba/measureobjectintensity.py index 558b51e..08c01db 100644 --- a/src/cp_measure/core/numba/measureobjectintensity.py +++ b/src/cp_measure/core/numba/measureobjectintensity.py @@ -14,7 +14,6 @@ import skimage.segmentation from numpy.typing import NDArray -from cp_measure._sanitize import sanitize_labels from cp_measure.core.measureobjectintensity import ( C_LOCATION, INTEGRATED_INTENSITY, @@ -50,7 +49,6 @@ ) -@sanitize_labels def get_intensity( masks: NDArray[numpy.integer], pixels: NDArray[numpy.floating], @@ -59,6 +57,10 @@ def get_intensity( ) -> dict[str, NDArray[numpy.floating]]: """masks is a labeled array where 0 are background. + Assumes labels are the contiguous integers ``1..N``; call via a + :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with + :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + ``legacy`` mirrors the numpy backend: False (default) uses ``numpy.percentile`` 'linear' quartiles + textbook median MAD; True reproduces the original CellProfiler ``n*q`` quartiles + ``(1/ndim)``-quantile MAD. diff --git a/src/cp_measure/featurizer.py b/src/cp_measure/featurizer.py index d6b2dd1..c10bc20 100644 --- a/src/cp_measure/featurizer.py +++ b/src/cp_measure/featurizer.py @@ -15,7 +15,6 @@ from __future__ import annotations -import functools import itertools import warnings @@ -27,13 +26,6 @@ _2D_ONLY = {"radial_distribution", "radial_zernikes", "zernike", "feret"} -def _unwrap(fn): - """Strip the @sanitize_labels wrapper; the featurizer sanitizes once up front.""" - if isinstance(fn, functools.partial): - return functools.partial(_unwrap(fn.func), *fn.args, **fn.keywords) - return getattr(fn, "__wrapped__", fn) - - def make_featurizer_config( channels: list[str] | None = None, *, @@ -228,17 +220,14 @@ def featurize( is_3d = image.ndim == 4 legacy = config.get("legacy", False) + # Each mask is sanitized once in the loop below, so fetch the raw + # implementations (sanitize=False) rather than re-sanitizing per call. core_funcs = ( - get_core_measurements_3d(legacy=legacy) + get_core_measurements_3d(legacy=legacy, sanitize=False) if is_3d - else get_core_measurements(legacy=legacy) + else get_core_measurements(legacy=legacy, sanitize=False) ) - corr_funcs = get_correlation_measurements() - - # Each mask is sanitized once below, so call the underlying implementations - # directly instead of re-sanitizing in every feature's decorator. - core_funcs = {name: _unwrap(fn) for name, fn in core_funcs.items()} - corr_funcs = {name: _unwrap(fn) for name, fn in corr_funcs.items()} + corr_funcs = get_correlation_measurements(sanitize=False) if is_3d: _warn_2d_only_in_3d(config) @@ -299,7 +288,7 @@ def featurize( block = np.column_stack([results[c] for c in columns]) all_blocks.append(block) - all_rows.extend((image_id, object_name, int(label)) for label in ids) + all_rows.extend((image_id, object_name, label) for label in ids) if not all_blocks: raise ValueError("all masks have no labels (all zeros)") diff --git a/test/test_sanitize.py b/test/test_sanitize.py index cec495d..1c266bd 100644 --- a/test/test_sanitize.py +++ b/test/test_sanitize.py @@ -3,13 +3,8 @@ import numpy as np import pytest -from cp_measure._sanitize import sanitize_masks -from cp_measure.bulk import ( - get_core_measurements, - get_correlation_measurements, - get_multimask_measurements, -) -from cp_measure.core.measurecolocalization import get_correlation_overlap +from cp_measure._sanitize import sanitize, sanitize_masks +from cp_measure.bulk import get_core_measurements from cp_measure.featurizer import featurize, make_featurizer_config @@ -56,31 +51,34 @@ def test_invalid_input_raises(bad): sanitize_masks(bad) -def _is_sanitized(fn): - return bool( - getattr(fn, "_sanitized", False) - or getattr(getattr(fn, "func", None), "_sanitized", False) - ) +def test_entry_point_sanitizes_by_default(): + # get_zernike assumes 1..N; the bulk entry point sanitizes by default, so + # gapped IDs yield the same values as their contiguous relabelling. + px = np.random.default_rng(0).random((64, 64)) + zernike = get_core_measurements()["zernike"] + gapped = zernike(_three_objects((1, 17, 5)), px) + contig = zernike(_three_objects((1, 3, 2)), px) + for key in gapped: + np.testing.assert_allclose(gapped[key], contig[key], equal_nan=True) -def test_public_funcs_sanitized_multimask_excluded(): - public = ( - *get_core_measurements().values(), - *get_core_measurements(legacy=True).values(), # partial-wrapped intensity - *get_correlation_measurements().values(), - get_correlation_overlap, # public, not in the registry - ) - assert all(_is_sanitized(fn) for fn in public) - # two-mask functions are out of scope and must stay unwrapped - assert not any(_is_sanitized(fn) for fn in get_multimask_measurements().values()) +def test_raw_function_does_not_sanitize(): + # sanitize=False returns the bare implementation: it assumes 1..N and breaks + # on gapped IDs (here zernike indexes past the relabelled range). + px = np.random.default_rng(0).random((64, 64)) + raw = get_core_measurements(sanitize=False)["zernike"] + raw(_three_objects((1, 3, 2)), px) # contiguous: fine + with pytest.raises(IndexError): + raw(_three_objects((1, 17, 5)), px) -def test_decorated_func_handles_gapped_labels(): - # get_zernike raises on gapped labels without sanitation; the decorator fixes it. +def test_sanitize_helper_wraps_raw_function(): + # The user-facing escape hatch: wrap a raw function to handle gapped IDs. px = np.random.default_rng(0).random((64, 64)) - zernike = get_core_measurements()["zernike"] - gapped = zernike(_three_objects((1, 17, 5)), px) - contig = zernike(_three_objects((1, 3, 2)), px) + raw = get_core_measurements(sanitize=False)["zernike"] + wrapped = sanitize(raw) + gapped = wrapped(_three_objects((1, 17, 5)), px) + contig = raw(_three_objects((1, 3, 2)), px) for key in gapped: np.testing.assert_allclose(gapped[key], contig[key], equal_nan=True) From c44c44dba6fdbb567339439f4627bd908bc40e9d Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Mon, 22 Jun 2026 16:42:49 +0200 Subject: [PATCH 10/11] fix(typing): satisfy numpy 2.5 stubs in costes + radial distribution numpy 2.5 tightened/regressed its type stubs, turning two pre-existing lines red under mypy on Python 3.11+ (CI installs numpy unpinned, so this surfaced on the first commit after the release; it reproduces on the parent commit too): - bisection/linear costes returned numpy floats from a tuple[float, float] signature -> cast to float(). - numpy.ma.masked_array is now a non-callable TypeAliasType in the stubs; call the MaskedArray class directly (identical at runtime). Co-Authored-By: Claude Opus 4.8 (1M context) --- src/cp_measure/core/measurecolocalization.py | 4 ++-- src/cp_measure/core/measureobjectintensitydistribution.py | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/cp_measure/core/measurecolocalization.py b/src/cp_measure/core/measurecolocalization.py index 00b289e..8e2a09e 100644 --- a/src/cp_measure/core/measurecolocalization.py +++ b/src/cp_measure/core/measurecolocalization.py @@ -444,7 +444,7 @@ def linear_costes( i -= i_step except ValueError: break - return thr_fi_c, thr_si_c + return float(thr_fi_c), float(thr_si_c) def bisection_costes( @@ -518,7 +518,7 @@ def bisection_costes( thr_fi_c = (valid - 1) / scale_max thr_si_c = (a * thr_fi_c) + b - return thr_fi_c, thr_si_c + return float(thr_fi_c), float(thr_si_c) # MODIFIED: This reproduces the behaviour of the block at diff --git a/src/cp_measure/core/measureobjectintensitydistribution.py b/src/cp_measure/core/measureobjectintensitydistribution.py index 57464c5..6ad5d44 100644 --- a/src/cp_measure/core/measureobjectintensitydistribution.py +++ b/src/cp_measure/core/measureobjectintensitydistribution.py @@ -285,7 +285,7 @@ def get_radial_distribution( mask = pixel_count == 0 - radial_means: numpy.ma.MaskedArray = numpy.ma.masked_array( + radial_means: numpy.ma.MaskedArray = numpy.ma.MaskedArray( radial_values / pixel_count, mask ) From 90f2eed9db73cc7676d1f285c721267be5984af6 Mon Sep 17 00:00:00 2001 From: Tim Treis Date: Mon, 22 Jun 2026 17:14:47 +0200 Subject: [PATCH 11/11] refactor(sanitize): trim duplicated docstrings, tighten sanitizer Code-review cleanup (no behavior change on supported paths): - Collapse the 14 repeated 3-line "assumes 1..N" docstring stanzas to one line each; trim the over-long README bullet and the _sanitize module + wrapper docstrings and the featurizer mask-param doc (the policy is stated once centrally). Drop 3 stray blank-line additions and a lazy import. - sanitize_masks: raise on negative labels in both branches (the unique branch previously dropped them silently), and return an integer dtype for a bool fast-path mask instead of leaking bool downstream. - radial_distribution: astype(numpy.integer) -> int64 (the abstract type is deprecated as a dtype). Net -42 LOC. 92 tests pass; ruff + mypy clean. Co-Authored-By: Claude Opus 4.8 (1M context) --- README.md | 2 +- src/cp_measure/_sanitize.py | 32 +++++++------------ src/cp_measure/bulk.py | 3 +- src/cp_measure/core/measurecolocalization.py | 21 +++--------- src/cp_measure/core/measuregranularity.py | 4 +-- src/cp_measure/core/measureobjectintensity.py | 5 +-- .../measureobjectintensitydistribution.py | 10 ++---- src/cp_measure/core/measureobjectsizeshape.py | 12 ++----- src/cp_measure/core/measuretexture.py | 5 +-- .../core/numba/measureobjectintensity.py | 4 +-- src/cp_measure/featurizer.py | 10 +++--- 11 files changed, 33 insertions(+), 75 deletions(-) diff --git a/README.md b/README.md index d6a637a..b4212f7 100644 --- a/README.md +++ b/README.md @@ -189,7 +189,7 @@ measurecolocalization.get_correlation_overlap ### Important notes -- **Labels**: Any positive integer labels work — non-contiguous IDs (e.g. `[1, 3, 4]`) are relabelled to `1..N` without modifying your array, and results are reported against your original IDs. This happens at the entry points: `featurize` always sanitizes, and the bulk registries (`get_core_measurements`, `get_core_measurements_3d`, `get_correlation_measurements`) sanitize by default (pass `sanitize=False` to opt out). The raw measurement functions assume contiguous `1..N` labels; if you import one directly and pass gapped IDs, wrap it with `cp_measure._sanitize.sanitize` first. +- **Labels**: Any positive integer labels work — non-contiguous IDs (e.g. `[1, 3, 4]`) are relabelled to `1..N` internally without modifying your array, and results are reported against your original IDs. `featurize` and the bulk `get_*` registries sanitize by default (`sanitize=False` to opt out); raw measurement functions assume contiguous `1..N`, so wrap them with `cp_measure._sanitize.sanitize` if you call them directly with gapped IDs. - **Fidelity**: If you need to match CellProfiler measurements 1:1, you must convert your image arrays to float values between 0 and 1. For instance, if you have an array of data type uint16, you must divide them all by 65535. This is important for radial distribution measurements. For the four intensity quantile features (`LowerQuartileIntensity`, `MedianIntensity`, `UpperQuartileIntensity`, `MADIntensity`) you additionally need `legacy=True` — see below. - **Speed**: The Granularity measurement is particularly slow (~80% of the compute time). Skip this one it if speed is of utmost importance. - **Legacy percentile convention**: `get_intensity` (numpy and numba backends), `get_core_measurements`, `get_core_measurements_3d`, and `make_featurizer_config` accept `legacy: bool = False`. The default uses `numpy.percentile` 'linear' (`(n-1)*q`) quartiles and the textbook `median(|x - median(x)|)` MAD. Pass `legacy=True` to reproduce the original cp_measure / CellProfiler behavior: `n*q` quartiles and the `(1/ndim)`-quantile MAD (which returns the 33rd percentile in 3D rather than the median). All other intensity features are identical either way. diff --git a/src/cp_measure/_sanitize.py b/src/cp_measure/_sanitize.py index f06ec26..7706eb7 100644 --- a/src/cp_measure/_sanitize.py +++ b/src/cp_measure/_sanitize.py @@ -1,10 +1,6 @@ -"""Central sanitation of non-contiguous mask label IDs. - -cp_measure's measurement functions assume labels are the contiguous integers -``1..N`` (see :func:`cp_measure.featurizer.featurize`). Real segmentations may -use gaps (``{1, 5, 17}``) or arbitrary values; this maps them to ``1..N`` before -any math runs and reports results against the original IDs, without mutating the -caller's array. +"""Relabel non-contiguous mask label IDs (gaps or arbitrary values) to contiguous +``1..N`` without mutating the caller's array, applied at the entry points. See +:func:`sanitize`. """ import functools @@ -29,33 +25,29 @@ def sanitize_masks(masks: NDArray) -> tuple[NDArray, NDArray[numpy.int64]]: """ if masks.dtype != bool and not numpy.issubdtype(masks.dtype, numpy.integer): raise ValueError(f"labels must be an integer array, got dtype {masks.dtype!r}") + if masks.min(initial=0) < 0: + raise ValueError("labels must be non-negative") mx = int(masks.max(initial=0)) if mx == 0: return masks, numpy.empty(0, dtype=numpy.int64) - # Cheap contiguity check (~4x faster than numpy.unique); bincount needs a - # bounded range, so fall back to unique for pathologically large labels. + # bincount is faster than unique but needs a bounded range; fall back for huge labels. if mx <= masks.size: ids = numpy.flatnonzero(numpy.bincount(masks.ravel(), minlength=mx + 1)) else: ids = numpy.unique(masks) ids = ids[ids > 0].astype(numpy.int64) if ids.size == mx: - return masks, ids # already 1..N: no copy + # already 1..N: no copy (cast only a bool mask up to an integer dtype) + return (masks if masks.dtype != bool else masks.astype(numpy.intp)), ids clean, _forward, _inverse = relabel_sequential(masks) return clean, ids def sanitize(func: Callable) -> Callable: - """Wrap a ``get_*`` function so its label argument (found by name in - :data:`_MASK_PARAMS`) is relabelled to ``1..N`` before the call; functions - with no such argument (e.g. the two-mask multimask functions) are returned - unchanged. - - Measurement functions are *not* sanitized by default — the bulk entry - points (:func:`cp_measure.bulk.get_core_measurements` and friends) apply - this for you, and :func:`cp_measure.featurizer.featurize` sanitizes once up - front. Apply it yourself only when calling a raw function directly with - gapped or arbitrary label IDs. + """Wrap a ``get_*`` function to relabel its label argument (named in + :data:`_MASK_PARAMS`) to ``1..N`` before the call; functions with no such + argument are returned unchanged. Use this only to call a raw measurement + function directly with gapped IDs — the entry points already sanitize. """ sig = inspect.signature(func) param = next((name for name in _MASK_PARAMS if name in sig.parameters), None) diff --git a/src/cp_measure/bulk.py b/src/cp_measure/bulk.py index d5444ce..c59266b 100644 --- a/src/cp_measure/bulk.py +++ b/src/cp_measure/bulk.py @@ -6,6 +6,7 @@ from functools import partial from typing import Callable +from cp_measure._sanitize import sanitize as _sanitize from cp_measure.core import ( measurecolocalization, measuregranularity, @@ -109,8 +110,6 @@ def _maybe_sanitize(funcs: dict[str, Callable], sanitize: bool) -> dict[str, Cal """Wrap each function to relabel arbitrary IDs to ``1..N``; no-op when False.""" if not sanitize: return funcs - from cp_measure._sanitize import sanitize as _sanitize - return {name: _sanitize(fn) for name, fn in funcs.items()} diff --git a/src/cp_measure/core/measurecolocalization.py b/src/cp_measure/core/measurecolocalization.py index 8e2a09e..2bd435e 100644 --- a/src/cp_measure/core/measurecolocalization.py +++ b/src/cp_measure/core/measurecolocalization.py @@ -86,7 +86,6 @@ import scipy.stats from scipy.linalg import lstsq - M_IMAGES = "Across entire image" M_OBJECTS = "Within objects" M_IMAGES_AND_OBJECTS = "Both" @@ -552,9 +551,7 @@ def get_correlation_pearson( ) -> dict[str, list[float]]: """Per-object Pearson correlation and slope between two channels. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ corrs: list[float] = [] slopes: list[float] = [] @@ -577,9 +574,7 @@ def get_correlation_manders_fold( ) -> dict[str, list[float]]: """Per-object Manders fold coefficients between two channels. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ m1_list: list[float] = [] m2_list: list[float] = [] @@ -602,9 +597,7 @@ def get_correlation_rwc( ) -> dict[str, list[float]]: """Per-object rank-weighted colocalization coefficients between two channels. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ r1: list[float] = [] r2: list[float] = [] @@ -628,9 +621,7 @@ def get_correlation_costes( ) -> dict[str, list[float]]: """Per-object Costes colocalization coefficients between two channels. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ scale = infer_scale(pixels_1) c1: list[float] = [] @@ -654,9 +645,7 @@ def get_correlation_overlap( ) -> dict[str, list[float]]: """Per-object overlap and k1/k2 colocalization coefficients between two channels. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ overlap_list: list[float] = [] k1_list: list[float] = [] diff --git a/src/cp_measure/core/measuregranularity.py b/src/cp_measure/core/measuregranularity.py index 235a9d6..1b26e26 100644 --- a/src/cp_measure/core/measuregranularity.py +++ b/src/cp_measure/core/measuregranularity.py @@ -68,9 +68,7 @@ def get_granularity( ) -> dict[str, NDArray[numpy.floating]]: """Per-object granularity spectrum features. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). 1. (Outcommented) Subsample image 2. Remove background pixels using a greyscale tophat filter diff --git a/src/cp_measure/core/measureobjectintensity.py b/src/cp_measure/core/measureobjectintensity.py index 6a3c0db..9913cd1 100644 --- a/src/cp_measure/core/measureobjectintensity.py +++ b/src/cp_measure/core/measureobjectintensity.py @@ -3,7 +3,6 @@ import skimage.segmentation from numpy.typing import NDArray - __doc__ = """ MeasureObjectIntensity ====================== @@ -145,9 +144,7 @@ def get_intensity( ) -> dict[str, NDArray[numpy.floating]]: """Per-object intensity features. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). Walks each object on its ``scipy.ndimage.find_objects`` bounding box rather than the full image; for each label the per-pixel reductions, quantiles, MAD, diff --git a/src/cp_measure/core/measureobjectintensitydistribution.py b/src/cp_measure/core/measureobjectintensitydistribution.py index 6ad5d44..4c781c4 100644 --- a/src/cp_measure/core/measureobjectintensitydistribution.py +++ b/src/cp_measure/core/measureobjectintensitydistribution.py @@ -99,9 +99,7 @@ def get_radial_distribution( """ Radial features (2D only) - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). zernike_degree : int Maximum zernike moment. @@ -142,7 +140,7 @@ def get_radial_distribution( return {} if labels.dtype == bool: - labels = labels.astype(numpy.integer) + labels = labels.astype(numpy.int64) nobjects = int(labels.max()) unique_labels = numpy.arange(1, nobjects + 1) @@ -315,9 +313,7 @@ def get_radial_zernikes( ) -> dict[str, NDArray[numpy.floating]]: """Per-object radial Zernike features (2D only). - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ if labels.ndim == 3: return {} diff --git a/src/cp_measure/core/measureobjectsizeshape.py b/src/cp_measure/core/measureobjectsizeshape.py index 4842773..7179604 100644 --- a/src/cp_measure/core/measureobjectsizeshape.py +++ b/src/cp_measure/core/measureobjectsizeshape.py @@ -575,9 +575,7 @@ def get_sizeshape( ) -> dict[str, NDArray[numpy.floating]]: """Compute the size/shape measurements for multiple object masks. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ # Properties available for both 2d and 3d desired_properties = [ @@ -1015,9 +1013,7 @@ def get_zernike( ) -> dict[str, NDArray[numpy.floating]]: """Per-object Zernike shape features (2D only). - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ if masks.ndim == 3: return {} @@ -1040,9 +1036,7 @@ def get_feret( ) -> dict[str, NDArray[numpy.floating]]: """Per-object Feret diameter features (2D only). - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). """ if masks.ndim == 3: return {} diff --git a/src/cp_measure/core/measuretexture.py b/src/cp_measure/core/measuretexture.py index 1612b47..d54b37c 100644 --- a/src/cp_measure/core/measuretexture.py +++ b/src/cp_measure/core/measuretexture.py @@ -149,7 +149,6 @@ from numpy.typing import NDArray import skimage.util - F_HARALICK = """AngularSecondMoment Contrast Correlation Variance InverseDifferenceMoment SumAverage SumVariance SumEntropy Entropy DifferenceVariance DifferenceEntropy InfoMeas1 InfoMeas2""".split() @@ -163,9 +162,7 @@ def get_texture( ) -> dict[str, NDArray[numpy.floating]]: """Per-object Haralick texture features. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). Parameters ---------- diff --git a/src/cp_measure/core/numba/measureobjectintensity.py b/src/cp_measure/core/numba/measureobjectintensity.py index 08c01db..11fd08c 100644 --- a/src/cp_measure/core/numba/measureobjectintensity.py +++ b/src/cp_measure/core/numba/measureobjectintensity.py @@ -57,9 +57,7 @@ def get_intensity( ) -> dict[str, NDArray[numpy.floating]]: """masks is a labeled array where 0 are background. - Assumes labels are the contiguous integers ``1..N``; call via a - :mod:`cp_measure.bulk` ``get_*`` entry point or wrap with - :func:`cp_measure._sanitize.sanitize` to handle gapped IDs. + Labels must be contiguous ``1..N`` (see :func:`cp_measure._sanitize.sanitize`). ``legacy`` mirrors the numpy backend: False (default) uses ``numpy.percentile`` 'linear' quartiles + textbook median MAD; True reproduces the original diff --git a/src/cp_measure/featurizer.py b/src/cp_measure/featurizer.py index c10bc20..433770d 100644 --- a/src/cp_measure/featurizer.py +++ b/src/cp_measure/featurizer.py @@ -186,10 +186,9 @@ def featurize( masks : numpy.ndarray Integer-labeled masks with shape ``(M, H, W)`` or ``(M, Z, H, W)``. Must have the same ``ndim`` as *image*. - Background is 0; labels are positive integers. Non-contiguous or - arbitrary IDs (e.g. ``{1, 5, 17}``) are sanitized to ``1..N`` - internally (see :mod:`cp_measure._sanitize`); the input array is not - modified and output rows are reported against the original IDs. + Background is 0; any positive integer labels (non-contiguous IDs are + relabelled internally, array untouched, original IDs reported in the + rows; see :mod:`cp_measure._sanitize`). config : dict, optional Configuration dictionary produced by :func:`make_featurizer_config`. If ``None``, all features are enabled with default parameters. @@ -220,8 +219,7 @@ def featurize( is_3d = image.ndim == 4 legacy = config.get("legacy", False) - # Each mask is sanitized once in the loop below, so fetch the raw - # implementations (sanitize=False) rather than re-sanitizing per call. + # Sanitize each mask once in the loop below, so fetch raw (unsanitized) funcs. core_funcs = ( get_core_measurements_3d(legacy=legacy, sanitize=False) if is_3d