Mapping neuro-vascular unit communications reveals distinct angiogenic programs across developing mouse brain regions.

Code ↔ Paper

4 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 4 matches

[1] § Methods › Analysis of postnatal brain angiogenesis › Analysis of the vessel density per brain area ↔ liom_toolkit/segmentation/stats.py, lines 148–206 · score 0.61 · vessel density, vessel area, Brain regions, voxel, regional, masks
[2] § Methods › Analysis of postnatal brain angiogenesis › Vessel segmentation ↔ liom_toolkit/segmentation/plane_segmentation.py, lines 9–24 · score 0.61 · Scikit Image, LIOM Toolkit, segmentation module
[3] § Methods › Analysis of postnatal brain angiogenesis › Vessel segmentation ↔ liom_toolkit/segmentation/stats.py, lines 534–584 · score 0.53 · vessel density, scikit-image, Heatmap, creation, module, segmented
[4] § Methods › Spatial transcriptomics › Data processing ↔ liom_toolkit/utils/io.py, lines 260–384 · score 0.50 · nearest neighbor, location, filtered

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The paper is loaded by your browser, with JavaScript: read it at doi.org or on Europe PMC.

The authors' code

Python · 725 lines · 27 KB · GPL-3.0 · 2 matches

"""Per-region vessel morphometric statistics and Allen atlas region filtering."""
from __future__ import annotations
import math
import tempfile
from pathlib import Path
from typing import Any
import dask.array as da
import imageio.v3 as iio
import numpy as np
import PIL.Image
from dask.distributed import Future
from numpy.typing import NDArray
from tqdm.auto import tqdm
# pandas + scipy + scikit-image are moved into the [seg]/[stats] extras
# (D-01/D-05). The upfront ImportError here is the honest signal on an
# io-only install. pandas is shared between [stats] and [antspy]; scipy is
# shared between [seg] and [stats] -- the message names both so the user
# picks the extra matching their workflow. The `from e` chain preserves the
# underlying error for debugging (AGENTS §2).
try:
import pandas as pd
import scipy.ndimage as ndi
from scipy.ndimage import distance_transform_edt
from skimage import measure
from skimage.color import gray2rgb
from skimage.draw import circle_perimeter
from skimage.measure import label
from skimage.measure._regionprops import RegionProperties
from skimage.morphology import skeletonize
from skimage.util import img_as_ubyte
except ImportError as e:
raise ImportError(
"Please install liom-toolkit[seg] or [stats] to use the segmentation stats module."
) from e
from liom_toolkit.utils.dask_client import dask_client_manager
PIL.Image.MAX_IMAGE_PIXELS = 2_000_000_000 # finite DoS-guard limit (not None — AGENTS §2)
# Precomputed branching-point detection kernel + valid signature set (PERF-01b).
# The inherited get_branching_points ran 5 structural elements x 4 rotations =
# 20 ndi.binary_hit_or_miss calls. This single-kernel replacement encodes each
# of the 20 valid 3x3 branching patterns as a unique integer signature: the 8
# neighbor positions carry distinct bit weights and the signature is the sum of
# weights at foreground positions (the center carries weight 0 -- it is always
# foreground in a branching point). One ndi.convolve pass with the bit-weight
# kernel produces the per-pixel neighbor signature; a foreground pixel whose
# signature is in the valid set is a branching point. The result is array_equal
# to the old 20-convolution result (numerical-equivalence regression test).
_BRANCHING_KERNEL = np.array([[1, 2, 4], [8, 0, 16], [32, 64, 128]], dtype=np.int32)
def _build_branching_signatures() -> np.ndarray:
"""Build the sorted array of valid branching-point neighbor signatures.
Returns
-------
np.ndarray
Sorted 1-D ``int32`` array of the 20 valid neighbor signatures.
"""
base_selems = [
np.array([[0, 1, 0], [1, 1, 1], [0, 0, 0]]),
np.array([[1, 0, 1], [0, 1, 0], [1, 0, 0]]),
np.array([[1, 0, 1], [0, 1, 0], [0, 1, 0]]),
np.array([[0, 1, 0], [1, 1, 0], [0, 0, 1]]),
np.array([[0, 0, 1], [1, 1, 1], [0, 1, 0]]),
]
selems = [np.rot90(s, k=j) for s in base_selems for j in range(4)]
sigs = {int((s * _BRANCHING_KERNEL).sum()) for s in selems}
return np.array(sorted(sigs), dtype=np.int32)
_BRANCHING_SIGNATURES = _build_branching_signatures()
def compute_slice_metrics(
output_dir: str,
image: str,
mask: NDArray[np.number],
vessel_mask: NDArray[np.number],
region_map: NDArray[np.number],
vessel_exclude: NDArray[np.number],
voxel_size: float = 0.65,
) -> None:
"""Compute the metrics for a brain slice and save the results to disk.
``image`` is a string label (filename/identifier) used in the output
DataFrame and progress messages, not an image array. The actual image
arrays are ``mask``, ``vessel_mask``, ``region_map``, and
``vessel_exclude``.
Vessel-free regions: regions with no vessels yield a row with vessel
density = 0.0, vessel area = 0.0, and branching points = 0, but the
'mean diameter (um)' entry is OMITTED (the row has no such column value)
because the mean diameter of an empty vessel set is undefined. The
omitted diameter row is itself a publishable 'no vessels detected'
signal, not a silent gap.
Parameters
----------
output_dir : str
The directory to save the output to.
image : str
The label (filename/identifier) of the brain slice, used in the
output DataFrame and progress messages.
mask : ArrayLike
The mask of the tissue in the brain slice.
vessel_mask : ArrayLike
The mask of the vessels in the brain slice.
region_map : ArrayLike
The map of the regions in the brain slice.
vessel_exclude : ArrayLike
The mask of the vessels to exclude from the analysis.
voxel_size : float
The size of the voxels in the image.
Notes
-----
``ValueError`` propagates from :func:`calculate_regional_density` when a
region has zero area (bad region mask, caller error).
"""
# Setup output directory (overwrite-safe via the symlink-aware
# create_directory helper from utils.zarr_writer: a second call into an
# existing output_dir shutil.rmtree's the directory then recreates it,
# eliminating the FileExistsError race on re-run. Function-scope import
# avoids a circular import with utils.zarr_writer at module load time,
# matching the conversion.py:save_zarr pattern.)
from liom_toolkit.utils.zarr_writer import create_directory
create_directory(Path(output_dir), overwrite=True)
df = pd.DataFrame(
columns=[
"image",
"region",
"vessel area (um2)",
"tissue area (um2)",
"vessel density (um2/um2)",
"branching points",
"mean diameter (um)",
]
)
# Get the different brain regions
regions, region_count = label(region_map, return_num=True)
props_list = measure.regionprops(regions)
full_vessel_mask = vessel_mask * mask
full_vessel_mask = full_vessel_mask * vessel_exclude
# Compute metrics per region
for i in tqdm(
range(region_count), desc="Computing metrics per region for " + image, leave=False
):
region = get_vessel_region(regions, i, full_vessel_mask)
# Calculate vessel density
vessel_area, total_area, density = calculate_regional_density(
region, i, props_list, output_dir, voxel_size
)
# Count branching points
branching_points_count, skeleton, branching_points = get_branching_point_count(
region, output_dir, filename=str(i) + "_skeleton.tif"
)
draw_branch_point_circles(
skeleton, branching_points, output_dir, filename=str(i) + "_skeleton_circled.png"
)
# Calculate average diameter. Vessel-free regions raise ValueError
# from compute_average_diameter (mean diameter of an empty set is
# undefined); the diameter row is OMITTED for those regions while
# the density=0.0 row above is kept. The omitted diameter row is
# itself a publishable 'no vessels detected' signal, not a silent
# gap.
try:
mean_diameter = compute_average_diameter(region, skeleton, voxel_size)
# Save data
entry = pd.DataFrame.from_dict(
{
"image": [image],
"region": [i],
"vessel area (um2)": [vessel_area],
"tissue area (um2)": [total_area],
"vessel density (um2/um2)": [density],
"branching points": [branching_points_count],
"mean diameter (um)": [mean_diameter],
}
)
except ValueError:
# Vessel-free region: density=0.0 row kept, diameter row omitted.
entry = pd.DataFrame.from_dict(
{
"image": [image],
"region": [i],
"vessel area (um2)": [vessel_area],
"tissue area (um2)": [total_area],
"vessel density (um2/um2)": [density],
"branching points": [branching_points_count],
}
)
df = pd.concat([df, entry])
# Compute metrics for the whole slice. Use full_vessel_mask
# (vessel_mask * mask * vessel_exclude) for the density calculation so
# the 'total' row is consistent with the per-region rows, which all
# derive their vessel area from full_vessel_mask via get_vessel_region.
# The pre-fix code passed the raw vessel_mask, so the total row's
# vessel_area and vessel_density included vessels outside the tissue
# mask and explicitly excluded vessels -- making the total density
# higher than the sum of per-region densities (a data inconsistency
# in the same output DataFrame).
tissue_area, vessel_area, vessel_density = calculate_density(full_vessel_mask, mask, voxel_size)
# Use full_vessel_mask (vessel_mask * mask * vessel_exclude) for the
# whole-slice branching and diameter calculations so the 'total' row is
# consistent with the per-region rows, which all derive their vessel area
# from full_vessel_mask via get_vessel_region. The pre-fix code passed
# the raw vessel_mask, so the total row's branching count and diameter
# included vessels outside the tissue mask and explicitly excluded
# vessels -- inconsistent with the per-region rows in the same DataFrame.
branching_points_count, skeleton, branching_points = get_branching_point_count(
full_vessel_mask, output_dir
)
draw_branch_point_circles(skeleton, branching_points, output_dir)
# Whole-slice diameter: wrap in try/except ValueError mirroring the
# per-region row-omission pattern above. A vessel-free slice (no vessels
# anywhere) hits the empty-vessel-set ValueError from
# compute_average_diameter (D-01 contract); without this wrap the whole
# metrics computation crashes and the user loses every per-region row
# computed up to this point. The 'total' row keeps the density=0.0 /
# branching-points / vessel-area entries and OMITS the mean-diameter
# entry, the same publishable 'no vessels detected' signal used for
# vessel-free regions.
try:
mean_diameter = compute_average_diameter(full_vessel_mask, skeleton, voxel_size)
total_entry = {
"image": [image],
"region": "total",
"vessel area (um2)": [vessel_area],
"tissue area (um2)": [tissue_area],
"vessel density (um2/um2)": [vessel_density],
"branching points": [branching_points_count],
"mean diameter (um)": [mean_diameter],
}
except ValueError:
# Vessel-free slice: density=0.0 row kept, diameter row omitted.
total_entry = {
"image": [image],
"region": "total",
"vessel area (um2)": [vessel_area],
"tissue area (um2)": [tissue_area],
"vessel density (um2/um2)": [vessel_density],
"branching points": [branching_points_count],
}
# Save intermediate results. The mask/region/skeleton arrays are
# small-integer label/mask arrays, so an explicit .astype(np.uint8)
# narrowing is used instead of skimage.util.img_as_ubyte: img_as_ubyte
# on an integer input whose max fits in uint8 emits a "Downcasting ...
# without scaling" UserWarning (skimage telling us it skipped the
# rescale), and the explicit astype is identical output with no warning
# and clearer intent (we are narrowing a small-integer mask, not
# rescaling a float image).
iio.imwrite(str(Path(output_dir) / "regions.png"), regions.astype(np.uint8))
iio.imwrite(str(Path(output_dir) / "vessel_exclude.png"), vessel_exclude.astype(np.uint8))
iio.imwrite(str(Path(output_dir) / "_complete_mask.png"), mask.astype(np.uint8))
iio.imwrite(str(Path(output_dir) / "vessels.png"), vessel_mask.astype(np.uint8))
# Save data
entry = pd.DataFrame.from_dict(total_entry)
df = pd.concat([df, entry])
df.to_excel(str(Path(output_dir) / "regions.xlsx"), index=False)
def get_vessel_region(
regions: NDArray[np.number], region_index: int, vessel_mask: NDArray[np.number]
) -> NDArray[np.number]:
"""Get the vessels in a region.
Parameters
----------
regions : ArrayLike
The regions of the tissue mask.
region_index : int
The index of the region.
vessel_mask : ArrayLike
The mask of the vessels.
Returns
-------
NDArray[np.generic]
The vessel within the masked region.
"""
region = regions == region_index + 1
return region * vessel_mask
def calculate_regional_density(
region: NDArray[np.number],
region_index: int,
props_list: list[RegionProperties],
output_dir: str,
voxel_size: float = 0.65,
) -> tuple[float, float, float]:
"""Calculate the density of vessels in a region.
Parameters
----------
region : ArrayLike
The region to calculate the density of.
region_index : int
The computational index of the region.
props_list : list[RegionProperties]
The list of properties of the regions.
output_dir : str
The directory to save the region mask to.
voxel_size : float
The size of the voxels in the image.
Returns
-------
tuple[float, float, float]
The area of the vessels, the area of the region, and the density of
the vessels in a specific region.
Raises
------
ValueError
If the region has zero area (bad region mask, caller error).
"""
vessel_area = float((region == 1).sum()) * math.pow(voxel_size, 2)
total_area = float(props_list[region_index].area) * math.pow(voxel_size, 2)
if total_area == 0:
raise ValueError("Empty region: regionprops area is 0 (bad region mask, caller error)")
iio.imwrite(str(Path(output_dir) / f"{region_index}.tif"), region.astype(np.uint8))
density = vessel_area / total_area
return vessel_area, total_area, density
def calculate_density(
vessel_mask: NDArray[np.number], mask: NDArray[np.number], voxel_size: float = 0.65
) -> tuple[float, float, float]:
"""Calculate the areas of the tissue and vessel to compute vessel density in a mask.
Empty-result contract:
* Empty vessel mask over positive tissue (``mask.sum() > 0``) returns
``(tissue_area, 0.0, 0.0)`` -- 0 vessels / positive tissue area is a
well-defined 0 density (the math is defined).
* Empty tissue mask (``mask.sum() == 0``) raises ``ValueError`` -- an
empty tissue mask is a bad region mask and a caller error, not a
valid result.
Parameters
----------
vessel_mask : ArrayLike
The mask of the vessels.
mask : ArrayLike
The mask of the tissue.
voxel_size : float
The size of the voxels in the image.
Returns
-------
tuple[float, float, float]
The area of the tissue, the area of the vessels, and the density of
the vessels.
Raises
------
ValueError
If the tissue mask is empty (``mask.sum() == 0``).
"""
tissue_area = float(mask.sum()) * math.pow(voxel_size, 2)
if tissue_area == 0:
raise ValueError("Empty tissue mask: mask.sum() == 0 (bad region mask, caller error)")
vessel_area = float(vessel_mask.sum()) * math.pow(voxel_size, 2)
vessel_density = vessel_area / tissue_area
return tissue_area, vessel_area, vessel_density
def get_branching_point_count(
vessel_mask: NDArray[np.number], output_dir: str, filename: str = "skeleton.tif"
) -> tuple[int, NDArray[np.bool_], NDArray[np.bool_]]:
"""Get the number of branching points in a vessel mask.
Parameters
----------
vessel_mask : ArrayLike
The mask of the vessels.
output_dir : str
The directory to save the skeleton to.
filename : str
The filename to save the skeleton to.
Returns
-------
tuple[int, NDArray[np.bool_], NDArray[np.bool_]]
The number of branching points in the vessel mask, the skeleton of
the vessel mask, and the location of the branching points.
"""
skeleton = skeletonize(vessel_mask)
branching_points = get_branching_points(skeleton)
points_count = branching_points.sum()
iio.imwrite(str(Path(output_dir) / filename), skeleton.astype(np.uint8))
return points_count, skeleton, branching_points
def get_branching_points(skeleton: NDArray[np.bool_]) -> NDArray[np.bool_]:
"""Get the branching points in a skeleton using a single convolution pass.
Replaces the inherited 20-convolution implementation (5 structural
elements x 4 rotations = 20 ``ndi.binary_hit_or_miss`` calls OR-ed
together) with a single ``ndi.convolve`` pass over a bit-packed 3x3
neighbor-signature kernel, followed by a vectorized membership check
against the 20 valid branching-point signatures.
Each of the 20 valid 3x3 patterns is encoded as a unique integer
signature: the 8 neighbor positions carry distinct bit weights
(1, 2, 4, 8, 16, 32, 64, 128) and the signature is the sum of weights at
foreground positions (the center is always foreground in a branching
point, so it carries weight 0). A single convolution with the bit-weight
kernel produces the per-pixel neighbor signature; a foreground pixel
whose signature is in the valid set is a branching point. The result is
``array_equal`` to the old 20-convolution result (gated by a
numerical-equivalence regression test asserting ``array_equal``, not
``allclose`` -- branching-point detection is a boolean topology operation).
Source:
https://stackoverflow.com/questions/43037692/how-to-find-branch-point-from-binary-skeletonize-image
Parameters
----------
skeleton : ArrayLike
The skeleton of the vessels.
Returns
-------
NDArray[np.bool_]
The branching points in the skeleton.
"""
conv = ndi.convolve(skeleton.astype(np.int32), _BRANCHING_KERNEL, mode="constant", cval=0)
return np.isin(conv, _BRANCHING_SIGNATURES) & skeleton.astype(bool)
def draw_branch_point_circles(
skeleton: NDArray[np.bool_],
branching_points: NDArray[np.bool_],
output_dir: str,
filename: str = "skeleton_circled.png",
) -> None:
"""Draw circles around the branching points in a skeleton and save to disk.
Parameters
----------
skeleton : ArrayLike
The skeleton of the vessels.
branching_points : ArrayLike
The location of the branching points.
output_dir : str
The directory to save the skeleton to.
filename : str
The filename to save the skeleton to.
"""
circled_skeleton = gray2rgb(skeleton.astype(np.uint8))
points_to_draw = np.argwhere(branching_points)
for point in points_to_draw:
circy, circx = circle_perimeter(point[0], point[1], 7, shape=skeleton.shape)
circled_skeleton[circy, circx] = (220, 20, 20)
iio.imwrite(str(Path(output_dir) / filename), circled_skeleton)
del circled_skeleton
def compute_average_diameter(
mask: NDArray[np.number], skeleton: NDArray[np.bool_], voxel_size: float = 0.65
) -> float:
"""Compute the average diameter of the vessels in a mask.
Empty-result contract:
* Empty vessel set (no positive radii in the skeleton) raises
``ValueError`` -- the mean diameter of an empty set is undefined;
returning 0.0 would imply zero-width vessels exist (a
plausible-shaped-but-wrong value), and returning ``NaN`` would let a
silent NaN escape into the published pandas DataFrame. The raise sits
BEFORE ``np.mean`` so no NaN + RuntimeWarning can escape.
* Empty tissue mask raises ``ValueError`` -- mean diameter is undefined
when there is no tissue.
Parameters
----------
mask : ArrayLike
The vessel mask.
skeleton : ArrayLike
The skeleton of the vessels.
voxel_size : float
The size of the voxels in the image.
Returns
-------
float
The average diameter of the vessels in the mask.
Raises
------
ValueError
If the vessel mask is empty or no positive radii are found in the
skeleton.
TypeError
If the computed average diameter is not a float.
"""
if mask.sum() == 0:
raise ValueError("Empty tissue mask: mean diameter is undefined")
distance = distance_transform_edt(mask.astype(np.float64))
radii = distance * skeleton.astype(bool)
positive_radii = radii[radii > 0]
if positive_radii.size == 0:
raise ValueError(
"Empty vessel set: mean diameter is undefined (no positive radii in skeleton)"
)
mean_radius = np.mean(positive_radii)
mean_diameter = 2 * mean_radius
result = mean_diameter * voxel_size
if not isinstance(result, float):
raise TypeError(f"Expected float, got {type(result)}")
return result
def create_heatmap(image: NDArray[np.number], output_dir: str, square_size: int = 150) -> None:
"""Create and save a heatmap of the vessel density in a brain slice.
Parameters
----------
image : ArrayLike
The image of the brain slice.
output_dir : str
The directory to save the heatmap to.
square_size : int
The size of the squares in the heatmap.
"""
# Overwrite-safe output directory creation (same create_directory helper
# as compute_slice_metrics above; function-scope import avoids a circular
# import with utils.zarr_writer at module load time).
from liom_toolkit.utils.zarr_writer import create_directory
create_directory(Path(output_dir), overwrite=True)
image = img_as_ubyte(image)
image = image / 255
image = image.astype(np.uint8)
heatmap = np.zeros_like(image, dtype=np.uint32)
# Compute the per-dimension block counts separately so non-square images
# iterate the correct number of squares per dimension (a single
# shape[0]/square_size count for both dims produces a staircase pattern
# on non-square or multi-row images).
n_x = int(image.shape[0] / square_size)
n_y = int(image.shape[1] / square_size)
x_start = 0
for _ in range(n_x):
y_start = 0 # reset at the start of each outer iteration
for _j in range(n_y):
heatmap[x_start : x_start + square_size, y_start : y_start + square_size] = image[
x_start : x_start + square_size, y_start : y_start + square_size
].sum()
y_start += square_size
x_start += square_size
# Set final square to max value to ensure same scaling across heatmaps
heatmap[-1, -1] = square_size**2
# heatmap is uint32 with max square_size**2 (e.g. 22500), which fits in
# uint16 without scaling. Use an explicit .astype(np.uint16) narrowing
# instead of skimage.util.img_as_uint: img_as_uint on an integer input
# whose max fits in uint16 emits a "Downcasting ... without scaling"
# UserWarning, and the explicit astype is identical output with no
# warning and clearer intent.
heatmap = heatmap.astype(np.uint16)
heatmap = heatmap.astype(float)
heatmap = heatmap / (square_size**2)
iio.imwrite(str(Path(output_dir) / "heatmap.tif"), heatmap)
def generate_itk_id_list_of_region(region: str, data_dir: str = "") -> list[int]:
"""Generate a list of itk ids for a given region.
Reconstructs the structure tree and gets the descendants contained
within the region.
Parameters
----------
region : str
The region to get the ids for.
data_dir : str
The directory where the atlas and structure tree are saved. Optional.
Returns
-------
list[int]
The list of itk ids for the region and its descendants.
Raises
------
TypeError
If the extracted itk ids are not a list.
"""
# Setup temporary directory if not given. Track whether WE created it
# so the cleanup runs unconditionally (the pre-fix code reassigned
# data_dir = temp_dir.name then re-tested ``if data_dir == ""``, which
# was always False after the reassignment, so temp_dir.cleanup() never
# ran and the temp directory leaked on every call where data_dir="" --
# the default).
use_temp = data_dir == ""
temp_dir: tempfile.TemporaryDirectory[str] | None = None
if use_temp:
temp_dir = tempfile.TemporaryDirectory()
data_dir = temp_dir.name
itk_ids: list[Any] = []
try:
# Construct reference space and get itk ids
from liom_toolkit.utils import construct_reference_space
rs = construct_reference_space(data_dir)
structure_tree = rs.structure_tree
_, labels = rs.export_itksnap_labels()
# Get the itk ids for the region
region_structures = structure_tree.get_structures_by_name([region])
region_id = region_structures[0]["id"]
region_sub = structure_tree.descendant_ids([region_id])
region_sub_acronyms = [
region["acronym"] for region in structure_tree.get_structures_by_id(region_sub[0])
]
itk_ids = labels.loc[labels["LABEL"].isin(region_sub_acronyms)]["IDX"].to_numpy().tolist()
finally:
if use_temp and temp_dir is not None:
temp_dir.cleanup()
if not isinstance(itk_ids, list):
raise TypeError(f"Expected list, got {type(itk_ids)}")
return [int(x) for x in itk_ids]
def create_filter_image(atlas: da.Array | Future[Any], region_ids: list[int]) -> da.Array:
"""Create a filter image based on the region ids.
Parameters
----------
atlas : da.Array | Future[Any]
The atlas containing the region ids.
region_ids : list[int]
The region ids to filter.
Returns
-------
da.Array
The filter image.
Raises
------
TypeError
If the gathered filter image is not a Dask array.
"""
client = dask_client_manager.get_client()
filter_image = client.submit(da.isin, atlas, region_ids)
result = client.gather(filter_image)
if not isinstance(result, da.Array):
raise TypeError(f"Expected dask Array, got {type(result)}")
return result
def filter_image_to_region(image_filter: da.Array, data: da.Array | Future[Any]) -> da.Array:
"""Filter an image to a region based on a filter.
Parameters
----------
image_filter : da.Array
The filter to apply.
data : da.Array | Future[Any]
The data to filter.
Returns
-------
da.Array
The filtered image.
Raises
------
TypeError
If the gathered filtered image is not a Dask array.
"""
client = dask_client_manager.get_client()
filtered_image = client.submit(da.where, image_filter, data, 0)
result = client.gather(filtered_image)
if not isinstance(result, da.Array):
raise TypeError(f"Expected dask Array, got {type(result)}")
return result
def compute_mask_area(mask: da.Array | Future[Any]) -> np.uint64:
"""Compute the area of a mask by summing the binary mask values.
Parameters
----------
mask : da.Array | Future[Any]
The mask to compute the area of.
Returns
-------
np.uint64
The area of the mask.
"""
client = dask_client_manager.get_client()
total_area = client.submit(da.sum, mask)
total_area = client.gather(total_area)
# BOUNDARY-REQUIRED materialization: client.gather returns a 0-dimensional
# dask.array.Array (a scalar Dask array), NOT a Python scalar, so .compute()
# is required to materialize the scalar the function promises to return.
# Removing it would return a Dask array (wrong type) -- KEEP this .compute().
result = total_area.compute()
return np.uint64(result)

