Tasks ¶

Multi-angle fusion ¶

timagetk.tasks.fusion.fusion_on_first(images, method='mean', init_trsfs=None, **kwargs)[source]¶

Fuse a list of image by registering them on the first of the list.

Parameters:

images (list of SpatialImage) – List of image to fuse together
method (str, optional) – Image fusion method to use, ‘mean’ by default, see AVERAGING_METHODS
init_trsfs (list of Trsf or None, optional) – List of transformations to use to initialize registration, must be of type rigid.
vectorfield (bool) – If True (default is False), also compute a non-linear deformation to match the images.
fuse_reference (bool) – If True (default), add the reference image to the list of images to fuse.
final_vxs (float) – If specified, the template image will be the isometric resampling of th reference image to this voxel-size float value
super_resolution (bool) – If True (default), compute a reference frame with minimum voxel-size.
global_averaging (bool) – If True (default is False), perform global averaging of each voxel by the total number of images. Else performs local averaging based on the number of defined imaged per voxels after registrations.

Returns:

timagetk.SpatialImage – The averaged images after registration.
list of Trsf – The list of initial transformation files.

Notes

Prior to averaging, masks are generated to average each voxel by the number of locally defined images and not the total number of images.

Examples

>>> import numpy as np
>>> from timagetk.tasks.fusion import fusion_on_first
>>> from timagetk.algorithms.pointmatching import pointmatching
>>> from timagetk.algorithms.resample import resample
>>> from timagetk.algorithms.trsf import inv_trsf
>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_data
>>> from timagetk.visu.stack import orthogonal_view

>>> # Using shared data as example, images are multi-angle view of the same floral meristem:
>>> fname = 'p58-t0-a{}.lsm'
>>> list_files = [shared_data(fname.format(a), 'p58') for a in range(3)]
>>> list_img = [imread(fname) for fname in list_files]
>>> # Load shared multi-angle landmarks for the first time point (t0) of 'p58' shared time-series
>>> ref_pts_01 = np.loadtxt(shared_data('p58_t0_reference_ldmk-01.txt', 'p58'))
>>> flo_pts_01 = np.loadtxt(shared_data('p58_t0_floating_ldmk-01.txt', 'p58'))
>>> ref_pts_02 = np.loadtxt(shared_data('p58_t0_reference_ldmk-02.txt', 'p58'))
>>> flo_pts_02 = np.loadtxt(shared_data('p58_t0_floating_ldmk-02.txt', 'p58'))
>>> # Creates manual initialization transformations with `pointmatching` algorithm:
>>> trsf_01 = pointmatching(flo_pts_01, ref_pts_01, template_img=list_img[0], method='rigid')
>>> trsf_02 = pointmatching(flo_pts_02, ref_pts_02, template_img=list_img[0], method='rigid')

>>> # Example 1 - Fuse the image, in super-resolution, after affine registration with manual initialization
>>> fused_img, _ = fusion_on_first(list_img, init_trsfs=[trsf_01, trsf_02], super_resolution=True)
>>> orthogonal_view(fused_img,suptitle="Fusion: centered local average image after affine registration with manual initialization")

>>> # Example 2 - Fuse the image, in super-resolution, after affine registration with manual initialization
>>> fused_img, _ = fusion_on_first(list_img, init_trsfs=[trsf_01, trsf_02], super_resolution=True, global_averaging=True)
>>> orthogonal_view(fused_img,suptitle="Fusion: centered global average image after affine registration with manual initialization")

>>> # Example 3 - Fuse the image, in super-resolution, after affine & non-linear registration with manual initialization
>>> fused_img, _ = fusion_on_first(list_img, init_trsfs=[trsf_01, trsf_02], super_resolution=True, vectorfield=True)
>>> orthogonal_view(fused_img,suptitle="Fusion: centered average image after affine & non-linear registration with manual initialization")

timagetk.tasks.fusion.get_ref_im_float_list(images, ref_img_id)[source]¶

