segment_images

Script to segment objects from images.

For fiber-like objects, binarize and skeletonize the image, then use skan to extract branches coordinates. For polygon-like objects, binarize the image and detect objects and extract contours coordinates. For points, treat that as polygons then extract the centroids instead of contours. Finally, export the coordinates as collections in geojson files, importable in QuPath. Supports any number of channel of interest within the same image. One file output file per channel will be created.

This script uses cuisto.seg. It is designed to work on probability maps generated from a pixel classifier in QuPath, but might work on raw images.

Usage : fill-in the Parameters section of the script and run it. A "geojson" folder will be created in the parent directory of IMAGES_DIR. To exclude objects near the edges of an ROI, specify the path to masks stored as images with the same names as probabilities images (without their suffix).

author : Guillaume Le Goc (g.legoc@posteo.org) @ NeuroPSI version : 2024.12.10

CHANNELS_PARAMS = [{'name': 'cy5', 'target_channel': 0, 'proba_threshold': 0.85, 'qp_class': 'Fibers: Cy5', 'qp_color': [164, 250, 120]}, {'name': 'dsred', 'target_channel': 1, 'proba_threshold': 0.65, 'qp_class': 'Fibers: DsRed', 'qp_color': [224, 153, 18]}, {'name': 'egfp', 'target_channel': 2, 'proba_threshold': 0.85, 'qp_class': 'Fibers: EGFP', 'qp_color': [135, 11, 191]}] module-attribute

This should be a list of dictionary (one per channel) with keys :

• name: str, used as suffix for output geojson files, not used if only one channel
• target_channel: int, index of the segmented channel of the image, 0-based
• proba_threshold: float < 1, probability cut-off for that channel
• qp_class: str, name of QuPath classification
• qp_color: list of RGB values, associated color

EDGE_DIST = 0 module-attribute

Distance to brain edge to ignore, in µm. 0 to disable.

FILTERS = {'length_low': 1.5, 'area_low': 10, 'area_high': 1000, 'ecc_low': 0.0, 'ecc_high': 0.9, 'dist_thresh': 30} module-attribute

Dictionary with keys :

• length_low: minimal length in microns - for lines
• area_low: minimal area in µm² - for polygons and points
• area_high: maximal area in µm² - for polygons and points
• ecc_low: minimal eccentricity - for polygons and points (0 = circle)
• ecc_high: maximal eccentricity - for polygons and points (1 = line)
• dist_thresh: maximal inter-point distance in µm - for points

IMAGES_DIR = '/path/to/images' module-attribute

Full path to the images to segment.

IMG_SUFFIX = '_Probabilities.tiff' module-attribute

Images suffix, including extension. Masks must be the same name without the suffix.

MASKS_DIR = 'path/to/corresponding/masks' module-attribute

Full path to the masks, to exclude objects near the brain edges (set to None or empty string to disable this feature).

MASKS_EXT = 'tiff' module-attribute

Masks files extension.

MAX_PIX_VALUE = 255 module-attribute

Maximum pixel possible value to adjust proba_threshold.

ORIGINAL_PIXELSIZE = 0.45 module-attribute

Original images pixel size in microns. This is in case the pixel classifier uses a lower resolution, yielding smaller probability maps, so output objects coordinates need to be rescaled to the full size images. The pixel size is written in the "Image" tab in QuPath.

QUPATH_TYPE = 'detection' module-attribute

QuPath object type.

SEGTYPE = 'boutons' module-attribute

Type of segmentation.

get_geojson_dir(images_dir)

Get the directory of geojson files, which will be in the parent directory of images_dir.

If the directory does not exist, create it.

Parameters:

Name	Type	Description	Default
images_dir	str		required

Returns:

Name	Type	Description
geojson_dir	str

Source code in scripts/segmentation/segment_images.py