stats.py at commit aadfa71, under GPL-3.0 · at the source

Overview

Authors: Mathilde Bizou¹, Elise Drapé^1,2, Gael Cagnone¹, Joel P Howard¹, Frans Irgolitsch³, Séverine Leclerc¹, Mei Xi Chen¹, Blanche Boisseau¹, Isabelle Robillard^4,5, Matthieu Ruiz^4,5,6, Fréderic Lesage^3,4, Jean-Sébastien Joyal^1,2,7, Gregor Andelfinger^1,7, Alexandre Dubrac^1,2,8,9

ORCID iDs: Gael Cagnone, Frans Irgolitsch, Mei Xi Chen, Matthieu Ruiz, Fréderic Lesage, Jean-Sébastien Joyal, Alexandre Dubrac

Centre de Recherche, CHU Sainte Justine, Montréal, QC Canada
Département de Pharmacologie et de Physiologie, Université de Montréal, Montréal, QC Canada
Laboratoire d’Imagerie optique et Moléculaire, Polytechnique Montréal, Montréal, QC Canada
Centre de Recherche, Montréal Heart Institute, Montréal, QC Canada
Metabolomics platform, Montréal Heart Institute, Montréal, QC Canada
Département de nutrition, Université de Montréal, Montréal, QC Canada
Département de Pédiatrie, Université de Montréal, Montréal, QC Canada
Département de Pathologie et Biologie Cellulaire, Université de Montréal, Montréal, QC Canada
Département d’Ophtalmologie, Université de Montréal, Montréal, QC Canada

