Skip to content

Commit

Permalink
- use lower resolution pixel classifier
Browse files Browse the repository at this point in the history
  • Loading branch information
GuillaumeLeGoc committed Dec 10, 2024
1 parent 73e2db6 commit cfd0966
Show file tree
Hide file tree
Showing 4 changed files with 87 additions and 25 deletions.
6 changes: 6 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
@@ -1,2 +1,8 @@
This file might no be up to date nor complete. Please check the Releases page for more information on versions.

## Version 2024.12.10
- Fixed offset writing geojson file in segmentation module.
- Added original_pixelsize parameters to rescale coordinates to match final image size. This allows the use of QuPath pixel classifier trained on resized image, eg. with a lower Resolution parameter (higher pixel size).

## Version 2024.11.19
- Initial public release.
57 changes: 45 additions & 12 deletions histoquant/seg.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,11 @@ def get_pixelsize(image_name: str) -> float:
"""

with tifffile.TiffFile(image_name) as tif:
return 1e6 / tif.pages[0].tags["XResolution"].value[0]
# XResolution is a tuple, numerator, denomitor. The inverse is the pixel size
return (
tif.pages[0].tags["XResolution"].value[1]
/ tif.pages[0].tags["XResolution"].value[0]
)


def convert_to_pixels(filters, pixelsize):
Expand Down Expand Up @@ -142,7 +146,7 @@ def get_image_skeleton(img: np.ndarray, minsize=0) -> np.ndarray:


def get_collection_from_skel(
skeleton: Skeleton, properties: dict
skeleton: Skeleton, properties: dict, rescale_factor: float = 1.0, offset=0.5
) -> geojson.FeatureCollection:
"""
Get the coordinates of each skeleton path as a GeoJSON Features in a
Expand All @@ -154,6 +158,11 @@ def get_collection_from_skel(
skeleton : skan.Skeleton
properties : dict
QuPatj objects' properties.
rescale_factor : float
Rescale output coordinates by this factor.
offset : float
Shift coordinates by this amount, typically to get pixel centers or edges.
Default is 0.5.
Returns
-------
Expand All @@ -171,7 +180,7 @@ def get_collection_from_skel(
collection.append(
geojson.Feature(
geometry=shapely.LineString(
skeleton.path_coordinates(ind)[:, ::-1]
(skeleton.path_coordinates(ind)[:, ::-1] + offset) * rescale_factor
), # shape object
properties=prop, # object properties
id=str(uuid.uuid4()), # object uuid
Expand All @@ -182,7 +191,7 @@ def get_collection_from_skel(


def get_collection_from_poly(
contours: list, properties: dict
contours: list, properties: dict, rescale_factor: float = 1.0, offset: float = 0.5
) -> geojson.FeatureCollection:
"""
Gather coordinates in the list and put them in GeoJSON format as Polygons.
Expand All @@ -195,17 +204,23 @@ def get_collection_from_poly(
contours : list
properties : dict
QuPatj objects' properties.
rescale_factor : float
Rescale output coordinates by this factor.
offset : float
Shift coordinates by this amount, typically to get pixel centers or edges.
Default is 0.5.
Returns
-------
collection : geojson.FeatureCollection
A FeatureCollection ready to be written as geojson.
"""

collection = [
geojson.Feature(
geometry=shapely.Polygon(np.fliplr(contour)), # shape object
geometry=shapely.Polygon(
np.fliplr((contour + offset) * rescale_factor)
), # shape object
properties=properties, # object properties
id=str(uuid.uuid4()), # object uuid
)
Expand All @@ -216,7 +231,7 @@ def get_collection_from_poly(


def get_collection_from_points(
coords: list, properties: dict
coords: list, properties: dict, rescale_factor: float = 1.0, offset: float = 0.5
) -> geojson.FeatureCollection:
"""
Gather coordinates from `coords` and put them in GeoJSON format.
Expand All @@ -228,6 +243,8 @@ def get_collection_from_points(
----------
coords : list
properties : dict
rescale_factor : float
Rescale output coordinates by this factor.
Returns
-------
Expand All @@ -237,7 +254,9 @@ def get_collection_from_points(

collection = [
geojson.Feature(
geometry=shapely.Point(np.flip(coord)), # shape object
geometry=shapely.Point(
np.flip((coord + offset) * rescale_factor)
), # shape object
properties=properties, # object properties
id=str(uuid.uuid4()), # object uuid
)
Expand All @@ -248,7 +267,7 @@ def get_collection_from_points(


def segment_lines(
img: np.ndarray, geojson_props: dict, minsize=0.0
img: np.ndarray, geojson_props: dict, minsize=0.0, rescale_factor=1.0
) -> geojson.FeatureCollection:
"""
Wraps skeleton analysis to get paths coordinates.
Expand All @@ -261,6 +280,8 @@ def segment_lines(
GeoJSON properties of objects.
minsize : float
Minimum size in pixels for an object.
rescale_factor : float
Rescale output coordinates by this factor.
Returns
-------
Expand All @@ -273,7 +294,9 @@ def segment_lines(

# get paths coordinates as FeatureCollection
skeleton = Skeleton(skel, keep_images=False)
return get_collection_from_skel(skeleton, geojson_props)
return get_collection_from_skel(
skeleton, geojson_props, rescale_factor=rescale_factor
)


def segment_polygons(
Expand All @@ -283,6 +306,7 @@ def segment_polygons(
area_max: float = np.inf,
ecc_min: float = 0.0,
ecc_max: float = 1.0,
rescale_factor: float = 1.0,
) -> geojson.FeatureCollection:
"""
Polygon segmentation.
Expand All @@ -297,6 +321,8 @@ def segment_polygons(
Minimum and maximum area in pixels for an object.
ecc_min, ecc_max : float
Minimum and maximum eccentricity for an object.
rescale_factor: float
Rescale output coordinates by this factor.
Returns
-------
Expand Down Expand Up @@ -328,7 +354,9 @@ def segment_polygons(
label_image = label_image > 0
contours = measure.find_contours(label_image)

return get_collection_from_poly(contours, geojson_props)
return get_collection_from_poly(
contours, geojson_props, rescale_factor=rescale_factor
)


def segment_points(
Expand All @@ -339,6 +367,7 @@ def segment_points(
ecc_min: float = 0,
ecc_max: float = 1,
dist_thresh: float = 0,
rescale_factor: float = 1,
) -> geojson.FeatureCollection:
"""
Point segmentation.
Expand All @@ -359,6 +388,8 @@ def segment_points(
dist_thresh : float
Maximal distance in pixels between objects before considering them as isolated and remove them.
0 disables it.
rescale_factor : float
Rescale output coordinates by this factor.
Returns
-------
Expand Down Expand Up @@ -410,4 +441,6 @@ def segment_points(
# get points coordinates
coords = np.argwhere(bw)

return get_collection_from_points(coords, geojson_props)
return get_collection_from_points(
coords, geojson_props, rescale_factor=rescale_factor
)
2 changes: 1 addition & 1 deletion pyproject.toml
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
[project]
name = "histoquant"
version = "2024.11.27"
version = "2024.12.10"
authors = [{ name = "Guillaume Le Goc", email = "[email protected]" }]
description = "Quantification of objects in histological slices"
readme = "README.md"
Expand Down
47 changes: 35 additions & 12 deletions scripts/segmentation/segment_images.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,8 +18,8 @@
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 ([email protected]) @ NeuroPSI
Version : 2024.11.27
author : Guillaume Le Goc ([email protected]) @ NeuroPSI
version : 2024.12.10
"""

Expand All @@ -37,17 +37,22 @@
pd.options.mode.copy_on_write = True # prepare for pandas 3

# --- Parameters
IMAGES_DIR = "/path/to/images/to/segment"
IMAGES_DIR = "/path/to/images"
"""Full path to the images to segment."""
MASKS_DIR = "path/to/corresponding/masks"
"""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"
"""Masks files extension."""
SEGTYPE = "fibers"
SEGTYPE = "boutons"
"""Type of segmentation."""
IMG_SUFFIX = "_Probabilities.tiff"
"""Images suffix, including extension. Masks must be the same name without the suffix."""
ORIGINAL_PIXELSIZE = 0.4500
"""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."""

CHANNELS_PARAMS = [
{
Expand All @@ -60,7 +65,7 @@
{
"name": "dsred",
"target_channel": 1,
"proba_threshold": 0.85,
"proba_threshold": 0.65,
"qp_class": "Fibers: DsRed",
"qp_color": [224, 153, 18],
},
Expand All @@ -80,13 +85,13 @@
- qp_class: str, name of QuPath classification
- qp_color: list of RGB values, associated color"""

EDGE_DIST = 50
EDGE_DIST = 0
"""Distance to brain edge to ignore, in µm. 0 to disable."""

FILTERS = {
"length_low": 1.5, # minimal length in microns - for lines
"area_low": 1.1, # minimal area in µm² - for polygons and points
"area_high": 10, # maximal area in µm² - for polygons and points
"area_low": 10, # minimal area in µm² - for polygons and points
"area_high": 1000, # maximal area in µm² - for polygons and points
"ecc_low": 0.0, # minimal eccentricity - for polygons and points (0 = circle)
"ecc_high": 0.9, # maximal eccentricity - for polygons and points (1 = line)
"dist_thresh": 30, # maximal inter-point distance in µm - for points
Expand Down Expand Up @@ -237,7 +242,9 @@ def parameters_as_dict(
}


def write_parameters(outfile: str, parameters: dict, filters: dict):
def write_parameters(
outfile: str, parameters: dict, filters: dict, original_pixelsize: float
):
"""
Write parameters to `outfile`.
Expand All @@ -252,12 +259,16 @@ def write_parameters(outfile: str, parameters: dict, filters: 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")

Expand All @@ -271,6 +282,7 @@ 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",
Expand All @@ -292,6 +304,8 @@ def process_directory(
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).
Expand Down Expand Up @@ -339,13 +353,16 @@ def process_directory(
if os.path.isfile(param_file):
raise FileExistsError("Parameters file already exists.")
else:
write_parameters(param_file, parameters, filters)
write_parameters(param_file, parameters, filters, original_pixelsize)

# convert parameters to pixels
# 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
Expand Down Expand Up @@ -385,7 +402,10 @@ def process_directory(

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

elif seg_method == "polygons":
Expand All @@ -396,6 +416,7 @@ def process_directory(
area_max=filters["area_high"],
ecc_min=filters["ecc_low"],
ecc_max=filters["ecc_high"],
rescale_factor=rescale_factor,
)

elif seg_method == "points":
Expand All @@ -407,6 +428,7 @@ def process_directory(
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
Expand Down Expand Up @@ -445,6 +467,7 @@ def make_suffix(s):
IMAGES_DIR,
img_suffix=IMG_SUFFIX,
segtype=SEGTYPE,
original_pixelsize=ORIGINAL_PIXELSIZE,
target_channel=param["target_channel"],
proba_threshold=param["proba_threshold"],
qupath_class=param["qp_class"],
Expand Down

0 comments on commit cfd0966

Please sign in to comment.