def get_geojson_dir(images_dir: str):
    """
    Get the directory of geojson files, which will be in the parent directory
    of `images_dir`.

    If the directory does not exist, create it.

    Parameters
    ----------
    images_dir : str

    Returns
    -------
    geojson_dir : str

    """

    geojson_dir = os.path.join(Path(images_dir).parent, "geojson")

    if not os.path.isdir(geojson_dir):
        os.mkdir(geojson_dir)

    return geojson_dir

get_geojson_properties(name, color, objtype='detection')

Return geojson objects properties as a dictionnary, ready to be used in geojson.Feature.

Parameters:

Name	Type	Description	Default
name	str	Classification name.	required
color	tuple or list	Classification color in RGB (3-elements vector).	required
objtype	str	Object type ("detection" or "annotation"). Default is "detection".	'detection'

Returns:

Name	Type	Description
props	dict

Source code in scripts/segmentation/segment_images.py

def get_geojson_properties(name: str, color: tuple | list, objtype: str = "detection"):
    """
    Return geojson objects properties as a dictionnary, ready to be used in geojson.Feature.

    Parameters
    ----------
    name : str
        Classification name.
    color : tuple or list
        Classification color in RGB (3-elements vector).
    objtype : str, optional
        Object type ("detection" or "annotation"). Default is "detection".

    Returns
    -------
    props : dict

    """

    return {
        "objectType": objtype,
        "classification": {"name": name, "color": color},
        "isLocked": "true",
    }

get_seg_method(segtype)

Determine what kind of segmentation is performed.

Segmentation kind are, for now, lines, polygons or points. We detect that based on hardcoded keywords.

Parameters:

Name	Type	Description	Default
segtype	str		required

Returns:

Name	Type	Description
seg_method	str

Source code in scripts/segmentation/segment_images.py

def get_seg_method(segtype: str):
    """
    Determine what kind of segmentation is performed.

    Segmentation kind are, for now, lines, polygons or points. We detect that based on
    hardcoded keywords.

    Parameters
    ----------
    segtype : str

    Returns
    -------
    seg_method : str

    """

    line_list = ["fibers", "axons", "fiber", "axon"]
    point_list = ["synapto", "synaptophysin", "syngfp", "boutons", "points"]
    polygon_list = ["cells", "polygon", "polygons", "polygon", "cell"]

    if segtype in line_list:
        seg_method = "lines"
    elif segtype in polygon_list:
        seg_method = "polygons"
    elif segtype in point_list:
        seg_method = "points"
    else:
        raise ValueError(
            f"Could not determine method to use based on segtype : {segtype}."
        )

    return seg_method

parameters_as_dict(images_dir, masks_dir, segtype, name, proba_threshold, edge_dist)

Get information as a dictionnary.

Parameters:

Name	Type	Description	Default
images_dir	str	Path to images to be segmented.	required
masks_dir	str	Path to images masks.	required
segtype	str	Segmentation type (eg. "fibers").	required
name	str	Name of the segmentation (eg. "green").	required
proba_threshold	float < 1	Probability threshold.	required
edge_dist	float	Distance in µm to the brain edge that is ignored.	required

Returns:

Name	Type	Description
params	dict

Source code in scripts/segmentation/segment_images.py

def parameters_as_dict(
    images_dir: str,
    masks_dir: str,
    segtype: str,
    name: str,
    proba_threshold: float,
    edge_dist: float,
):
    """
    Get information as a dictionnary.

    Parameters
    ----------
    images_dir : str
        Path to images to be segmented.
    masks_dir : str
        Path to images masks.
    segtype : str
        Segmentation type (eg. "fibers").
    name : str
        Name of the segmentation (eg. "green").
    proba_threshold : float < 1
        Probability threshold.
    edge_dist : float
        Distance in µm to the brain edge that is ignored.

    Returns
    -------
    params : dict

    """

    return {
        "images_location": images_dir,
        "masks_location": masks_dir,
        "type": segtype,
        "probability threshold": proba_threshold,
        "name": name,
        "edge distance": edge_dist,
    }