Return the reference image and the list of floating images from given index.

Parameters:

images (list(SpatialImages)) – The list of SpatialImage to consider.
ref_img_id (int) – The index of the reference image in the initial images list.

Returns:

timagetk.SpatialImage – The reference SpatialImage.
list(SpatialImages) – The list of floating SpatialImages.

Examples

>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_data
>>> from timagetk.tasks.fusion import get_ref_im_float_list
>>> # Using shared data as example, images are multi-angle view of the same floral meristem:
>>> fname = 'p58-t0-a{}.lsm'
>>> list_files = [shared_data(fname.format(a), 'p58') for a in range(3)]
>>> list_img = [imread(fname) for fname in list_files]
>>> # Example 1 - Use the first image of the list as reference:
>>> ref_img, float_imgs = get_ref_im_float_list(list_img, 0)  # index starts at 0!
>>> print(f"Reference image: {ref_img.filename}")
>>> print(f"Floating images: {', '.join([flo_img.filename for flo_img in float_imgs])}")
Reference image: p58-t0-a0.lsm
Floating images: p58-t0-a1.lsm, p58-t0-a2.lsm
>>> # Example 2 - Use the third image of the list as reference:
>>> ref_img, float_imgs = get_ref_im_float_list(list_img, 2)  # index starts at 0!
>>> print(f"Reference image: {ref_img.filename}")
>>> print(f"Floating images: {', '.join([flo_img.filename for flo_img in float_imgs])}")
Reference image: p58-t0-a2.lsm
Floating images: p58-t0-a0.lsm, p58-t0-a1.lsm

timagetk.tasks.fusion.iterative_fusion(images, method='mean', init_trsfs=None, n_iter=3, **kwargs)[source]¶

Iterative fusion of a list of multi-angle images.

Parameters:

images (list of SpatialImage) – List of image to fuse together
method (str, optional) – Image fusion method to use, ‘mean’ by default, see AVERAGING_METHODS
init_trsfs (list of Trsf, optional) – List of transformations to use to initialize registration, must be of type rigid.
n_iter (int, optional) – Number of iterations to perform, including the initial registration on first image.
fuse_reference (bool) – If True (default), add the reference image to the list of images to fuse.
global_averaging (bool) – If True (default is False), perform global averaging of each voxel by the total number of images. Else performs local averaging based on the number of defined image per voxels after registrations.
super_resolution (bool) – Defines the voxelsize of the returned fused image. If True (default is False), use the smallest voxelsize of the reference image as isometric voxelsize. Else, use the largest voxelsize of the reference image as isometric voxelsize.
vectorfield_at_last (bool) – If True (default is False), perform non-linear registration at last iteration step. Else, do not performs non-linear registration, only rigid and affine.
final_vxs (float) – The voxel-size of the final isometric image, override the voxelsize computed by super_resolution.
interpolation ({'linear', 'cspline'}) – The type of interpolation to performs when resampling and applying transformation to intensity images. Defaults to ‘linear’.

Returns:

timagetk.SpatialImage – The fused image after iterative registration.

Notes

Only the last iteration will use the maximum resolution. The previous steps are performed on 2x down-sampled images.

Examples

>>> import numpy as np
>>> from timagetk.tasks.fusion import iterative_fusion
>>> from timagetk.algorithms.pointmatching import pointmatching
>>> from timagetk.algorithms.resample import resample
>>> from timagetk.algorithms.trsf import inv_trsf
>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_data
>>> from timagetk.visu.stack import stack_browser
>>> from timagetk.visu.stack import orthogonal_view
>>> # Using shared data as example, images are multi-angle view of the same floral meristem:
>>> fname = 'p58-t0-a{}.lsm'
>>> list_files = [shared_data(fname.format(a), 'p58') for a in range(3)]
>>> list_img = [imread(fname) for fname in list_files]
>>> ref_pts_01 = np.loadtxt(shared_data('p58_t0_reference_ldmk-01.txt', 'p58'))
>>> flo_pts_01 = np.loadtxt(shared_data('p58_t0_floating_ldmk-01.txt', 'p58'))
>>> ref_pts_02 = np.loadtxt(shared_data('p58_t0_reference_ldmk-02.txt', 'p58'))
>>> flo_pts_02 = np.loadtxt(shared_data('p58_t0_floating_ldmk-02.txt', 'p58'))
>>> # Creates manual initialization transformations with `pointmatching` algorithm:
>>> trsf_01 = pointmatching(flo_pts_01, ref_pts_01, template_img=list_img[0], method='rigid')
>>> trsf_02 = pointmatching(flo_pts_02, ref_pts_02, template_img=list_img[0], method='rigid')