Institutions: Centre Hospitalier Universitaire Sainte-Justine (Canada); CHU Sainte-Justine Research Center (Canada); Université de Montréal (Canada); Polytechnique Montréal (Canada); Montreal Heart Institute (Canada)

Journal: Nature communications, volume 17, issue 1, article 6746

Dates: received 18 August 2023; accepted 28 April 2026; published online 22 May 2026

Type: Research article · Language: English

License: CC BY-NC-ND

Identifiers: DOI 10.1038/s41467-026-73373-w · PMID 42173924 · PMCID PMC13385802 · OpenAlex W7162147048

Open access: gold, a free copy (OpenAlex)

Status: code verified

Categories: genetics / omics (modality), mouse (organism)

Methods: Spectral & time-frequency, Statistics, Smoothing, state filtering, decompositions, Connectivity, Machine learning

Keywords: Angiogenesis, Neuro-vascular interactions

MeSH: Angiogenesis*, Brain*, Neovascularization, Physiologic*, Animals, Endothelial Cells, Female, Male, Mice, Mice, Inbred C57BL, Neurodevelopment, Neuroglia, Neurons, Signal Transduction, Spatial Transcriptomics, Thalamus, TOR Serine-Threonine Kinases (* major topic)

Topic: Axon Guidance and Neuronal Signaling (Cellular and Molecular Neuroscience, Neuroscience), according to OpenAlex

