Segmentation

The SIMA package implements multiple segmentation strategies, which share a common interface defined by the abstract class SegmentationStrategy. Any strategy with this interface can be used in conjunction with the ImagingDataset.segment() method to segment an ImagingDataset object.

ImagingDataset.segment(strategy, label=None, planes=None)

Segment an ImagingDataset to generate ROIs.

Parameters:
  • strategy (sima.segment.SegmentationStrategy) – The strategy for segmentation.
  • label (str, optional) – Label to be associated with the segmented set of ROIs.
  • planes (list of int) – List of the planes that are to be segmented.
Returns:

ROIs – The segmented regions of interest.

Return type:

sima.ROI.ROIList

Segmentation strategies

Note

Please consider contributing additional methods to the SIMA project.

class sima.segment.SegmentationStrategy

Abstract class implementing the interface for segmentation strategies.

This class can be subclassed to create a concrete segmentation strategy by implementing a _segment() method. As well, any existing segmentation strategy can be extended by adding a PostProcessingStep with the append() method.

_segment(dataset)

Implementation of the segmentation strategy.

This abstract method must be implemented by any subclass of SegmentationStrategy.

Parameters:dataset (sima.ImagingDataset) – The dataset to be segmented.
Returns:ROIs – The segmented regions of interest.
Return type:sima.ROI.ROIList
append(post_processing_step)

Add a post processing step.

Parameters:post_processing_step (PostProcessingStep) –
segment(dataset)

Apply the segmentation method to a dataset.

Parameters:dataset (ImagingDataset) – The dataset whose affinity matrix is being calculated. return self._segment(dataset)
Returns:rois – A list of sima.ROI ROI objects.
Return type:sima.ROI.ROIList

The specific segmentation strategies that have been implemented are documented below. Once initialized as documented, these strategies all share the above interface.

Plane-Wise Segmentation

class sima.segment.PlaneWiseSegmentation(plane_strategy)

Bases: sima.segment.segment.SegmentationStrategy

Segmentation approach with each plane segmented separately.

Parameters:plane_strategy (SegmentationStrategy or list of SegmentationStrategy) – The strategies to be applied to each plane.

Examples

One use is to apply the same segmentation method to multiple layers. For example, the user may wish to apply the PlaneCA1PC strategy separately to multiple well-separated planes in which distinct cell bodies are imaged:

>>> from sima.segment import PlaneWiseSegmentation, PlaneCA1PC
>>> layer_strategy = PlaneCA1PC()
>>> strategy = PlaneWiseSegmentation(layer_strategy)

Alternatively, one may wish to use different segmentation strategies on each plane. For example, to segment one plane of dendrites with stICA and one plane of cell bodies with the CA1PC strategy, a PlaneWiseSegmentation strategy can be created as follows:

>>> from sima.segment import PlaneWiseSegmentation, PlaneCA1PC, STICA
>>> strategy = PlaneWiseSegmentation([STICA(), PlaneCA1PC()])

Spatialtemporal Independent Component Analysis

class sima.segment.STICA(channel=0, mu=0.01, components=75, verbose=False)

Bases: sima.segment.segment.SegmentationStrategy

Segmentation using spatiotemporal independent component analysis (stICA).

Parameters:
  • channel (string or int, optional) – Channel containing the signal to be segmented, either an integer index or a channel name. Default: 0.
  • mu (float, optional) – Weighting parameter for the trade off between spatial and temporal information. Must be between 0 and 1. Low values give higher weight to temporal information. Default: 0.01
  • components (int or list, optional) – Number of principal components to use. If list is given, then use only the principal components indexed by the list Default: 75
  • verbose (bool, optional) – Whether to print progress updates.

Notes

Spatiotemporal (stICA) [1] is a procedure which applies ICA to extracted PCA components in a process that takes into consideration both the spatial and temporal character of these components. This method has been used to segment calcium imaging data [2], and can be used to segment cell bodies, dendrites, and axons.

In order to implement spatio and temporal ICA, temporal components from PCA are concatenated to the spatial ones. The following spatiotemporal variable \(y_i\) and the resulting ICA components \(z_i\) are defined by:

\[\begin{split}y_i &= \begin{cases} \mu U_{ki} & i \leq N_x \\ (1-\mu)V_{ki} & N_x < i \leq (N_x+N_t) \end{cases} \\ z_i^{k} &= \sum_j \mathbf{W}_{i,j}^{(n)} y^{(j)},\end{split}\]

where \(U\) corresponds to the spatio PCA component matrix with dimensions \(N_x\), pixels, by \(k\) principal components and \(V\) corresponds to the \(N_t\), time frames, by \(k\) temporal PCA component matrix. \(\mu\) is a weighting parameter to balance the trade-off between the spatial and temporal information with low values of \(\mu\) giving higher weight to the signals temporal components. ICA is performed on \(y_i\) to extract the independent components \(z_i\).

References