>>> # Example #1 - Fast iterative fusion: non-linear registration only at last iteration
>>> vf_fused_img = iterative_fusion(list_img, init_trsfs=[trsf_01, trsf_02], vectorfield_at_last=True)
>>> orthogonal_view(vf_fused_img,suptitle="p58-t0 multi-angle iterative fusion in super-resolution from AFFINE registration with manual initialization")

>>> # Example #2 - Slower iterative fusion: non-linear registration only at every iteration
>>> fused_img = iterative_fusion(list_img, init_trsfs=[trsf_01, trsf_02])
>>> orthogonal_view(fused_img,suptitle="p58-t0 multi-angle iterative fusion in super-resolution from AFFINE registration with manual initialization")

>>> from timagetk.components.multi_channel import combine_channels
>>> from timagetk.visu.mplt import grayscale_imshow
>>> cuts = [1 / 6., 1 / 4., 1 / 2., 2 / 3., 3 / 4.]
>>> z_sh = fused_img.get_shape('z')
>>> z_slices = [int(round(c * z_sh, 0)) for c in cuts]
>>> blend = combine_channels([fused_img, vf_fused_img], colors=["red", "green"])
>>> thumbs = [blend.get_slice(zsl) for zsl in z_slices]
>>> titles = [f"z-slice {zsl}/{z_sh}" for zsl in z_slices]
>>> grayscale_imshow(thumbs, titles=titles)

Intensity image segmentation ¶

timagetk.tasks.segmentation.background_detection(image, threshold=None, sigma=None, equalize_hist=True, contrast_stretch=False, **kwargs)[source]¶

Detect the background position, i.e. not part of the object.

Parameters:

image (timagetk.SpatialImage or timagetk.MultiChannelImage) – The intensity image to segment.
threshold (int, optional) – The threshold to use for bacground detection, if not defined use the li_threshold method to estimate it.
sigma (float, optional) – Sigma value used with Gaussian smoothing
equalize_hist (bool, optional) – If True (default), performs adaptive histogram equalization to image prior to any other step.
contrast_stretch (bool, optional) – If True (default is False), performs contrast stretching to image prior to any other step.
sigma – Sigma value, in real units, used to smooth the input image, use 0 to skip. If None (default), automatically set to 1.5 * max(image.voxelsize).
channel (str) – If a MultiChannelImage is used as input image, select the channel to use with this algorithm.

Returns:

timagetk.components.spatial_image.SpatialImage – Mask image (binary) where 1 indicate background position.

See also

timagetk.array_util.guess_intensity_threshold, scipy.ndimage.binary_fill_holes

Notes

This strategy assumes the background is where “there is no signal”. We fill holes in the detected structure, slice by slice, moving along the z-axis.

Example

>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_dataset
>>> from timagetk.tasks.segmentation import background_detection
>>> img_path = shared_dataset("p58")[0]
>>> image = imread(img_path)
>>> bkgd_img, _ = background_detection(image)
>>> from timagetk.visu.stack import stack_browser
>>> b = stack_browser(bkgd_img, val_range=[0,1])

timagetk.tasks.segmentation.seed_detection(intensity_image, h_min, **kwargs)[source]¶

Performs a seed detection an intensity image.

