Re-written more specific question
For the details of my segmentation algorithm, see the original question below.
I have managed to improve the segmentation considerably by using a different strategy to generate markers. Since I know that my objects of interest are round and quite big, instead of finding peaks of the distance to background I simply use erosion in each plane to keep the "middle" of the objects as markers. And then I apply random walker in 3D. This works pretty well, except I get a strange effect of mislabelling small parts of one object close to the other object. I can remove these small parts by simply removing small or elongated objects in each plane, but there is no way to attach them to the real object, I think. And there is no way to fix the same issue when this small part is connected to the other big object. Here is an example:
These are two z-planes from the same 3D image, on the left there was a small part labelled as yellow that I removed by size/shape, on the right you can see how yellow object "invaded" the blue object. (The transparent red spot is the marker used for randomwalker)
So the question is, how do I prevent these "invasions"? Since I am segmenting a binary image, I think the beta in the randomwalker doesn't affect anything - I tried raising it to 10000 and don't think anything was different then...
I might try doing randomwalker in 2D for each slice, but it gets a bit tricky since the same object can get a different label in different slices (especially if there is only one object in one slice, but two in the next) and then it's difficult to merge them in 3D.
Original question:
This is the first question I am asking here, please let me know if there are any details I forgot to mention.
I need to segment a set of anisotropic 3D images - confocal images of DAPI staining of zygotes. I am using scikit-image. I have been struggling with it for a long time, trying to improve the success rate, but whatever I do, I might improve segmentation of some images, but segmentation of others gets worse...
I know that in each image in the end I have only two objects of interest (rarely one, but in that case I don't have to analyse the image) - which are nuclei. However, often they are close together and simple thresholding creates just one object where there should be two. And sometimes there is one or two small and much brighter objects in the image (polar bodies) which I need to exclude, but they interfere with thresholding since Otsu often considers my objects of interest background if these are present; I try to tackle them by simply clipping the brightest 1% of pixels in the image, and often it works.
My workflow currently involves first processing each plane in the z-stack separately: I do gaussian smoothing, Otsu thresholding, some filtering of small artefactual objects and misshapen holes (I know that only holes I want to keep are round-ish).
Here is an example of a relatively easy case where just thresholding works fine - transparent yellow-ish and blue-ish colours are overlaid on the two nuclei (and you can see a small very bright object in the 3-21 planes) which is luckily avoided using the described workflow):
But many images are more complicated (and I would prefer the same code to work for all images, and not tailor the segmentation method for each image, since it defeats the purpose of doing it automatically). I tried watershed or random walker segmentation using peaks of distance to background as markers (following http://www.scipy-lectures.org/packages/scikit-image/#marker-based-methods). This is tricky, since defining peaks can be done very differently. The distance to background can be calculated in 2D for each plane, or in 3D; number of peaks (or peaks per label) can be different; peaks can be defined using footprint of different size (and in 2D or in 3D). And, finally, the watershed or random-walker can be applied either to the whole image, or to each plane separately (however for randomwalker with each plane separately I think there is a bug that prevents me from using this approach https://github.com/scikit-image/scikit-image/issues/2783)
Playing with these parameters showed me that they can affect the segmentation quality a lot, and it seems ideal parameters for one image might be suboptimal for another; moreover, sometimes applying watershed or randomwalker even reduces the quality of initial segmentation. So I guess my question is: is there anything else I could try that would be more robust, or is there a smart way to choose the parameters for the methods I have tried?
The full code is quite long, I put the current version of it onto github: https://github.com/Phlya/ImageScripts/blob/master/zygote_segmentation_improved.py