Funding: Gouvernement du Canada | Instituts de Recherche en Santé du Canada | CIHR Skin Research Training Centre (Skin Research Training Centre) (202203PJT-183658, 202403PJT-517269, 2019PJT-165871); Canadian Network for Research and Innovation in Machining Technology, Natural Sciences and Engineering Research Council of Canada (RGPIN-2022-04726)

Citations: cited by 1 paper (Europe PMC); 66 references in the paper

Abstract

The abstract is not reproduced here: the paper's license (CC BY-NC-ND) does not allow it. Read it in the paper, at the publisher or on Europe PMC.

Repository

Its files are read in the Code ↔ Paper reader above, with 4 matches between paragraphs and lines of code.

LIOMLab/liom-toolkit

License: GPL-3.0

State: the link answers, verified on 28 September 2026

Evidence: files inventoried

Commit: aadfa715e14ba44154053d8c1cd3cd15e0fa4c45, 28 August 2026

Languages: Python (88), Jupyter (6)

Size: 149 files, 94 scripts

Software Heritage: not archived

Found in: the text, “Vessel segmentation”

Holds: README, license file, environment (pyproject.toml, uv.lock), tests, continuous integration, documentation, 6 notebooks

Not found: CITATION.cff

Tools: NumPy (45 files), imageio (16 files), ANTs (11 files), PyTorch (11 files), scikit-image (10 files), pandas (6 files), Pillow (4 files), SciPy (4 files), h5py (3 files), NiBabel (2 files), tifffile (2 files), Matplotlib (1 file), OpenCV (1 file), scikit-learn (1 file), SimpleITK (1 file)