Parameters:

intensity_image (timagetk.SpatialImage or timagetk.MultiChannelImage) – The intensity image to segment
h_min (int) – Height-minima value to use for seed detection
channel (str) – If a MultiChannelImage is used as input image, select the channel to use with this algorithm.
process (dict) – A dictionary keeping track of the algorithms and their parameters.

Returns:

LabelledImage – The seeds image.
dict – Dictionary summarizing the process with used parameters

timagetk.tasks.segmentation.seeded_watershed(intensity_image, h_min, **kwargs)[source]¶

Performs a seeded watershed of an intensity image.

Parameters:

intensity_image (timagetk.SpatialImage or timagetk.MultiChannelImage) – The intensity image to segment.
h_min (int) – Height-minima value to use for seed detection.
channel (str) – If a MultiChannelImage is used as input image, select the channel to use with this algorithm.
labelchoice ({"first", "min", "most"}) – How to deal with “labels conflicts”, i.e. where several labels meet.
process (dict) – A dictionary keeping track of the algorithms and their parameters.

Returns:

LabelledImage – The segmented image
LabelledImage – The seeds image
dict – Dictionary summarizing the process with used parameters

Notes

Explanations of available labelchoice keyword argument:

first: the first label wins (default);

min: the less represented label wins;

most: the most represented label wins;

Examples

>>> from timagetk import TissueImage3D
>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_dataset
>>> from timagetk.tasks.segmentation import seeded_watershed
>>> from timagetk.visu.stack import stack_browser
>>> img_path = shared_dataset("p58")[0]
>>> image = imread(img_path)

>>> # Example 1 - Basic example:
>>> seg_img, _, _ = seeded_watershed(image, 15)
>>> b = stack_browser(TissueImage3D(seg_img, background=1), cmap='glasbey')

>>> # Example 2 - Add a pre-processing step:
>>> from timagetk.tasks.segmentation import watershed_preprocessing
>>> proc_image, _ = watershed_preprocessing(image,0.5,equalize_hist=True)
>>> seg_img, _, _ = seeded_watershed(proc_image, 25)
>>> b = stack_browser(TissueImage3D(seg_img, background=1), cmap='glasbey')

timagetk.tasks.segmentation.watershed_postprocessing(intensity_image, segmented_image, seeds_image, min_size=None, max_size=None, cell_sigma=1, **kwargs)[source]¶

Performs the post-processing step after watershed segmentation of an intensity image.

Parameters:

intensity_image (timagetk.components.spatial_image.SpatialImage) – The intensity image to segment
segmented_image (timagetk.LabelledImage) – The segmented image from the seeded watershed step
seeds_image (timagetk.LabelledImage) – The seeds image from the seeded watershed step
min_size (float, optional) – A minimal area (if 2D) or volume (if 3D), in real units, to accept a label, None by default. If equal or lower, remove it from the seed image and re-run the watershed.
max_size (float, optional) – A maximal area (if 2D) or volume (if 3D), in real units, to accept a label, None by default. If equal or superior, remove it from the seed image and re-run the watershed.
cell_sigma (int, optional) – Apply a “cell smoothing” with a sigma in voxel units, 1 by default. Set to 0 to avoid this step. Warnings, even small values (like 2) can drastically change the topology!
sizes_data (dict) – Dictionary of size data to use.
labelchoice ({"first", "min", "most"}) – How to deal with “labels conflicts”, i.e. where several labels meet.
process (dict) – A dictionary keeping track of the algorithms and their parameters.

Returns:

LabelledImage – The segmented image
LabelledImage – The seeds image
dict – Dictionary summarizing the process with used parameters

Notes

We recommend performing this step to remove too small cells (un-realistic) and smooth the cells interfaces.

Explanations of available labelchoice keyword argument:

first: the first label wins (default);

min: the less represented label wins;

most: the most represented label wins;

Examples

