Image segmentation¶

Segment an image one field of view at a time.

Segment an OME-Zarr image with ngio and skimage, one field of view at a time, and write the result back as a label. The second half repeats the segmentation inside a mask, so it only runs where you want it to.

Step 1: set up¶

Start with a function that segments an image, using skimage to do the work.

# Setup a simple segmentation function
import numpy as np
import skimage


def otsu_threshold_segmentation(image: np.ndarray, max_label: int) -> np.ndarray:
    """Simple segmentation using Otsu thresholding."""
    threshold = skimage.filters.threshold_otsu(image)
    binary = image > threshold
    label_image = skimage.measure.label(binary)
    label_image += max_label
    label_image = np.where(binary, label_image, 0)
    return label_image.astype(np.uint32)

Step 2: open the OME-Zarr container¶

from pathlib import Path

from ngio import open_ome_zarr_container
from ngio.utils import download_ome_zarr_dataset

# Download the dataset
download_dir = Path("./data").absolute()
hcs_path = download_ome_zarr_dataset("CardiomyocyteTiny", download_dir=download_dir)
image_path = hcs_path / "B" / "03" / "0"

# Open the OME-Zarr container
ome_zarr = open_ome_zarr_container(image_path)

Step 3: segment the image¶

Rather than segmenting the image all at once, iterate over its FOVs and segment them one by one.

from ngio.iterators import SegmentationIterator

# Take the image to read from, and the FOV table naming the regions to walk
image = ome_zarr.get_image()
roi_table = ome_zarr.get_roi_table("FOV_ROI_table")

# Derive an empty label image to write the segmentation into
label = ome_zarr.derive_label("new_label", overwrite=True)

# Setup the segmentation iterator
seg_iterator = SegmentationIterator(
    input_image=image,
    output_label=label,
    channel_selection="DAPI",
    axes_order=["z", "y", "x"],
)
seg_iterator = seg_iterator.product(roi_table)

# Split any remaining time axis, so each step yields one whole ZYX volume
seg_iterator = seg_iterator.by_zyx()

max_label = 0  # Carried across regions so the label ids never collide
for image_data, label_writer in seg_iterator.iter_as_numpy():
    roi_segmentation = otsu_threshold_segmentation(
        image_data, max_label
    )  # Segment the image

    max_label = roi_segmentation.max()  # Get the max label for the next iteration

    label_writer(patch=roi_segmentation)  # Write the segmentation back to the label

# No need to consolidate, the iterator does it automatically after the last write

Plot the segmentation¶

rand_cmap = random_label_cmap()
original = image.get_as_numpy(c=0, z=1, axes_order=["y", "x"])

# The data does not fill its uint16 range, so window it on percentiles.
vmin, vmax = np.percentile(original, (1, 99.8))

fig, axs = plt.subplots(2, 1, figsize=(8, 6))
axs[0].set_title("Original image")
axs[0].imshow(original, cmap="gray", vmin=vmin, vmax=vmax)
axs[1].set_title("Final segmentation")
axs[1].imshow(
    label.get_as_numpy(z=1, axes_order=["y", "x"]),
    cmap=rand_cmap,
    interpolation="nearest",
)
for ax in axs:
    ax.axis("off")
fig.tight_layout()
print(figure_html(fig))

Step 4: masked image segmentation¶

Now use a mask to restrict the segmentation to certain areas of the image. Here you create the mask by hand for illustration, but in a real pipeline it would usually come from another segmentation.

# Create a basic mask for illustration purposes
mask = ome_zarr.derive_label("mask", overwrite=True)
mask_data = mask.get_as_numpy(axes_order=["z", "y", "x"])
mask_data[:, 200:-200, 500:2000] = 1
mask_data[:, 200:-200, 3000:-500] = 2
mask_data[:, 600:-600, 1200:-1000] = 0
mask_data[:, 700:-700, 1600:-1500] = 3
mask.set_array(mask_data, axes_order=["z", "y", "x"])
mask.consolidate()

fig, axs = plt.subplots(2, 1, figsize=(8, 6))
axs[0].set_title("Original image")
axs[0].imshow(original, cmap="gray", vmin=vmin, vmax=vmax)
axs[1].set_title("Mask")
axs[1].imshow(
    mask.get_as_numpy(z=1, axes_order=["y", "x"]),
    cmap=rand_cmap,
    interpolation="nearest",
)
for ax in axs:
    ax.axis("off")
fig.tight_layout()
print(figure_html(fig))

Note that the next step rebinds image to the masked image, so the plot below shows the masked image rather than the original one.

from ngio.iterators import MaskedSegmentationIterator

# Take a masked image, which carries its masking ROI table with it
image = ome_zarr.get_masked_image(masking_label_name="mask")

# Derive an empty label image to write the segmentation into
label = ome_zarr.derive_label("masked_new_label", overwrite=True)

# Setup the masked segmentation iterator
seg_iterator = MaskedSegmentationIterator(
    input_image=image,
    output_label=label,
    channel_selection="DAPI",
    axes_order=["z", "y", "x"],
)

# Split any remaining time axis, so each step yields one whole ZYX volume
seg_iterator = seg_iterator.by_zyx()

max_label = 0  # Carried across regions so the label ids never collide
for image_data, label_writer in seg_iterator.iter_as_numpy():
    roi_segmentation = otsu_threshold_segmentation(
        image_data, max_label
    )  # Segment the image

    max_label = roi_segmentation.max()  # Get the max label for the next iteration

    label_writer(patch=roi_segmentation)  # Write the segmentation back to the label

# No need to consolidate, the iterator does it automatically after the last write

fig, axs = plt.subplots(2, 1, figsize=(8, 6))
axs[0].set_title("Original image")
axs[0].imshow(
    image.get_as_numpy(c=0, z=1, axes_order=["y", "x"]),
    cmap="gray",
    vmin=vmin,
    vmax=vmax,
)
axs[1].set_title("Final segmentation")
axs[1].imshow(
    label.get_as_numpy(z=1, axes_order=["y", "x"]),
    cmap=rand_cmap,
    interpolation="nearest",
)
for ax in axs:
    ax.axis("off")
fig.tight_layout()
print(figure_html(fig))

Next steps¶

Feature extraction — measure the objects you just segmented.
Masked images and labels — read data object-by-object.

Beyond the tutorials¶

The ngio workshop has hands-on marimo notebooks covering containers, images, labels and tables, and the processing iterators. Run them locally with uv, in the browser via molab, or read them as static pages.