[1]Stone JV, Porrill J, Porter NR, Wilkinson ID. Spatiotemporal independent component analysis of event-related fMRI data using skewed probability density functions. Neuroimage. 2002 Feb;15(2):407-21.
[2]Mukamel EA, Nimmerjahn A, Schnitzer MJ. Automated analysis of cellular signals from large-scale calcium imaging data. Neuron. 2009 Sep 24;63(6):747-60.

Warning

In version 1.0, this method currently only works on datasets with a single plane, or in conjunction with sima.segment.PlaneWiseSegmentation.

Normalized cuts

class sima.segment.PlaneNormalizedCuts(affinity_method=None, cut_max_pen=0.01, cut_min_size=40, cut_max_size=200)

Bases: sima.segment.segment.SegmentationStrategy

Segment image by iteratively performing normalized cuts.

Parameters:
  • affinity_method (AffinityMatrixMethod) – The method used to calculate the affinity matrix.
  • max_pen (float) – Iterative cutting will continue as long as the cut cost is less than max_pen.
  • cut_max_size (cut_min_size,) – Regardless of the cut cost, iterative cutting will not be performed on regions with fewer pixels than min_size and will always be performed on regions larger than max_size.

Notes

The normalized cut procedure [3] is iteratively applied, first to the entire image, and then to each cut made from the previous application of the procedure.

References

[3]Jianbo Shi and Jitendra Malik. Normalized Cuts and Image Segmentation. IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. 22, NO. 8, AUGUST 2000.

Warning

In version 1.0, this method currently only works on datasets with a single plane, or in conjunction with sima.segment.PlaneWiseSegmentation.

Affinity Matrix Methods

class sima.segment.AffinityMatrixMethod

Method for calculating the affinity matrix

calculate(dataset)

Calculate the affinity matrix for a dataset.

Parameters:dataset (sima.ImagingDataset) – The dataset for which the affinity matrix is to be calculated.
Returns:affinities – The affinities between the image pixels.
Return type:scipy.sparse.coo_matrix
class sima.segment.BasicAffinityMatrix(channel=0, max_dist=None, spatial_decay=None, num_pcs=75, verbose=False)

Bases: sima.segment.normcut.AffinityMatrixMethod

Return a sparse affinity matrix for use with normalized cuts.

The affinity \(A_{ij}\) between each pair of pixels \(i,j\) is a function of the correlation \(c_{i,j}\) of the pixel-intensity time series, and the relative locations (\(\mathbf X_i,\mathbf X_j\)) of the pixels:

\[A_{ij} = e^{k_cc_{ij}} \cdot e^{-\frac{|\mathbf X_i-\mathbf X_j|_2^2}{\sigma_{\mathbf X}^2}},\]

with \(k_c\) and \(\sigma_{\mathbf X}^2\) being automatically determined constants.

Parameters:
  • channel (string or int, optional) – Channel containing the signal to be processed, either an integer index or a channel name. Default: 0.
  • max_dist (tuple of int, optional) – Defaults to (2, 2).
  • spatial_decay (tuple of int, optional) – Defaults to (2, 2).
  • num_pcs (int, optional) – The number of principal component to use. Default: 75.
  • verbose (bool, optional) – Whether to print progress status. Default: False.
class sima.segment.AffinityMatrixCA1PC(channel=0, max_dist=None, spatial_decay=None, num_pcs=75, x_diameter=10, y_diameter=10, verbose=False)

Bases: sima.segment.normcut.BasicAffinityMatrix

CA1 pyramidal cells

class sima.segment.PlaneCA1PC(channel=0, num_pcs=75, max_dist=None, spatial_decay=None, cut_max_pen=0.01, cut_min_size=40, cut_max_size=200, x_diameter=10, y_diameter=10, circularity_threhold=0.5, min_roi_size=20, min_cut_size=30, verbose=False)

Bases: sima.segment.segment.SegmentationStrategy

Segmentation method designed for finding CA1 pyramidal cell somata.

Parameters:
  • channel (int, optional) – Channel containing the signal to be segmented, either an integer index or a channel name. Default: 0.
  • num_pcs (int, optional) – The number of principle components to be used in the calculations. Default: 75.
  • max_dist (tuple of int, optional) – Defaults to (2, 2).
  • spatial_decay (tuple of int, optional) – Defaults to (2, 2).
  • max_pen (float) – Iterative cutting will continue as long as the cut cost is less than max_pen.
  • cut_max_size (cut_min_size,) – Regardless of the cut cost, iterative cutting will not be performed on regions with fewer pixels than min_size and will always be performed on regions larger than max_size.
  • circularity_threhold (float) – ROIs with circularity below threshold are discarded. Default: 0.5.
  • min_roi_size (int, optional) – ROIs with fewer than min_roi_size pixels are discarded. Default: 20.
  • min_cut_size (int, optional) – No ROIs are made from cuts with fewer than min_cut_size pixels. Default: 30.
  • x_diameter (int, optional) – The estimated x-diameter of the nuclei in pixels. Default: 8
  • y_diameter (int, optional) – The estimated x-diameter of the nuclei in pixels. Default: 8

Notes

This method begins with the normalized cuts procedure. The affinities used for normalized cuts also depend on the time-averaged image, which is used to increase the affinity between pixels within the same dark patch (putative nucleus).