>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_dataset
>>> from timagetk.tasks.segmentation import watershed_preprocessing
>>> from timagetk.tasks.segmentation import seeded_watershed
>>> from timagetk.tasks.segmentation import watershed_postprocessing
>>> from timagetk.visu.stack import stack_browser
>>> img_path = shared_dataset("p58")[0]
>>> image = imread(img_path)
>>> proc_image, _ = watershed_preprocessing(image,0.5,equalize_hist=True)
>>> seg_img, seeds, _ = seeded_watershed(proc_image, 25)
>>> seg_img, _, _ = watershed_postprocessing(proc_image,seg_img,seeds,min_size=50,cell_sigma=1)
>>> b = stack_browser(seg_img, cmap='glasbey')

timagetk.tasks.segmentation.watershed_preprocessing(image, sigma=None, equalize_hist=True, contrast_stretch=False, **kwargs)[source]¶

Performs the pre-processing step prior to watershed segmentation of an intensity image.

A Gaussian smoothing is applied prior to local minima detection to avoid detecting to many small local minima due to noise.

Parameters:

image (timagetk.SpatialImage or timagetk.MultiChannelImage) – The intensity image to segment.
sigma (float, optional) – Sigma value, in real units, used to smooth the input image, use 0 to skip. If None (default), automatically set to 1.5 * max(image.voxelsize).
equalize_hist (bool, optional) – If True (default), performs adaptive histogram equalization to image prior to any other step.
contrast_stretch (bool, optional) – If True (default is False), performs contrast stretching to image prior to any other step.
channel (str) – If a MultiChannelImage is used as input image, select the channel to use with this algorithm.
process (dict) – A dictionary keeping track of the algorithms and their parameters.

Returns:

timagetk.SpatialImage – The segmented image
dict – Dictionary summarizing the process with used parameters

Notes

We recommend performing automatic global contrast stretching OR adaptive histogram equalization to improve overall segmentation quality. It is done prior to the Gaussian smoothing step.

Examples

>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_dataset
>>> from timagetk.tasks.segmentation import watershed_preprocessing
>>> from timagetk.visu.stack import stack_browser
>>> img_path = shared_dataset("p58")[0]
>>> image = imread(img_path)
>>> proc_image, _ = watershed_preprocessing(image, 1.0, equalize_hist=True)
>>> stack_browser(proc_image, cmap='gray')

>>> from timagetk import SpatialImage
>>> from skimage.morphology import white_tophat, ball
>>> noise_img = SpatialImage(white_tophat(image.get_array(), ball(3)), voxelsize=image.voxelsize)
>>> stack_browser(noise_img, cmap='gray')
>>> denoise_img = SpatialImage(image.get_array() - noise_img, voxelsize=image.voxelsize)
>>> stack_browser(denoise_img, cmap='gray')
>>> proc_image, _ = watershed_preprocessing(denoise_img, 1.0, equalize_hist=True)
>>> stack_browser(proc_image, cmap='gray')

>>> from timagetk import MultiChannelImage
>>> im_url = "https://zenodo.org/record/3737795/files/qDII-CLV3-PIN1-PI-E35-LD-SAM4.czi"
>>> ch_names = ['DII-VENUS-N7', 'pPIN1:PIN1-GFP', 'Propidium Iodide', 'pRPS5a:TagBFP-SV40', 'pCLV3:mCherry-N7']
>>> colors = ['yellow', 'green', 'red', 'blue', 'purple']
>>> image = imread(im_url, channel_names=ch_names)
>>> proc_image, _ = watershed_preprocessing(image, 1.0, equalize_hist=True, channel='Propidium Iodide')

timagetk.tasks.segmentation.watershed_segmentation(intensity_image, h_min, sigma=None, equalize_hist=True, contrast_stretch=False, min_size=50, max_size=None, cell_sigma=1)[source]¶

Performs watershed segmentation of an intensity image after detection of connexe local minima.

A Gaussian smoothing is applied prior to local minima detection to avoid detecting to many small local minima due to noise. The (stretched/equalized) smoothed intensity image is used for seed detection (connexe local minima) AND watershed.

