cell segmentation algorithms

Cell segmentation algorithms

Brightfield image segmentation

This mode was developed to segment brightfield images captured at 2 um BELOW the focal plane (cells apear surrounded by a light halo).

Selecting this mode

To select this mode, set the SegMode property of the CyTracker object to 'brightfield':

CT.SegMode = 'brightfield';

Description

1. Median-filter image using a 3x3 kernel.
2. Background subtract the image using tophat subtraction: a. Morphologically open the image with a disk-shaped structuring element of radius 20 pixels. b. Subtract the original image from the opened image.
3. Get a value for intensity thresholding by fitting the image intensity histogram to a Gaussian (fits the background intensity). greyscaleThreshold = B * ThresholdLevel * C where B is the offset of the Gaussian and C is the FWHM.
4. Compute the mask: a. Initial mask is obtained by intensity thresholding the image. b. Holes in the mask are filled. c. Open the mask with a disk-shaped structuring element, radius 3 px. d. Close the mask with a disk-shaped structuring element, radius 3 px. e. Erode the mask with a single-pixel structuring element.
5. Label the cells: a. Compute the inverse distance-transform of the cell mask. b. Set the minimum depth to MaxCellMinDepth. c. Watershed to label cells. d. Final mask is thickened by 8 pixels.
6. Refine cell labels a. Remove labels which are below/above to specified min/max cell areas (CellAreaLim)

Fluorescence image segmentation

This mode was created for whole-cell fluorescence (i.e. GFP or mOrange rather than Cy5).

Selecting this mode

To select this mode, set the SegMode property of the CyTracker object to 'fluorescence':

CT.SegMode = 'fluorescence';

Description

1. The image is intensity-normalized (i.e. the intensity range is set to 0 to 1 where 1 is the maximum pixel value of the image).
2. Get a value for intensity thresholding by fitting the image intensity histogram to a Gaussian (fits the background intensity). The initial fitting parameters are determined by using the highest count of the histogram. greyscaleThreshold = B * ThresholdLevel * C where B is the offset of the Gaussian and C is the FWHM.
3. Compute the mask: a. Initial mask is obtained by intensity thresholding the image. b. Open the mask with a disk-shaped structuring element, radius 2 px. c. Regions touching image edges are cleared. d. Active contour algorithm is used to fit the initial masks to the cell image. e. Small regions (<100 pixels in area) are removed. f. The mask is openen using a disk-shaped structuring element, radius 2 px. g. Holes in the mask are filled.
4. Markers for each cell is created (for marker-based watershedding): a. The original cell image is median filtered using a 10x10 kernel. b. The regional maxima are found. c. The mask is tidied up using dilation (disk, 6px) and erosion disk, 3 px).
5. Refine the cell markers. The median of the mean intensities from each of the objects is calculated. Any marker whichfalls below 20% of the median value is removed.
6. Label the cells: a. The original fluorescence image is median filtered using a 4x4 kernel. a. Compute the distance-transform of the color-complemented median-filtered image. b. Minima are imposed on the distance-transform using the markers and mask (everything outside the masks are set to a minimum). c. Watershed to label cells. d. Final mask is opened using a disk-shaped structuring element of 6 pixels.
7. Refine cell labels a. Remove labels which are below/above to specified min/max cell areas (CellAreaLim)

Experimental (In progress)

The segmentation progresses similar to the 'brightfield' algorithm with two differences: (1) A step is introduced after the initial watershed (Step 5) to split cell masks which are too large, and (2) cylinders are drawn over the cell masks.

Selecting this mode

To select this mode, set the SegMode property of the CyTracker object to 'fluorescence':

CT.SegMode = 'experimental';

Description

(This follows Step 5 of brightfield)

Clump splitting

1. Identify cells which are "too large" using the absolute deviation around the median 1. The absolute deviation around the median (MAD) is calculated as follows: a. The median of all cell areas M is calculated. b. The absolute deviation is calculated as: D = 1.4826 * median(|Area - M|). The 1.4826 factor is applied, assuming that the distribution is normal. c. Cells are too large if their area is greater than M + 2 * D.
2. For each identified clump of cells: a. Calculate a new threshold using the 40-th percentile of the cell intensity. b. A mask is generated by thresholding the cell image. c. The mask is tidied up by opening, then closing using a disk-shaped structuring element of 3 px radius. d. The mask is eroded by 1-pixel. e. The distance transform is calculated using the inverse of this maks. Minimas less than the MaxCellMinDepth in the distance transform are surpressed. f. Watershed. g. Cells which are smaller than the minimum specified are removed.

Cylinder drawing

After splitting all large cell clumps, cylinders are drawn over the cell masks:

1. The centroid, length, axis, and orientation for each cell mask is measured using regionprops.
2. Cells which are too large/small are removed (using the limits set by CellAreaLim).
3. A capsule is drawn over each mask using these parameters.

The capsule algorithm is:

    function imgOut = drawCapsule(imgSize, props)
        %DRAWCAPSULE  Draws a capsule at the specified location
        %
        %  L = CyTracker.DRAWCAPSULE(I, RP) will draw capsules given
        %  the Centroid, MajorAxisLength, MinorAxisLength, and
        %  Orientation in RP. RP should be a struct, e.g. generated
        %  using REGIONPROPS.
        
        if numel(imgSize) ~= 2
            error('CyTracker:drawCapsule:InvalidImageSize',...
                'Image size should be a 1x2 vector.');                
        end
        
        imgOut = zeros(imgSize);
        
        xx = 1:size(imgOut, 2);
        yy = 1:size(imgOut, 1);
        
        [xx, yy] = meshgrid(xx,yy);
        
        for ii = 1:numel(props)
            
            center = props(ii).Centroid;
            theta = props(ii).Orientation/180 * pi;
            cellLen = props(ii).MajorAxisLength;
            cellWidth = props(ii).MinorAxisLength;
            
            rotX = (xx - center(1)) * cos(theta) - (yy - center(2)) * sin(theta);
            rotY = (xx - center(1)) * sin(theta) + (yy - center(2)) * cos(theta);
            
            %Plot the rectangle
            imgOut(abs(rotX) < (cellLen/2 - cellWidth/2) & abs(rotY) < cellWidth/2) = ii;
            
            % %Plot circles on either end
            imgOut(((rotX-(cellLen/2 - cellWidth/2)).^2 + rotY.^2) < (cellWidth/2)^2 ) = ii;
            imgOut(((rotX+(cellLen/2 - cellWidth/2)).^2 + rotY.^2) < (cellWidth/2)^2 ) = ii;
            
        end
        
    end

References

Back to user documentation