Availability: 1 check, the latest on 28 September 2026: the link answers

28 September 2026: the link answers

96 files

docs/source/conf.py, Python, 214 lines
docs/source/notebooks/save_masks.ipynb, Jupyter, 37 lines
docs/source/notebooks/segment_vessels.ipynb, Jupyter, 99 lines
docs/source/notebooks/templating.ipynb, Jupyter, 52 lines
docs/source/notebooks/vessel_density.ipynb, Jupyter, 63 lines
docs/source/notebooks/vseg_full.ipynb, Jupyter, 384 lines
docs/source/notebooks/zarr_conversion.ipynb, Jupyter, 58 lines
liom_toolkit/__init__.py, Python, 25 lines
liom_toolkit/_logging.py, Python, 80 lines
liom_toolkit/conversion/__init__.py, Python, 23 lines
liom_toolkit/conversion/conversion.py, Python, 552 lines
liom_toolkit/registration/__init__.py, Python, 19 lines
liom_toolkit/registration/register.py, Python, 539 lines
liom_toolkit/registration/templating.py, Python, 769 lines
liom_toolkit/scripts/__init__.py, Python, 7 lines
liom_toolkit/scripts/_common.py, Python, 56 lines
liom_toolkit/scripts/liom_align_annotations.py, Python, 112 lines
liom_toolkit/scripts/liom_build_template.py, Python, 118 lines
liom_toolkit/scripts/liom_compute_slice_metrics.py, Python, 92 lines
liom_toolkit/scripts/liom_convert_hdf5_to_zarr.py, Python, 84 lines
liom_toolkit/scripts/liom_create_mask.py, Python, 85 lines
liom_toolkit/scripts/liom_segment_2d.py, Python, 105 lines
liom_toolkit/scripts/liom_train_model.py, Python, 142 lines
liom_toolkit/segmentation/__init__.py, Python, 23 lines
liom_toolkit/segmentation/plane_segmentation.py, Python, 263 lines, 1 match
liom_toolkit/segmentation/stats.py, Python, 725 lines, 2 matches
liom_toolkit/segmentation/volume_segmentation.py, Python, 123 lines
liom_toolkit/segmentation/vseg/__init__.py, Python, 7 lines
liom_toolkit/segmentation/vseg/cldice.py, Python, 58 lines
liom_toolkit/segmentation/vseg/dataset.py, Python, 694 lines
liom_toolkit/segmentation/vseg/loss.py, Python, 82 lines
liom_toolkit/segmentation/vseg/model.py, Python, 214 lines
liom_toolkit/segmentation/vseg/prediction.py, Python, 293 lines
liom_toolkit/segmentation/vseg/training.py, Python, 602 lines
liom_toolkit/segmentation/vseg/utils.py, Python, 354 lines
liom_toolkit/segmentation/vseg/validation.py, Python, 158 lines
liom_toolkit/utils/__init__.py, Python, 61 lines
liom_toolkit/utils/allen_sdk.py, Python, 795 lines
liom_toolkit/utils/ants.py, Python, 109 lines
liom_toolkit/utils/checkpoint.py, Python, 454 lines
liom_toolkit/utils/dask_client.py, Python, 167 lines
liom_toolkit/utils/io.py, Python, 662 lines, 1 match
liom_toolkit/utils/utils.py, Python, 66 lines
liom_toolkit/utils/zarr_writer.py, Python, 750 lines
liom_toolkit/visualization/__init__.py, Python, 17 lines
liom_toolkit/visualization/slice_extraction.py, Python, 251 lines
tests/__init__.py, Python, 1 line
tests/conftest.py, Python, 173 lines
tests/test_api_surface.py, Python, 147 lines
tests/test_canary.py, Python, 28 lines
tests/test_conversion/__init__.py, Python, 1 line
tests/test_conversion/test_conversion.py, Python, 602 lines
tests/test_conversion/test_use_custom_atlas.py, Python, 194 lines
tests/test_docs/__init__.py, Python, 1 line
tests/test_docs/test_docs_structure.py, Python, 361 lines
tests/test_imports.py, Python, 237 lines
tests/test_package_metadata.py, Python, 361 lines
tests/test_pytest_config.py, Python, 320 lines
tests/test_registration/__init__.py, Python, 6 lines
tests/test_registration/test_register.py, Python, 444 lines
tests/test_registration/test_templating.py, Python, 648 lines
tests/test_release_config.py, Python, 528 lines
tests/test_scripts/__init__.py, Python, 1 line
tests/test_scripts/test_cli_common.py, Python, 108 lines
tests/test_scripts/test_liom_align_annotations.py, Python, 123 lines
tests/test_scripts/test_liom_build_template.py, Python, 99 lines
tests/test_scripts/test_liom_compute_slice_metrics.py, Python, 100 lines
tests/test_scripts/test_liom_convert_hdf5_to_zarr.py, Python, 95 lines
tests/test_scripts/test_liom_create_mask.py, Python, 70 lines
tests/test_scripts/test_liom_segment_2d.py, Python, 51 lines
tests/test_scripts/test_liom_train_model.py, Python, 116 lines
tests/test_segmentation/__init__.py, Python, 1 line
tests/test_segmentation/test_plane_segmentation.py, Python, 137 lines
tests/test_segmentation/test_stats.py, Python, 739 lines
tests/test_segmentation/test_volume_segmentation.py, Python, 102 lines
tests/test_segmentation/vseg/__init__.py, Python, 1 line
tests/test_segmentation/vseg/conftest.py, Python, 1 line
tests/test_segmentation/vseg/test_cldice.py, Python, 13 lines
tests/test_segmentation/vseg/test_dataset.py, Python, 598 lines
tests/test_segmentation/vseg/test_model.py, Python, 370 lines
tests/test_segmentation/vseg/test_prediction.py, Python, 149 lines
tests/test_segmentation/vseg/test_training.py, Python, 321 lines
tests/test_segmentation/vseg/test_utils.py, Python, 171 lines
tests/test_utils/__init__.py, Python, 1 line
tests/test_utils/fixtures/allen_itksnap_100um/regenerate.py, Python, 104 lines
tests/test_utils/test_allen_sdk.py, Python, 786 lines
tests/test_utils/test_checkpoint.py, Python, 363 lines
tests/test_utils/test_dask_client.py, Python, 216 lines
tests/test_utils/test_io.py, Python, 357 lines
tests/test_utils/test_logging.py, Python, 158 lines
tests/test_utils/test_utils.py, Python, 44 lines
tests/test_utils/test_zarr_writer.py, Python, 764 lines
tests/test_visualization/__init__.py, Python, 1 line
tests/test_visualization/test_slice_extraction.py, Python, 120 lines
LICENSE, License, 683 lines
README.md, Text, 190 lines