process_directory(images_dir, img_suffix='', segtype='', original_pixelsize=1.0, target_channel=0, proba_threshold=0.0, qupath_class='Object', qupath_color=[0, 0, 0], channel_suffix='', edge_dist=0.0, filters={}, masks_dir='', masks_ext='')

Main function, processes the .ome.tiff files in the input directory.

Parameters:

Name	Type	Description	Default
images_dir	str	Animal ID to process.	required
img_suffix	str	Images suffix, including extension.	''
segtype	str	Segmentation type.	''
original_pixelsize	float	Original images pixel size in microns.	1.0
target_channel	int	Index of the channel containning the objects of interest (eg. not the background), in the probability map (not the original images channels).	0
proba_threshold	float < 1	Probability below this value will be discarded (multiplied by MAX_PIXEL_VALUE)	0.0
qupath_class	str	Name of the QuPath classification.	'Object'
qupath_color	list of three elements	Color associated to that classification in RGB.	[0, 0, 0]
channel_suffix	str	Channel name, will be used as a suffix in output geojson files.	''
edge_dist	float	Distance to the edge of the brain masks that will be ignored, in microns. Set to 0 to disable this feature.	0.0
filters	dict	Filters values to include or excludes objects. See the top of the script.	{}
masks_dir	str	Path to images masks, to exclude objects found near the edges. The masks must be with the same name as the corresponding image to be segmented, without its suffix. Default is "", which disables this feature.	''
masks_ext	str	Masks files extension, without leading ".". Default is ""	''

Source code in scripts/segmentation/segment_images.py

