Oversplitting by watershed

Tue Nov 13 06:57:33 EST 2012

Hi,

just found a way to get my desired result: applying a gaussian filter to
the distance map allows be to adjust the number of local maxima found and
thereby controlling the sensitivity of the following watershed. Maybe not
the best option, but it serves the purpose. Here the code if you are
interested to check it yourself.

# packages needed to perform image processing and analysis
import numpy as np
import scipy as scp
import matplotlib.pyplot as plt
import matplotlib.image as mpimg
import scipy.ndimage as nd
import skimage
from skimage import io
from skimage.morphology import watershed, is_local_maximum
from skimage.segmentation import find_boundaries, visualize_boundaries
from skimage.color import gray2rgb

#read files jpeg file
image = mpimg.imread('c:\\test.jpg')
image_thresh = image > 140
labels = nd.label(image_thresh)[0]

distance = nd.distance_transform_edt(image_thresh)
#apply Gaussian filter to the distance map to merge several local maxima
into one
distance=nd.gaussian_filter(distance,3)
local_maxi = is_local_maximum(distance, labels=labels,
footprint=np.ones((9, 9)))
markers = nd.label(local_maxi)[0]
labelled_image = watershed(-distance, markers, mask=image_thresh)

#find outline of objects for plotting
boundaries = find_boundaries(labelled_image)
img_rgb = gray2rgb(image)
overlay = np.flipud(visualize_boundaries(img_rgb,boundaries))
imshow(overlay)

Cheers,

    Frank

On Tue, Nov 13, 2012 at 10:54 AM, Frank Pennekamp <
pennekampster at googlemail.com> wrote:

> Hi Tony,
>
> thanks for helping me out on this again. Your solution produces a nice
> segmentation of the image, but the particles that need to be split remain
> touching (the diving cell left of the big blob in the middle; the two cells
> in the lower left quarter that touch on their tips). I think it is the same
> result as just using the global threshold.
>
> I agree with you that the problem seem to be the markers. I have about
> three times more markers than actual objects, so that's not corresponding
> to the actual number of objects at all. On the other extreme, replacing the
> regional maxima with the centroid of the thresholded blobs is not splitting
> the touching objects, because there is only one centroid per object.
>
> I had some success in splitting objects with the watershed algorithm
> implemented in ImageJ, maybe there is a way of translating their approach
> into Python. Their description is the follwing:
>
> *Watershed segmentation* of the Euclidian distance map (EDM) is a way of
> automatically separating or cutting apart particles that touch (Watershed
> separation of a grayscale image is available via the Find Maxima...<http://imagejdocu.tudor.lu/doku.php?id=gui:process:binary#find_maxima>command). The Watershed command requires a binary image containing black
> particles on a white background. It first calculates the Euclidian distance
> map and finds the ultimate eroded points (UEPs). It then dilates each of
> the UEPs (the peaks or local maxima of the EDM) as far as possible - either
> until the edge of the particle is reached, or the edge of the region of
> another (growing) UEP. Watershed segmentation works best for smooth convex
> objects that don't overlap too much.
>
> [image: watershed example]<http://imagejdocu.tudor.lu/doku.php?id=gui:process:start>
> *Ultimate points:* generates the ultimate eroded points (UEPs) of the
> EDM. Requires a binary image as input. The UEPs represent the centers of
> particles that would be separated by segmentation. The UEP's gray value is
> equal to the radius of the inscribed circle of the corresponding particle.
> Use Process>Binary>Options<http://imagejdocu.tudor.lu/doku.php?id=gui:process:binary#options>to set the background color (black or white) and the output type.
>
> How could i get the ultimate eroded points in scikit image? There seems no
> function to do so for the moment, but may you have a suggestion how to
> tackle this problem?
>
> Many thanks in any case for your help already!
>
> Best,
>
>    Frank
>
>
> On Mon, Nov 12, 2012 at 10:56 PM, Tony Yu <tsyu80 at gmail.com> wrote:
>
>>
>>
>> On Mon, Nov 12, 2012 at 7:43 AM, Frank <pennekampster at googlemail.com>wrote:
>>
>>> Dear group,
>>>
>>> I have some issues with the watershed algorithm implemented in scikits
>>> image. I use a global threshold to segment cells from background, but some
>>> cells touch and I want them to be split. Watershed seems the appropriate
>>> way to deal with my problem, however my particles are split in too many
>>> pieces. Is there a way to adjust the sensitivity of the watershed method?
>>>
>>> Many thanks for any suggestion!
>>>
>>> The code that I use looks like below. An example image that I want to
>>> process can be downloaded here:
>>> https://dl.dropbox.com/u/10373933/test.jpg
>>>
>>> # packages needed to perform image processing and analysis
>>> import numpy as np
>>> import scipy as scp
>>> import matplotlib.pyplot as plt
>>> import matplotlib.image as mpimg
>>> import scipy.ndimage as nd
>>> import skimage
>>> from skimage import io
>>> from skimage.morphology import watershed, is_local_maximum
>>> from skimage.segmentation import find_boundaries, visualize_boundaries
>>> from skimage.color import gray2rgb
>>>
>>> #read files jpeg file
>>> image = mpimg.imread('c:\\test.jpg')
>>> image_thresh = image > 140
>>> labels = nd.label(image_thresh)[0]
>>> distance = nd.distance_transform_edt(image_thresh)
>>> local_maxi = is_local_maximum(distance, labels=labels,
>>> footprint=np.ones((9, 9)))
>>> markers = nd.label(local_maxi)[0]
>>> labelled_image = watershed(-distance, markers, mask=image_thresh)
>>>
>>> #find outline of objects for plotting
>>> boundaries = find_boundaries(labelled_image)
>>> img_rgb = gray2rgb(image)
>>> overlay = np.flipud(visualize_boundaries(img_rgb,boundaries))
>>> imshow(overlay)
>>
>>
>> Hi Frank,
>>
>> Actually, I don't think the issue is in the watershed segmentation.
>> Instead, I think the problem is in the marker specification: Using local
>> maxima creates too many marker points when a blob deviates greatly from a
>> circle. (BTW, does anyone know if there are any differences between
>> `is_local_maximum` and `peak_local_max`? Maybe the former should be
>> deprecated.)
>>
>> Using the centroids of blobs gives cleaner results. See slightly-modified
>> example below.
>>
>> Best,
>> -Tony
>>
>> # packages needed to perform image processing and analysis
>> import numpy as np
>> import matplotlib.pyplot as plt
>> import scipy.ndimage as nd
>>
>> from skimage import io
>> from skimage import measure
>> from skimage.morphology import watershed
>> from skimage.segmentation import find_boundaries, visualize_boundaries
>> from skimage.color import gray2rgb
>>
>> #read files jpeg file
>> image = io.imread('test.jpg')
>>
>> image_thresh = image > 140
>> labels = nd.label(image_thresh)[0]
>> distance = nd.distance_transform_edt(image_thresh)
>>
>> props = measure.regionprops(labels, ['Centroid'])
>> coords = np.array([np.round(p['Centroid']) for p in props], dtype=int)
>> # Create marker image where blob centroids are marked True
>> markers = np.zeros(image.shape, dtype=bool)
>> markers[tuple(np.transpose(coords))] = True
>>
>> labelled_image = watershed(-distance, markers, mask=image_thresh)
>>
>> #find outline of objects for plotting
>> boundaries = find_boundaries(labelled_image)
>> img_rgb = gray2rgb(image)
>> overlay = visualize_boundaries(img_rgb, boundaries, color=(1, 0, 0))
>>
>> plt.imshow(overlay)
>> plt.show()
>>
>> --
>>
>>
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scikit-image/attachments/20121113/b1461638/attachment.html>