At the source: github.com/LIOMLab/liom-toolkit

Code availability statement

The paper has a code availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

it says that the code is available on request

Read it in the paper: doi.org/10.1038/s41467-026-73373-w.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

1 repository of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
94 scripts, each with its path and the digest of its content;
4 matches between paragraphs of the paper and lines of the code (method lexical-v1);
neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

Datasets cited

geo:GSE242215, at NCBI GEO; found in “Data availability”
urldefense.com/v3/__https:, at urldefense.com; found in “Data availability”

Code and data availability statement

The paper has a code and data availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

it points to 2 datasets: NCBI GEO GSE242215, urldefense.com/v3/__https:
it says that the data are available on request
it says that the code is available on request

Read it in the paper: doi.org/10.1038/s41467-026-73373-w.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 28 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 14 authors, 2 keywords, 16 MeSH terms, 2 funders, 66 references.

Cite

This paper

Bizou, M., Drapé, E., Cagnone, G., Howard, J. P., Irgolitsch, F., Leclerc, S., Chen, M. X., Boisseau, B., Robillard, I., Ruiz, M., Lesage, F., Joyal, J.-S., Andelfinger, G., & Dubrac, A. (2026). Mapping neuro-vascular unit communications reveals distinct angiogenic programs across developing mouse brain regions. Nature communications, 17(1), 6746. https://doi.org/10.1038/s41467-026-73373-w