def process_directory(
    images_dir: str,
    img_suffix: str = "",
    segtype: str = "",
    original_pixelsize: float = 1.0,
    target_channel: int = 0,
    proba_threshold: float = 0.0,
    qupath_class: str = "Object",
    qupath_color: list = [0, 0, 0],
    channel_suffix: str = "",
    edge_dist: float = 0.0,
    filters: dict = {},
    masks_dir: str = "",
    masks_ext: str = "",
):
    """
    Main function, processes the .ome.tiff files in the input directory.

    Parameters
    ----------
    images_dir : str
        Animal ID to process.
    img_suffix : str
        Images suffix, including extension.
    segtype : str
        Segmentation type.
    original_pixelsize : float
        Original images pixel size in microns.
    target_channel : int
        Index of the channel containning the objects of interest (eg. not the
        background), in the probability map (*not* the original images channels).
    proba_threshold : float < 1
        Probability below this value will be discarded (multiplied by `MAX_PIXEL_VALUE`)
    qupath_class : str
        Name of the QuPath classification.
    qupath_color : list of three elements
        Color associated to that classification in RGB.
    channel_suffix : str
        Channel name, will be used as a suffix in output geojson files.
    edge_dist : float
        Distance to the edge of the brain masks that will be ignored, in microns. Set to
        0 to disable this feature.
    filters : dict
        Filters values to include or excludes objects. See the top of the script.
    masks_dir : str, optional
        Path to images masks, to exclude objects found near the edges. The masks must be
        with the same name as the corresponding image to be segmented, without its
        suffix. Default is "", which disables this feature.
    masks_ext : str, optional
        Masks files extension, without leading ".". Default is ""

    """

    # -- Preparation
    # get segmentation type
    seg_method = get_seg_method(segtype)

    # get output directory path
    geojson_dir = get_geojson_dir(images_dir)

    # get images list
    images_list = [
        os.path.join(images_dir, filename)
        for filename in os.listdir(images_dir)
        if filename.endswith(img_suffix)
    ]

    # write parameters
    parameters = parameters_as_dict(
        images_dir, masks_dir, segtype, channel_suffix, proba_threshold, edge_dist
    )
    param_file = os.path.join(geojson_dir, "parameters" + channel_suffix + ".txt")
    if os.path.isfile(param_file):
        raise FileExistsError("Parameters file already exists.")
    else:
        write_parameters(param_file, parameters, filters, original_pixelsize)

    # convert parameters to pixels in probability map
    pixelsize = hq.seg.get_pixelsize(images_list[0])  # get pixel size
    edge_dist = int(edge_dist / pixelsize)
    filters = hq.seg.convert_to_pixels(filters, pixelsize)

    # get rescaling factor
    rescale_factor = pixelsize / original_pixelsize

    # get GeoJSON properties
    geojson_props = get_geojson_properties(
        qupath_class, qupath_color, objtype=QUPATH_TYPE
    )

    # -- Processing
    pbar = tqdm(images_list)
    for imgpath in pbar:
        # build file names
        imgname = os.path.basename(imgpath)
        geoname = imgname.replace(img_suffix, "")
        geojson_file = os.path.join(
            geojson_dir, geoname + "_segmentation" + channel_suffix + ".geojson"
        )

        # checks if output file already exists
        if os.path.isfile(geojson_file):
            continue

        # read images
        pbar.set_description(f"{geoname}: Loading...")
        img = tifffile.imread(imgpath, key=target_channel)
        if (edge_dist > 0) & (len(masks_dir) != 0):
            mask = tifffile.imread(os.path.join(masks_dir, geoname + "." + masks_ext))
            mask = hq.seg.pad_image(mask, img.shape)  # resize mask
            # apply mask, eroding from the edges
            img = img * hq.seg.erode_mask(mask, edge_dist)

        # image processing
        pbar.set_description(f"{geoname}: IP...")

        # threshold probability and binarization
        img = img >= proba_threshold * MAX_PIX_VALUE

        # segmentation
        pbar.set_description(f"{geoname}: Segmenting...")

        if seg_method == "lines":
            collection = hq.seg.segment_lines(
                img,
                geojson_props,
                minsize=filters["length_low"],
                rescale_factor=rescale_factor,
            )

        elif seg_method == "polygons":
            collection = hq.seg.segment_polygons(
                img,
                geojson_props,
                area_min=filters["area_low"],
                area_max=filters["area_high"],
                ecc_min=filters["ecc_low"],
                ecc_max=filters["ecc_high"],
                rescale_factor=rescale_factor,
            )

        elif seg_method == "points":
            collection = hq.seg.segment_points(
                img,
                geojson_props,
                area_min=filters["area_low"],
                area_max=filters["area_high"],
                ecc_min=filters["ecc_low"],
                ecc_max=filters["ecc_high"],
                dist_thresh=filters["dist_thresh"],
                rescale_factor=rescale_factor,
            )
        else:
            # we already printed an error message
            return

        # save geojson
        pbar.set_description(f"{geoname}: Saving...")
        with open(geojson_file, "w") as fid:
            fid.write(geojson.dumps(collection))

write_parameters(outfile, parameters, filters, original_pixelsize)

Write parameters to outfile.

A timestamp will be added. Parameters are written as key = value, and a [filters] is added before filters parameters.

Parameters:

Name	Type	Description	Default
outfile	str	Full path to the output file.	required
parameters	dict	General parameters.	required
filters	dict	Filters parameters.	required
original_pixelsize	float	Size of pixels in original image.	required

Source code in scripts/segmentation/segment_images.py

def write_parameters(
    outfile: str, parameters: dict, filters: dict, original_pixelsize: float
):
    """
    Write parameters to `outfile`.

    A timestamp will be added. Parameters are written as key = value,
    and a [filters] is added before filters parameters.

    Parameters
    ----------
    outfile : str
        Full path to the output file.
    parameters : dict
        General parameters.
    filters : dict
        Filters parameters.
    original_pixelsize : float
        Size of pixels in original image.

    """

    with open(outfile, "w") as fid:
        fid.writelines(f"date = {datetime.now().strftime('%d-%B-%Y %H:%M:%S')}\n")

        fid.writelines(f"original_pixelsize = {original_pixelsize}\n")

        for key, value in parameters.items():
            fid.writelines(f"{key} = {value}\n")

        fid.writelines("[filters]\n")

        for key, value in filters.items():
            fid.writelines(f"{key} = {value}\n")