Parameters:

intensity_image (timagetk.SpatialImage or timagetk.MultiChannelImage) – The intensity image to segment.
h_min (int) – Height-minima value to use for seed detection
sigma (float, optional) – Sigma value used with Gaussian smoothing
equalize_hist (bool, optional) – If True (default), performs adaptive histogram equalization to image prior to any other step.
contrast_stretch (bool, optional) – If True (default is False), performs contrast stretching to image prior to any other step.
min_size (float, optional) – A minimal volume or area, in real units, to accept a label, 50 µm³ by default. If equal or lower, remove it from the seed image and re-run the watershed.
max_size (float, optional) – A maximal volume or area, in real units, to accept a label, None by default. If equal or superior, remove it from the seed image and re-run the watershed.
cell_sigma (int, optional) – Apply a “cell smoothing” with a sigma in voxel units, 1 by default. Set to 0 to avoid this step. Warnings, even small values (like 2) can drastically change the topology!
channel (str) – If a MultiChannelImage is used as input image, select the channel to use with this algorithm.

Returns:

LabelledImage – The segmented image.
LabelledImage – The seeds image.
dict – Dictionary summarizing the process with used parameters.

Notes

We recommend performing automatic global contrast stretching OR adaptive histogram equalization to improve overall segmentation quality. It is done prior to the Gaussian smoothing step.

Examples

>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_dataset
>>> from timagetk.tasks.segmentation import watershed_segmentation
>>> from timagetk.visu.stack import stack_browser
>>> img_path = shared_dataset("p58")[0]
>>> image = imread(img_path)
>>> # Example 1 - Basic example:
>>> seg_img, seeds, _ = watershed_segmentation(image, 5)
>>> from timagetk.visu.util import greedy_colormap
>>> b = stack_browser(seg_img, cmap=greedy_colormap(seg_img))

>>> # Example 2 - Effect of cell sigma smoothing:
>>> from timagetk.algorithms.reconstruction import project_segmentation
>>> from timagetk.visu.mplt import grayscale_imshow
>>> seg_img_0, _, _ = watershed_segmentation(image, 30, cell_sigma=0)
>>> seg_img_1, _, _ = watershed_segmentation(image, 30, cell_sigma=1)
>>> seg_img_2, _, _ = watershed_segmentation(image, 30, cell_sigma=2)
>>> seg_img_3, _, _ = watershed_segmentation(image, 30, cell_sigma=3)
>>> z_proj = [project_segmentation(seg_img_0)]
>>> z_proj += [project_segmentation(seg_img_1)]
>>> z_proj += [project_segmentation(seg_img_2)]
>>> z_proj += [project_segmentation(seg_img_3)]
>>> grayscale_imshow(z_proj, title=["cell_sigma=0", "cell_sigma=1", "cell_sigma=2", "cell_sigma=3"], cmap='glasbey', val_range='auto')

>>> # Example 3 - Show volumes distribution (in real units)
>>> from timagetk import TissueImage3D
>>> seg_img, _, _ = watershed_segmentation(image, 15, min_vol=0)
>>> tissue = TissueImage3D(seg_img, background=1, not_a_label=0)
>>> volumes = tissue.cells.volume(real=True)
>>> volumes.pop(1)
>>> import matplotlib.pyplot as plt
>>> fig = plt.figure()
>>> plt.hist(volumes.values(), range=(0, 500), bins=100)  # crop the range for better output
>>> plt.xlabel("Volumes (µm³)")
>>> plt.ylabel("Frequency")
>>> plt.show()

Table of Contents

Tasks ¶

Multi-angle fusion ¶

Intensity image segmentation ¶

Clustering ¶

Table of Contents

Tasks¶

Multi-angle fusion¶

Intensity image segmentation¶

Clustering¶

Tasks ¶

Multi-angle fusion ¶

Intensity image segmentation ¶

Clustering ¶