BibTeX

@article{bizou2026mapping, author = {Bizou, Mathilde and Drapé, Elise and Cagnone, Gael and Howard, Joel P and Irgolitsch, Frans and Leclerc, Séverine and Chen, Mei Xi and Boisseau, Blanche and Robillard, Isabelle and Ruiz, Matthieu and Lesage, Fréderic and Joyal, Jean-Sébastien and Andelfinger, Gregor and Dubrac, Alexandre}, title = {{Mapping neuro-vascular unit communications reveals distinct angiogenic programs across developing mouse brain regions}}, journal = {Nature communications}, year = {2026}, month = may, volume = {17}, number = {1}, pages = {6746}, publisher = {Nature Publishing Group}, issn = {2041-1723}, doi = {10.1038/s41467-026-73373-w}, url = {https://doi.org/10.1038/s41467-026-73373-w}, pmid = {42173924}, pmcid = {PMC13385802} }

RIS

TY - JOUR AU - Bizou, Mathilde AU - Drapé, Elise AU - Cagnone, Gael AU - Howard, Joel P AU - Irgolitsch, Frans AU - Leclerc, Séverine AU - Chen, Mei Xi AU - Boisseau, Blanche AU - Robillard, Isabelle AU - Ruiz, Matthieu AU - Lesage, Fréderic AU - Joyal, Jean-Sébastien AU - Andelfinger, Gregor AU - Dubrac, Alexandre TI - Mapping neuro-vascular unit communications reveals distinct angiogenic programs across developing mouse brain regions T2 - Nature communications J2 - Nat Commun PY - 2026 DA - 2026/05/22 VL - 17 IS - 1 SP - 6746 SN - 2041-1723 PB - Nature Publishing Group DO - 10.1038/s41467-026-73373-w UR - https://doi.org/10.1038/s41467-026-73373-w LA - en ER -

CSL-JSON

{ "id": "10.1038/s41467-026-73373-w", "type": "article-journal", "title": "Mapping neuro-vascular unit communications reveals distinct angiogenic programs across developing mouse brain regions", "container-title": "Nature communications", "author": [ { "family": "Bizou", "given": "Mathilde" }, { "family": "Drapé", "given": "Elise" }, { "family": "Cagnone", "given": "Gael" }, { "family": "Howard", "given": "Joel P" }, { "family": "Irgolitsch", "given": "Frans" }, { "family": "Leclerc", "given": "Séverine" }, { "family": "Chen", "given": "Mei Xi" }, { "family": "Boisseau", "given": "Blanche" }, { "family": "Robillard", "given": "Isabelle" }, { "family": "Ruiz", "given": "Matthieu" }, { "family": "Lesage", "given": "Fréderic" }, { "family": "Joyal", "given": "Jean-Sébastien" }, { "family": "Andelfinger", "given": "Gregor" }, { "family": "Dubrac", "given": "Alexandre" } ], "container-title-short": "Nat Commun", "volume": "17", "issue": "1", "page": "6746", "DOI": "10.1038/s41467-026-73373-w", "PMID": "42173924", "PMCID": "PMC13385802", "ISSN": "2041-1723", "publisher": "Nature Publishing Group", "URL": "https://doi.org/10.1038/s41467-026-73373-w", "language": "en", "issued": { "date-parts": [ [ 2026, 5, 22 ] ] } }

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1364/boe.605322 [code]: Generalized plaque digitization framework for multi-dimensional mesoscopic images.
Journal: Biomedical optics express
In common: imageio, SimpleITK, tifffile, 10 other tools, mouse
[2] doi:10.1038/s41598-026-57519-w [code]: Automated segmentation of neurons and spinal cord structures in immunofluorescence images using SpineDL.
Journal: Scientific reports
In common: imageio, tifffile, OpenCV, 10 other tools, mouse
[3] doi:10.1016/j.xcrm.2026.102766 [code]: A longitudinal single-cell and spatial multiomic atlas of pediatric high-grade glioma.
Journal: Cell reports. Medicine
In common: imageio, tifffile, OpenCV, 8 other tools, genetics / omics, 2 references
[4] doi:10.3389/frai.2026.1771088 [code]: Few-shot deployment of pretrained MRI transformers in brain imaging tasks.
Journal: Frontiers in artificial intelligence
In common: imageio, SimpleITK, OpenCV, 10 other tools
[5] doi:10.1038/s41597-026-07248-6 [code]: A large-scale fMRI dataset for vision-language semantic association.
Journal: Scientific data
In common: imageio, ANTs, OpenCV, 10 other tools
[6] doi:10.3389/fnins.2026.1870124 [code]: An end-to-end pipeline for automated fetal brain segmentation and biometry from 3D SSFP MRI.
Journal: Frontiers in neuroscience
In common: imageio, SimpleITK, tifffile, 9 other tools
[7] doi:10.1162/imag.a.1326 [code]: RAVEN: Robust, generalizable, multi-resolution structural MRI upsampling using autoencoders.
Journal: Imaging neuroscience (Cambridge, Mass.)
In common: imageio, SimpleITK, OpenCV, 9 other tools
[8] doi:10.1002/alz.71649 [code]: Postmortem brain MRI reveals differential associations of subcortical and limbic volumes with cortical thinning and neurodegenerative pathologies.
Journal: Alzheimer's & dementia : the journal of the Alzheimer's Association
In common: SimpleITK, ANTs, OpenCV, 9 other tools
[9] doi:10.1038/s41586-026-10679-1 [code]: Cortical development dynamics across autism spectrum disorder mouse models.
Journal: Nature
In common: imageio, tifffile, OpenCV, 7 other tools, genetics / omics, mouse, 1 reference
[10] doi:10.1371/journal.pcbi.1014263 [code]: MIRAGE: Robust multi-modal architectures translate fMRI-to-image models from vision to mental imagery.
Journal: PLoS computational biology
In common: imageio, OpenCV, scikit-image, 9 other tools

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.