From 6b35422c0b98a4b3ccc2a44ba7f8b24f3af9ddb9 Mon Sep 17 00:00:00 2001 From: Yohsuke Fukai Date: Tue, 16 Jun 2026 21:46:09 -0700 Subject: [PATCH 1/4] fix? --- src/tracksdata/nodes/_regionprops.py | 47 +++++++++- .../nodes/_test/test_regionprops.py | 88 +++++++++++++++++++ 2 files changed, 133 insertions(+), 2 deletions(-) diff --git a/src/tracksdata/nodes/_regionprops.py b/src/tracksdata/nodes/_regionprops.py index 11ad7261..3cbb1ba4 100644 --- a/src/tracksdata/nodes/_regionprops.py +++ b/src/tracksdata/nodes/_regionprops.py @@ -37,6 +37,13 @@ class RegionPropsNodes(BaseNodesOperator): Physical spacing between pixels. If provided, affects distance-based measurements. Should be (row_spacing, col_spacing) for 2D or (depth_spacing, row_spacing, col_spacing) for 3D. + separate_arrays : bool, optional + If True, array-like properties (e.g. ``inertia_tensor`` or multi-channel + ``intensity_mean``) are flattened into multiple scalar attributes instead + of being stored as a single array attribute. The new attributes are named + ``_`` (e.g. ``intensity_mean_0``, ``intensity_mean_1``) using + the same convention as ``node_attrs(unpack=True)``. This makes individual + components filterable. Defaults to False. Attributes ---------- @@ -44,6 +51,8 @@ class RegionPropsNodes(BaseNodesOperator): List of additional properties to compute. _spacing : tuple[float, float] | None Physical spacing between pixels. + _separate_arrays : bool + Whether array-like properties are flattened into scalar attributes. Examples -------- @@ -92,6 +101,7 @@ def __init__( self, extra_properties: list[str | Callable[[RegionProperties], Any]] | None = None, spacing: tuple[float, float] | None = None, + separate_arrays: bool = False, ): super().__init__() self._extra_properties = extra_properties or [] @@ -102,6 +112,7 @@ def __init__( if "bbox" in self._extra_properties: raise ValueError("`bbox` is not supported as an extra property. It's already included by default.") self._spacing = spacing + self._separate_arrays = separate_arrays def _axis_names(self, labels: NDArray[np.integer]) -> list[str]: """ @@ -124,6 +135,37 @@ def _axis_names(self, labels: NDArray[np.integer]) -> list[str]: else: raise ValueError(f"`labels` must be 't + 2D' or 't + 3D', got '{labels.ndim}' dimensions.") + def _attr_items(self, key: str, value: Any) -> list[tuple[str, Any]]: + """ + Normalize a single property value into one or more node attribute items. + + Tuple/list/array-like numeric values are converted to numpy arrays so they + are stored consistently as fixed-shape array attributes. When + ``separate_arrays`` is enabled, such values are instead flattened into + scalar attributes named ``_`` (row-major), matching the + ``node_attrs(unpack=True)`` naming convention. + + Parameters + ---------- + key : str + The base attribute name. + value : Any + The property value returned by regionprops or a custom callable. + + Returns + ------- + list[tuple[str, Any]] + The (name, value) pairs to add to the node attributes. + """ + if isinstance(value, np.ndarray | tuple | list): + arr = np.asarray(value) + if arr.dtype.kind in "biufc" and arr.ndim >= 1: + if self._separate_arrays: + return [("_".join([key, *map(str, idx)]), arr[idx]) for idx in np.ndindex(arr.shape)] + return [(key, arr)] + + return [(key, value)] + def _init_node_attrs(self, graph: BaseGraph, node_attrs: dict[str, Any]) -> None: """ Initialize the node attributes for the graph. @@ -302,9 +344,10 @@ def _nodes_per_time( for prop in self._extra_properties: if callable(prop): - attrs[prop.__name__] = prop(obj) + key, value = prop.__name__, prop(obj) else: - attrs[prop] = getattr(obj, prop) + key, value = prop, getattr(obj, prop) + attrs.update(self._attr_items(key, value)) attrs[DEFAULT_ATTR_KEYS.MASK] = Mask(obj.image, obj.bbox) attrs[DEFAULT_ATTR_KEYS.BBOX] = np.asarray(obj.bbox, dtype=int) diff --git a/src/tracksdata/nodes/_test/test_regionprops.py b/src/tracksdata/nodes/_test/test_regionprops.py index 567c62e0..1fcca449 100644 --- a/src/tracksdata/nodes/_test/test_regionprops.py +++ b/src/tracksdata/nodes/_test/test_regionprops.py @@ -1,7 +1,9 @@ import numpy as np +import polars as pl import pytest from skimage.measure._regionprops import RegionProperties +from tracksdata.attrs import NodeAttr from tracksdata.constants import DEFAULT_ATTR_KEYS from tracksdata.graph import RustWorkXGraph from tracksdata.nodes import Mask, RegionPropsNodes @@ -334,3 +336,89 @@ def test_regionprops_multiprocessing_isolation() -> None: """Test that multiprocessing options don't affect subsequent tests.""" # Verify default n_workers is 1 assert get_options().n_workers == 1 + + +def test_regionprops_multichannel_intensity_array() -> None: + """Multi-channel intensity props are stored as fixed-shape array attributes (#195).""" + graph = RustWorkXGraph() + + labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) + intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) + intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] + intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] + + operator = RegionPropsNodes(extra_properties=["intensity_max"]) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + nodes_df = graph.node_attrs() + assert isinstance(nodes_df.schema["intensity_max"], pl.Array) + assert nodes_df["intensity_max"].dtype.shape == (2,) + + +def test_regionprops_tuple_property_stored_as_array() -> None: + """Tuple-returning props (e.g. centroid_weighted) are normalized to array attributes (#191).""" + graph = RustWorkXGraph() + + labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) + intensity = np.array([[[10, 20, 0], [30, 0, 40], [0, 50, 60]]], dtype=np.float32) + + operator = RegionPropsNodes(extra_properties=["centroid_weighted"]) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + nodes_df = graph.node_attrs() + # tuple props must become fixed-shape arrays (not pl.List) so they are unpackable + assert isinstance(nodes_df.schema["centroid_weighted"], pl.Array) + unpacked = graph.node_attrs(unpack=True) + assert "centroid_weighted_0" in unpacked.columns + assert "centroid_weighted_1" in unpacked.columns + + +def test_regionprops_separate_arrays() -> None: + """`separate_arrays=True` flattens array props into filterable scalar columns (#269).""" + graph = RustWorkXGraph() + + labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) + intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) + intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] + intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] + + operator = RegionPropsNodes(extra_properties=["intensity_max", "inertia_tensor"], separate_arrays=True) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + nodes_df = graph.node_attrs() + # 1D property -> single index suffix; 2D property -> row-major double index suffix + for col in ["intensity_max_0", "intensity_max_1", "inertia_tensor_0_0", "inertia_tensor_1_1"]: + assert col in nodes_df.columns + assert nodes_df[col].dtype == pl.Float64 + + # the array column itself must not exist when separated + assert "intensity_max" not in nodes_df.columns + + # separated columns are now filterable + subgraph = graph.filter(NodeAttr("intensity_max_0") > 30) + filtered = subgraph.node_attrs() + assert len(filtered) == 1 + assert filtered["intensity_max_0"][0] == 60.0 + + +def test_regionprops_separate_arrays_matches_unpack() -> None: + """`separate_arrays=True` column names match `node_attrs(unpack=True)`.""" + labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) + intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) + intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] + intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] + + extra = ["intensity_max", "inertia_tensor"] + + sep_graph = RustWorkXGraph() + RegionPropsNodes(extra_properties=extra, separate_arrays=True).add_nodes( + sep_graph, labels=labels, intensity_image=intensity + ) + + packed_graph = RustWorkXGraph() + RegionPropsNodes(extra_properties=extra).add_nodes(packed_graph, labels=labels, intensity_image=intensity) + + def _prop_cols(df: pl.DataFrame) -> set[str]: + return {c for c in df.columns if c.startswith(("intensity_max", "inertia_tensor"))} + + assert _prop_cols(sep_graph.node_attrs()) == _prop_cols(packed_graph.node_attrs(unpack=True)) From 3ee0f06ce88c382b23afb0da42b4d634e0f7adfa Mon Sep 17 00:00:00 2001 From: Yohsuke Fukai Date: Tue, 28 Jul 2026 13:35:51 +0900 Subject: [PATCH 2/4] refactor: drop separate_arrays from RegionPropsNodes `node_attrs(unpack=True)` already expands fixed-shape array columns into `_` scalars, so the operator-level flag was redundant for reading. The only capability it added was filtering on individual components, which belongs to #269 and is better served by struct-typed attributes (#268). Keeps the tuple/list -> ndarray normalization, which is the actual fix for #195 and #191. Co-Authored-By: Claude Opus 5 --- src/tracksdata/nodes/_regionprops.py | 19 +------ .../nodes/_test/test_regionprops.py | 53 ++++--------------- 2 files changed, 11 insertions(+), 61 deletions(-) diff --git a/src/tracksdata/nodes/_regionprops.py b/src/tracksdata/nodes/_regionprops.py index 3cbb1ba4..ffc97298 100644 --- a/src/tracksdata/nodes/_regionprops.py +++ b/src/tracksdata/nodes/_regionprops.py @@ -37,13 +37,6 @@ class RegionPropsNodes(BaseNodesOperator): Physical spacing between pixels. If provided, affects distance-based measurements. Should be (row_spacing, col_spacing) for 2D or (depth_spacing, row_spacing, col_spacing) for 3D. - separate_arrays : bool, optional - If True, array-like properties (e.g. ``inertia_tensor`` or multi-channel - ``intensity_mean``) are flattened into multiple scalar attributes instead - of being stored as a single array attribute. The new attributes are named - ``_`` (e.g. ``intensity_mean_0``, ``intensity_mean_1``) using - the same convention as ``node_attrs(unpack=True)``. This makes individual - components filterable. Defaults to False. Attributes ---------- @@ -51,8 +44,6 @@ class RegionPropsNodes(BaseNodesOperator): List of additional properties to compute. _spacing : tuple[float, float] | None Physical spacing between pixels. - _separate_arrays : bool - Whether array-like properties are flattened into scalar attributes. Examples -------- @@ -101,7 +92,6 @@ def __init__( self, extra_properties: list[str | Callable[[RegionProperties], Any]] | None = None, spacing: tuple[float, float] | None = None, - separate_arrays: bool = False, ): super().__init__() self._extra_properties = extra_properties or [] @@ -112,7 +102,6 @@ def __init__( if "bbox" in self._extra_properties: raise ValueError("`bbox` is not supported as an extra property. It's already included by default.") self._spacing = spacing - self._separate_arrays = separate_arrays def _axis_names(self, labels: NDArray[np.integer]) -> list[str]: """ @@ -140,10 +129,8 @@ def _attr_items(self, key: str, value: Any) -> list[tuple[str, Any]]: Normalize a single property value into one or more node attribute items. Tuple/list/array-like numeric values are converted to numpy arrays so they - are stored consistently as fixed-shape array attributes. When - ``separate_arrays`` is enabled, such values are instead flattened into - scalar attributes named ``_`` (row-major), matching the - ``node_attrs(unpack=True)`` naming convention. + are stored consistently as fixed-shape array attributes, which + ``node_attrs(unpack=True)`` can expand into ``_`` columns. Parameters ---------- @@ -160,8 +147,6 @@ def _attr_items(self, key: str, value: Any) -> list[tuple[str, Any]]: if isinstance(value, np.ndarray | tuple | list): arr = np.asarray(value) if arr.dtype.kind in "biufc" and arr.ndim >= 1: - if self._separate_arrays: - return [("_".join([key, *map(str, idx)]), arr[idx]) for idx in np.ndindex(arr.shape)] return [(key, arr)] return [(key, value)] diff --git a/src/tracksdata/nodes/_test/test_regionprops.py b/src/tracksdata/nodes/_test/test_regionprops.py index 1fcca449..0dc4eb01 100644 --- a/src/tracksdata/nodes/_test/test_regionprops.py +++ b/src/tracksdata/nodes/_test/test_regionprops.py @@ -3,7 +3,6 @@ import pytest from skimage.measure._regionprops import RegionProperties -from tracksdata.attrs import NodeAttr from tracksdata.constants import DEFAULT_ATTR_KEYS from tracksdata.graph import RustWorkXGraph from tracksdata.nodes import Mask, RegionPropsNodes @@ -373,52 +372,18 @@ def test_regionprops_tuple_property_stored_as_array() -> None: assert "centroid_weighted_1" in unpacked.columns -def test_regionprops_separate_arrays() -> None: - """`separate_arrays=True` flattens array props into filterable scalar columns (#269).""" +def test_regionprops_multidim_array_property_unpacks() -> None: + """2D array props (e.g. inertia_tensor) unpack into row-major scalar columns.""" graph = RustWorkXGraph() labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) - intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) - intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] - intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] - - operator = RegionPropsNodes(extra_properties=["intensity_max", "inertia_tensor"], separate_arrays=True) - operator.add_nodes(graph, labels=labels, intensity_image=intensity) - - nodes_df = graph.node_attrs() - # 1D property -> single index suffix; 2D property -> row-major double index suffix - for col in ["intensity_max_0", "intensity_max_1", "inertia_tensor_0_0", "inertia_tensor_1_1"]: - assert col in nodes_df.columns - assert nodes_df[col].dtype == pl.Float64 - - # the array column itself must not exist when separated - assert "intensity_max" not in nodes_df.columns - - # separated columns are now filterable - subgraph = graph.filter(NodeAttr("intensity_max_0") > 30) - filtered = subgraph.node_attrs() - assert len(filtered) == 1 - assert filtered["intensity_max_0"][0] == 60.0 - -def test_regionprops_separate_arrays_matches_unpack() -> None: - """`separate_arrays=True` column names match `node_attrs(unpack=True)`.""" - labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) - intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) - intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] - intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] - - extra = ["intensity_max", "inertia_tensor"] - - sep_graph = RustWorkXGraph() - RegionPropsNodes(extra_properties=extra, separate_arrays=True).add_nodes( - sep_graph, labels=labels, intensity_image=intensity - ) - - packed_graph = RustWorkXGraph() - RegionPropsNodes(extra_properties=extra).add_nodes(packed_graph, labels=labels, intensity_image=intensity) + operator = RegionPropsNodes(extra_properties=["inertia_tensor"]) + operator.add_nodes(graph, labels=labels) - def _prop_cols(df: pl.DataFrame) -> set[str]: - return {c for c in df.columns if c.startswith(("intensity_max", "inertia_tensor"))} + assert isinstance(graph.node_attrs().schema["inertia_tensor"], pl.Array) - assert _prop_cols(sep_graph.node_attrs()) == _prop_cols(packed_graph.node_attrs(unpack=True)) + unpacked = graph.node_attrs(unpack=True) + for col in ["inertia_tensor_0_0", "inertia_tensor_1_1"]: + assert col in unpacked.columns + assert unpacked[col].dtype == pl.Float64 From bd65af4ea060ad5bb69d6221bc93472b00910438 Mon Sep 17 00:00:00 2001 From: Yohsuke Fukai Date: Tue, 28 Jul 2026 16:56:48 +0900 Subject: [PATCH 3/4] revert --- src/tracksdata/nodes/_regionprops.py | 19 ++++++- .../nodes/_test/test_regionprops.py | 53 +++++++++++++++---- 2 files changed, 61 insertions(+), 11 deletions(-) diff --git a/src/tracksdata/nodes/_regionprops.py b/src/tracksdata/nodes/_regionprops.py index ffc97298..3cbb1ba4 100644 --- a/src/tracksdata/nodes/_regionprops.py +++ b/src/tracksdata/nodes/_regionprops.py @@ -37,6 +37,13 @@ class RegionPropsNodes(BaseNodesOperator): Physical spacing between pixels. If provided, affects distance-based measurements. Should be (row_spacing, col_spacing) for 2D or (depth_spacing, row_spacing, col_spacing) for 3D. + separate_arrays : bool, optional + If True, array-like properties (e.g. ``inertia_tensor`` or multi-channel + ``intensity_mean``) are flattened into multiple scalar attributes instead + of being stored as a single array attribute. The new attributes are named + ``_`` (e.g. ``intensity_mean_0``, ``intensity_mean_1``) using + the same convention as ``node_attrs(unpack=True)``. This makes individual + components filterable. Defaults to False. Attributes ---------- @@ -44,6 +51,8 @@ class RegionPropsNodes(BaseNodesOperator): List of additional properties to compute. _spacing : tuple[float, float] | None Physical spacing between pixels. + _separate_arrays : bool + Whether array-like properties are flattened into scalar attributes. Examples -------- @@ -92,6 +101,7 @@ def __init__( self, extra_properties: list[str | Callable[[RegionProperties], Any]] | None = None, spacing: tuple[float, float] | None = None, + separate_arrays: bool = False, ): super().__init__() self._extra_properties = extra_properties or [] @@ -102,6 +112,7 @@ def __init__( if "bbox" in self._extra_properties: raise ValueError("`bbox` is not supported as an extra property. It's already included by default.") self._spacing = spacing + self._separate_arrays = separate_arrays def _axis_names(self, labels: NDArray[np.integer]) -> list[str]: """ @@ -129,8 +140,10 @@ def _attr_items(self, key: str, value: Any) -> list[tuple[str, Any]]: Normalize a single property value into one or more node attribute items. Tuple/list/array-like numeric values are converted to numpy arrays so they - are stored consistently as fixed-shape array attributes, which - ``node_attrs(unpack=True)`` can expand into ``_`` columns. + are stored consistently as fixed-shape array attributes. When + ``separate_arrays`` is enabled, such values are instead flattened into + scalar attributes named ``_`` (row-major), matching the + ``node_attrs(unpack=True)`` naming convention. Parameters ---------- @@ -147,6 +160,8 @@ def _attr_items(self, key: str, value: Any) -> list[tuple[str, Any]]: if isinstance(value, np.ndarray | tuple | list): arr = np.asarray(value) if arr.dtype.kind in "biufc" and arr.ndim >= 1: + if self._separate_arrays: + return [("_".join([key, *map(str, idx)]), arr[idx]) for idx in np.ndindex(arr.shape)] return [(key, arr)] return [(key, value)] diff --git a/src/tracksdata/nodes/_test/test_regionprops.py b/src/tracksdata/nodes/_test/test_regionprops.py index 0dc4eb01..1fcca449 100644 --- a/src/tracksdata/nodes/_test/test_regionprops.py +++ b/src/tracksdata/nodes/_test/test_regionprops.py @@ -3,6 +3,7 @@ import pytest from skimage.measure._regionprops import RegionProperties +from tracksdata.attrs import NodeAttr from tracksdata.constants import DEFAULT_ATTR_KEYS from tracksdata.graph import RustWorkXGraph from tracksdata.nodes import Mask, RegionPropsNodes @@ -372,18 +373,52 @@ def test_regionprops_tuple_property_stored_as_array() -> None: assert "centroid_weighted_1" in unpacked.columns -def test_regionprops_multidim_array_property_unpacks() -> None: - """2D array props (e.g. inertia_tensor) unpack into row-major scalar columns.""" +def test_regionprops_separate_arrays() -> None: + """`separate_arrays=True` flattens array props into filterable scalar columns (#269).""" graph = RustWorkXGraph() labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) + intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) + intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] + intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] - operator = RegionPropsNodes(extra_properties=["inertia_tensor"]) - operator.add_nodes(graph, labels=labels) + operator = RegionPropsNodes(extra_properties=["intensity_max", "inertia_tensor"], separate_arrays=True) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) - assert isinstance(graph.node_attrs().schema["inertia_tensor"], pl.Array) + nodes_df = graph.node_attrs() + # 1D property -> single index suffix; 2D property -> row-major double index suffix + for col in ["intensity_max_0", "intensity_max_1", "inertia_tensor_0_0", "inertia_tensor_1_1"]: + assert col in nodes_df.columns + assert nodes_df[col].dtype == pl.Float64 - unpacked = graph.node_attrs(unpack=True) - for col in ["inertia_tensor_0_0", "inertia_tensor_1_1"]: - assert col in unpacked.columns - assert unpacked[col].dtype == pl.Float64 + # the array column itself must not exist when separated + assert "intensity_max" not in nodes_df.columns + + # separated columns are now filterable + subgraph = graph.filter(NodeAttr("intensity_max_0") > 30) + filtered = subgraph.node_attrs() + assert len(filtered) == 1 + assert filtered["intensity_max_0"][0] == 60.0 + + +def test_regionprops_separate_arrays_matches_unpack() -> None: + """`separate_arrays=True` column names match `node_attrs(unpack=True)`.""" + labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) + intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) + intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] + intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] + + extra = ["intensity_max", "inertia_tensor"] + + sep_graph = RustWorkXGraph() + RegionPropsNodes(extra_properties=extra, separate_arrays=True).add_nodes( + sep_graph, labels=labels, intensity_image=intensity + ) + + packed_graph = RustWorkXGraph() + RegionPropsNodes(extra_properties=extra).add_nodes(packed_graph, labels=labels, intensity_image=intensity) + + def _prop_cols(df: pl.DataFrame) -> set[str]: + return {c for c in df.columns if c.startswith(("intensity_max", "inertia_tensor"))} + + assert _prop_cols(sep_graph.node_attrs()) == _prop_cols(packed_graph.node_attrs(unpack=True)) From 1a7bd1211ee35046ebbe24e00c42603bae9642da Mon Sep 17 00:00:00 2001 From: Yohsuke Fukai Date: Thu, 30 Jul 2026 17:27:50 +0900 Subject: [PATCH 4/4] splitting channels --- src/tracksdata/nodes/_regionprops.py | 142 +++++++++++++++++- .../nodes/_test/test_regionprops.py | 127 ++++++++++++++++ 2 files changed, 262 insertions(+), 7 deletions(-) diff --git a/src/tracksdata/nodes/_regionprops.py b/src/tracksdata/nodes/_regionprops.py index 3cbb1ba4..201f764c 100644 --- a/src/tracksdata/nodes/_regionprops.py +++ b/src/tracksdata/nodes/_regionprops.py @@ -1,4 +1,4 @@ -from collections.abc import Callable +from collections.abc import Callable, Sequence from functools import partial from typing import Any @@ -6,7 +6,12 @@ import polars as pl from numpy.typing import NDArray from polars.datatypes import numpy_char_code_to_dtype -from skimage.measure._regionprops import RegionProperties, regionprops +from skimage.measure._regionprops import ( + PROPS, + RegionProperties, + _require_intensity_image, + regionprops, +) from typing_extensions import override from tracksdata.constants import DEFAULT_ATTR_KEYS @@ -44,6 +49,14 @@ class RegionPropsNodes(BaseNodesOperator): ``_`` (e.g. ``intensity_mean_0``, ``intensity_mean_1``) using the same convention as ``node_attrs(unpack=True)``. This makes individual components filterable. Defaults to False. + channel_names : Sequence[str] | None, optional + Names of the channels of a multi-channel `intensity_image`, whose channels + are expected on its last axis. When provided, intensity-dependent properties + (e.g. ``intensity_mean``, ``centroid_weighted``) are split along the channel + axis into one attribute per channel, named ``_`` + (e.g. ``intensity_mean_dapi``) instead of a single array attribute. + Custom callables in `extra_properties` are never split, since their output + layout is unknown. Defaults to None (no splitting). Attributes ---------- @@ -53,6 +66,8 @@ class RegionPropsNodes(BaseNodesOperator): Physical spacing between pixels. _separate_arrays : bool Whether array-like properties are flattened into scalar attributes. + _channel_names : list[str] | None + Names used to split intensity properties along the channel axis. Examples -------- @@ -95,6 +110,18 @@ def custom_property(region): labels_series = np.random.randint(0, 10, (10, 100, 100)) node_op.add_nodes(graph, labels=labels_series) ``` + + Name the channels of a multi-channel intensity image: + + ```python + node_op = RegionPropsNodes( + extra_properties=["intensity_mean"], + channel_names=["dapi", "gfp"], + ) + # intensity_image has shape (t, y, x, 2) + node_op.add_nodes(graph, labels=labels, intensity_image=intensity) + # nodes have `intensity_mean_dapi` and `intensity_mean_gfp` attributes + ``` """ def __init__( @@ -102,6 +129,7 @@ def __init__( extra_properties: list[str | Callable[[RegionProperties], Any]] | None = None, spacing: tuple[float, float] | None = None, separate_arrays: bool = False, + channel_names: Sequence[str] | None = None, ): super().__init__() self._extra_properties = extra_properties or [] @@ -113,6 +141,53 @@ def __init__( raise ValueError("`bbox` is not supported as an extra property. It's already included by default.") self._spacing = spacing self._separate_arrays = separate_arrays + self._channel_names = self._validate_channel_names(channel_names) + + @staticmethod + def _validate_channel_names(channel_names: Sequence[str] | None) -> list[str] | None: + """ + Validate and normalize the `channel_names` argument. + + Parameters + ---------- + channel_names : Sequence[str] | None + The channel names provided by the user. + + Returns + ------- + list[str] | None + The channel names as a list, or None when not provided. + + Raises + ------ + ValueError + If the names are empty, not strings, or not unique. + """ + if channel_names is None: + return None + + channel_names = list(channel_names) + if len(channel_names) == 0: + raise ValueError("`channel_names` must not be empty, use `None` for single-channel intensity images.") + + non_str = [name for name in channel_names if not isinstance(name, str)] + if non_str: + raise ValueError(f"`channel_names` must be strings, got {non_str}.") + + if len(set(channel_names)) != len(channel_names): + raise ValueError(f"`channel_names` must be unique, got {channel_names}.") + + return channel_names + + @staticmethod + def _is_intensity_prop(prop: str) -> bool: + """ + Whether a regionprops property name is computed from the intensity image. + + Such properties gain a trailing channel axis for multi-channel intensity + images, hence they are the ones split by `channel_names`. + """ + return PROPS.get(prop, prop) in _require_intensity_image def _axis_names(self, labels: NDArray[np.integer]) -> list[str]: """ @@ -135,7 +210,7 @@ def _axis_names(self, labels: NDArray[np.integer]) -> list[str]: else: raise ValueError(f"`labels` must be 't + 2D' or 't + 3D', got '{labels.ndim}' dimensions.") - def _attr_items(self, key: str, value: Any) -> list[tuple[str, Any]]: + def _attr_items(self, key: str, value: Any, split_channels: bool = False) -> list[tuple[str, Any]]: """ Normalize a single property value into one or more node attribute items. @@ -145,21 +220,45 @@ def _attr_items(self, key: str, value: Any) -> list[tuple[str, Any]]: scalar attributes named ``_`` (row-major), matching the ``node_attrs(unpack=True)`` naming convention. + When ``split_channels`` is set, the last axis is first split into one item + per channel, named ``_``; the remainder of each channel's + value is then normalized as above. + Parameters ---------- key : str The base attribute name. value : Any The property value returned by regionprops or a custom callable. + split_channels : bool + Whether the last axis of `value` is a channel axis to be split using + `channel_names`. Returns ------- list[tuple[str, Any]] The (name, value) pairs to add to the node attributes. + + Raises + ------ + ValueError + If the channel axis length does not match the number of channel names. """ if isinstance(value, np.ndarray | tuple | list): arr = np.asarray(value) if arr.dtype.kind in "biufc" and arr.ndim >= 1: + if split_channels: + if arr.shape[-1] != len(self._channel_names): + raise ValueError( + f"Property '{key}' has {arr.shape[-1]} channels, " + f"but {len(self._channel_names)} `channel_names` were provided: {self._channel_names}." + ) + items = [] + for i, name in enumerate(self._channel_names): + channel_value = arr[..., i] + # 0-dim arrays are stored as scalars, not as shape-() array attributes + items.extend(self._attr_items(f"{key}_{name}", channel_value[()])) + return items if self._separate_arrays: return [("_".join([key, *map(str, idx)]), arr[idx]) for idx in np.ndindex(arr.shape)] return [(key, arr)] @@ -192,6 +291,9 @@ def attr_keys(self) -> list[str]: Returns only the keys for extra_properties. The centroid coordinates (x, y, z) and mask are always included but not listed here. + When `channel_names` is set, intensity-dependent properties are listed + once per channel. Note that `separate_arrays` further splits array-valued + properties into per-index keys that are not listed here. Returns ------- @@ -206,7 +308,15 @@ def attr_keys(self) -> list[str]: print(keys) # ['area', 'perimeter'] ``` """ - return [prop.__name__ if callable(prop) else prop for prop in self._extra_properties] + keys = [] + for prop in self._extra_properties: + if callable(prop): + keys.append(prop.__name__) + elif self._channel_names is not None and self._is_intensity_prop(prop): + keys.extend(f"{prop}_{name}" for name in self._channel_names) + else: + keys.append(prop) + return keys @override def add_nodes( @@ -242,7 +352,14 @@ def add_nodes( intensity_image : NDArray | None, optional Intensity image(s) corresponding to the labels. Used for computing intensity-based properties. Must have the same shape as labels - (excluding the label values). + (excluding the label values), plus a trailing channel axis when + `channel_names` is used. + + Raises + ------ + ValueError + If `channel_names` was provided but `intensity_image` is missing or + does not have a matching trailing channel axis. Examples -------- @@ -271,6 +388,16 @@ def add_nodes( node_op.add_nodes(graph, labels=labels, t=0, intensity_image=fluorescence_image) ``` """ + if self._channel_names is not None: + if intensity_image is None: + raise ValueError("`channel_names` was provided but `intensity_image` is None.") + if intensity_image.ndim != labels.ndim + 1 or intensity_image.shape[-1] != len(self._channel_names): + raise ValueError( + f"`intensity_image` must have shape '{(*labels.shape, len(self._channel_names))}' to match " + f"`labels` plus the {len(self._channel_names)} channels of " + f"`channel_names` {self._channel_names}, got '{intensity_image.shape}'." + ) + if "shape" not in graph.metadata: graph.metadata.update(shape=labels.shape) @@ -344,10 +471,11 @@ def _nodes_per_time( for prop in self._extra_properties: if callable(prop): - key, value = prop.__name__, prop(obj) + key, value, split_channels = prop.__name__, prop(obj), False else: key, value = prop, getattr(obj, prop) - attrs.update(self._attr_items(key, value)) + split_channels = self._channel_names is not None and self._is_intensity_prop(prop) + attrs.update(self._attr_items(key, value, split_channels=split_channels)) attrs[DEFAULT_ATTR_KEYS.MASK] = Mask(obj.image, obj.bbox) attrs[DEFAULT_ATTR_KEYS.BBOX] = np.asarray(obj.bbox, dtype=int) diff --git a/src/tracksdata/nodes/_test/test_regionprops.py b/src/tracksdata/nodes/_test/test_regionprops.py index 1fcca449..31073c6b 100644 --- a/src/tracksdata/nodes/_test/test_regionprops.py +++ b/src/tracksdata/nodes/_test/test_regionprops.py @@ -401,6 +401,133 @@ def test_regionprops_separate_arrays() -> None: assert filtered["intensity_max_0"][0] == 60.0 +def _multichannel_data() -> tuple[np.ndarray, np.ndarray]: + labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32) + intensity = np.zeros((1, 3, 3, 2), dtype=np.float32) + intensity[..., 0] = [[10, 20, 0], [30, 0, 40], [0, 50, 60]] + intensity[..., 1] = [[1, 2, 0], [3, 0, 4], [0, 5, 6]] + return labels, intensity + + +def test_regionprops_channel_names() -> None: + """`channel_names` splits intensity props into one attribute per named channel.""" + graph = RustWorkXGraph() + labels, intensity = _multichannel_data() + + operator = RegionPropsNodes(extra_properties=["intensity_max", "area"], channel_names=["dapi", "gfp"]) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + nodes_df = graph.node_attrs() + # the packed array column is replaced by one scalar column per channel + assert "intensity_max" not in nodes_df.columns + assert nodes_df["intensity_max_dapi"].dtype == pl.Float64 + assert nodes_df["intensity_max_gfp"].dtype == pl.Float64 + # non-intensity properties are untouched + assert "area" in nodes_df.columns + + assert sorted(nodes_df["intensity_max_dapi"]) == [30.0, 60.0] + assert sorted(nodes_df["intensity_max_gfp"]) == [3.0, 6.0] + + # per-channel columns are filterable + filtered = graph.filter(NodeAttr("intensity_max_gfp") > 4).node_attrs() + assert len(filtered) == 1 + assert filtered["intensity_max_dapi"][0] == 60.0 + + +def test_regionprops_channel_names_multi_axis_property() -> None: + """Channel splitting keeps the remaining axes of multi-axis intensity props as arrays.""" + graph = RustWorkXGraph() + labels, intensity = _multichannel_data() + + operator = RegionPropsNodes(extra_properties=["centroid_weighted"], channel_names=["dapi", "gfp"]) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + nodes_df = graph.node_attrs() + for name in ["dapi", "gfp"]: + col = f"centroid_weighted_{name}" + assert isinstance(nodes_df.schema[col], pl.Array) + # (y, x) is kept as a fixed-shape array + assert nodes_df.schema[col].shape == (2,) + + +def test_regionprops_channel_names_with_separate_arrays() -> None: + """`separate_arrays` further splits the per-channel values into scalars.""" + graph = RustWorkXGraph() + labels, intensity = _multichannel_data() + + operator = RegionPropsNodes( + extra_properties=["intensity_max", "centroid_weighted"], + channel_names=["dapi", "gfp"], + separate_arrays=True, + ) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + nodes_df = graph.node_attrs() + for col in [ + "intensity_max_dapi", + "intensity_max_gfp", + "centroid_weighted_dapi_0", + "centroid_weighted_dapi_1", + "centroid_weighted_gfp_0", + "centroid_weighted_gfp_1", + ]: + assert col in nodes_df.columns + assert nodes_df[col].dtype == pl.Float64 + + +def test_regionprops_channel_names_custom_property_not_split() -> None: + """Custom callables keep their own layout, they are not split per channel.""" + graph = RustWorkXGraph() + labels, intensity = _multichannel_data() + + def max_per_channel(region: RegionProperties) -> np.ndarray: + return np.asarray(region.intensity_max) + + operator = RegionPropsNodes(extra_properties=[max_per_channel], channel_names=["dapi", "gfp"]) + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + nodes_df = graph.node_attrs() + assert isinstance(nodes_df.schema["max_per_channel"], pl.Array) + assert nodes_df.schema["max_per_channel"].shape == (2,) + + +def test_regionprops_channel_names_attr_keys() -> None: + """attr_keys lists intensity properties once per channel.""" + operator = RegionPropsNodes(extra_properties=["area", "intensity_mean"], channel_names=["dapi", "gfp"]) + assert operator.attr_keys() == ["area", "intensity_mean_dapi", "intensity_mean_gfp"] + + +@pytest.mark.parametrize( + "channel_names,match", + [ + ([], "must not be empty"), + (["dapi", "dapi"], "must be unique"), + (["dapi", 0], "must be strings"), + ], +) +def test_regionprops_channel_names_invalid(channel_names: list, match: str) -> None: + """Invalid `channel_names` are rejected at construction time.""" + with pytest.raises(ValueError, match=match): + RegionPropsNodes(channel_names=channel_names) + + +def test_regionprops_channel_names_shape_mismatch() -> None: + """The intensity image must have a trailing axis matching `channel_names`.""" + graph = RustWorkXGraph() + labels, intensity = _multichannel_data() + + operator = RegionPropsNodes(extra_properties=["intensity_max"], channel_names=["dapi", "gfp", "rfp"]) + with pytest.raises(ValueError, match="must have shape"): + operator.add_nodes(graph, labels=labels, intensity_image=intensity) + + # single-channel intensity image, no trailing channel axis + with pytest.raises(ValueError, match="must have shape"): + operator.add_nodes(graph, labels=labels, intensity_image=intensity[..., 0]) + + with pytest.raises(ValueError, match="`intensity_image` is None"): + operator.add_nodes(graph, labels=labels) + + def test_regionprops_separate_arrays_matches_unpack() -> None: """`separate_arrays=True` column names match `node_attrs(unpack=True)`.""" labels = np.array([[[1, 1, 0], [1, 0, 2], [0, 2, 2]]], dtype=np.int32)