Following the normalized cut procedure, the algorithm attempts to identify a nucleus within each cut.

  • Otsu thresholding of the time-averaged image after processing with CLAHE and unsharp mask procedures.
  • Binary opening and closing of the resulting ROIs.
  • Rejection of ROIs not passing the circularity and size requirements.

See also

sima.segment.normcut

Warning

In version 1.0, this method currently only works on datasets with a single plane, or in conjunction with sima.segment.PlaneWiseSegmentation.

Post-Processing Steps

Any number of post-processing steps can be added to a segmentation method using the SegmentationStrategy.append() method. These appended post-processing steps must have the interface defined by PostProcessingStep below. The appended steps can be selected from those documented in this section, or can be created by the user by subclassing any of the classes listed below.

class sima.segment.PostProcessingStep

Abstract class representing the interface for post processing steps that can be appended to a segmentation method to modify the the segmented ROIs.

Examples

To apply a binary opening to all ROIs, a subclass of PostProcessingStep can be created as follows:

>>> from scipy import ndimage
>>> import sima.segment
>>> class BinaryOpening(sima.segment.PostProcessingStep):
...     def apply(self, rois, dataset=None):
...         for r in rois:
...             r.mask = ndimage.binary_opening(r.mask)
...         return rois

We can then append an instance of this class to any segmentation strategy.

>>> strategy = sima.segment.STICA()
>>> strategy.append(BinaryOpening())
apply(rois, dataset=None)

Apply the post-processing step to rois from a dataset.

Parameters:
Returns:

ROIs – The post-processed ROIs.

Return type:

sima.ROI.ROIList

class sima.segment.ROIFilter(func)

Bases: sima.segment.segment.PostProcessingStep

Post-processing step for generic filtering of ROIs.

ROIs produced by the segmentation are filtered to retain only the ROIs that cause the specified function to evaluate to be True.

Parameters:func (function) – A boolean-valued function taking arguments (rois, dataset) that is used to filter the ROIs.

Example

To select ROIs with a size (i.e. number of non-zero pixels) of at least 20 pixels and no more than 50 pixels, the following ROIFilter can be created:

>>> from sima.segment import ROIFilter
>>> size_filter = ROIFilter(lambda roi: roi.size >= 20 and roi.size <= 50)
class sima.segment.CircularityFilter(circularity_threhold=0.5)

Bases: sima.segment.segment.ROIFilter

Post-processing step to filter ROIs based on circularity.

Parameters:circularity_threhold (float, optional) – ROIs with circularity below threshold are discarded. Default: 0.5. Range: 0 to 1.
class sima.segment.MergeOverlapping(threshold)

Bases: sima.segment.segment.PostProcessingStep

Post-processing step to merge overlapping ROIs.

Parameters:threshold (float) – Minimum percent of the smaller ROIs total area which must be covered in order for the ROIs to be evaluated as overlapping.
apply(rois, dataset=None)

Remove overlapping ROIs

Parameters:
  • rois (list) – list of sima.ROI ROIs
  • percent_overlap (float) – percent of the smaller ROIs total area which must be covered in order for the ROIs to be evaluated as overlapping
Returns:

rois – A list of sima.ROI ROI objects with the overlapping ROIs combined

Return type:

list

class sima.segment.SmoothROIBoundaries(tolerance=0.5, min_verts=8, n_processes=1)

Bases: sima.segment.segment.PostProcessingStep

Post-processing step to smooth ROI boundaries.

Reduces the number of points in the ROI polygons.

Parameters:
  • tolerance (float) – Maximum distance from original points of polygon to approximated polygonal chain. If tolerance is 0, the original roi coordinates array are returned. (See skimage.measure.approximate_polygon). Default is 0.5.
  • min_verts (int) – Minimum number of verticies an ROI polygon can have prior to attemting smoothing. Defaults is 8.
  • n_processes (int, optional) – Number of processes to farm out the roi smoothing across. Should be at least 1 and at most one less than the number of CPUs in the computer. Defaults to 1.
class sima.segment.SparseROIsFromMasks(min_size=50, static_threshold=0.5, smooth_size=4, sign_split=True, n_processes=1)

Bases: sima.segment.segment.PostProcessingStep

PostProcessingStep to extract sparse ROIs from masks.

Parameters:
  • min_size (int, optional) – The minimum size in number of pixels that an ROI can be. Default: 50.
  • static_threshold (float, optional) – threshold on the static allowable in an ICA components, eliminating high scoring components speeds the ROI extraction and may improve the results. Default: 0.5
  • smooth_size (int, optional) – number of iterations of static removal and Gaussian blur to perform on each stICA component. 0 provides no Gaussian blur, larger values produce stICA components with less static but the ROIs lose definition. Default: 4
  • sign_split (bool, optional) – Whether to split each mask into its positive and negative components.
  • n_processes (int, optional) – Number of processes to farm out the roi smoothing across. Should be at least 1 and at most one less than the number of CPUs in the computer. Defaults to 1.