0

I'm using the watershed algorithm to try and segment touching nuclei. A typical image may look like:enter image description here or this: enter image description here

I'm trying to apply the watershed algorithm with this code:

show(RGB_img)


%Convert to grayscale image
I = rgb2gray(RGB_img);

%Take structuring element of a disk of size 10, for the morphological transformations
%Attempt to subtract the background from the image: top hat is the
%subtraction of the open image from the original


%Morphological transformation to subtract background noise from the image
%Tophat is the subtraction of an opened image from the original. Remove all
%images smaller than the structuring element of 10
I1 = imtophat(I, strel('disk', 10));

%Increases contrast
I2 = imadjust(I1);
%show(I2,'contrast')
%Assume we have background and foreground and assess thresh as such 
level = graythresh(I2);
%Convert to binary image based on graythreshold
BW = im2bw(I2,level);
show(BW,'C');



BW = bwareaopen(BW,8);
show(BW,'C2');

BW = bwdist(BW) <= 1;
show(BW,'joined');
%Complement because we want image to be black and background white
C = ~BW;
%Use distance tranform to find nearest nonzero values from every pixel
D = -bwdist(C);

%Assign Minus infinity values to the values of C inside of the D image
%   Modify the image so that the background pixels and the extended maxima
%   pixels are forced to be the only local minima in the image (So you could
%   hypothetically fill in water on the image

D(C) = -Inf;

%Gets 0 for all watershed lines and integers for each object (basins)
L = watershed(D);
show(L,'L');

%Takes the labels and converts to an RGB (Using hot colormap)
fin = label2rgb(L,'hot','w');

% show(fin,'fin');
im = I;

%Superimpose ridgelines,L has all of them as 0 -> so mark these as 0(black)
im(L==0)=0;

clean_img = L;
show(clean_img)

After C = ~BW; the whole image goes dark. I believe this is because the image pixels are all -inf or some smaller negative number. This is there a way around this and if so what could I change in my code to get this algorithm working? I've experimented a ton and I don't really know what's happening. Any help would be great!

muZero
  • 948
  • 9
  • 22
  • what does the `show` command do? I don't recognise it. Is it a 2016a thing? – Tasos Papastylianou Aug 22 '16 at 16:45
  • @tasosPapastylianou show is a small function I wrote- it really only encapsulates `figure`, `imshow`, and adding a label to the figure. It should be the same display though as `imshow` – muZero Aug 22 '16 at 16:47

1 Answers1

2

The problem is with your show command. As you said in the comments this uses imshow under the hood. If you try imshow directly you'll see you also get a black image. However, if you call it with appropriate limits:

imshow(clean_img,[min(clean_img(:)), max(clean_img(:))])

you'll see everything you expect to see.

In general I usually prefer imagesc for that reason. imshow makes arbitrary judgements as to what range to represent, and I usually can't be bothered to keep up with it. I think in your case, your end image is uint16 so imshow chooses to represent the range [1, 65025]. Since all your pixel values are below 400, they look black to the naked eye for that range.

Tasos Papastylianou
  • 21,371
  • 2
  • 28
  • 57
  • 1
    You can also just do `imshow(clean_img,[]);`, which will automatically compute the min and max values and scale to `[0,1]` respectively. – rayryeng Aug 22 '16 at 17:19
  • Thanks a lot, this helped. However, I'm now using another algorithm in which `L = watershed(I_mod);` comes after: `I_eq_c = imcomplement(I_eq);` and `I_mod = imimposemin(I_eq_c, ~bw4 | mask_em);` After `L` the image's pixels are all set to 1 and is therefore all black again, any ideas? – muZero Aug 22 '16 at 17:28
  • Also is there a more general purpose method to segment these? I seem to be getting spotty results if results at all- even though the tutorial code used just as complicated images. – muZero Aug 22 '16 at 17:34
  • There's no single generally applicable segmentation method. There's still people doing PhDs on cell segmentation. Try `imopen` instead of `imtophat`, you'll get a slightly more relevant result in your case, perhaps you can take it from there. One approach is to crudely isolate the shapes first using a crude threshold, and only apply the watershed algorithm on 'joined' cells; you could identify those as being structures that fail some sort of "convexity" measure threshold (e.g. area / surface) – Tasos Papastylianou Aug 22 '16 at 18:18
  • thanks @rayryeng , I tend to avoid `imshow` so I keep forgetting that! :p – Tasos Papastylianou Aug 22 '16 at 18:19
  • @Sam also in general I'd avoid working on indexed integer images unless there's a particular reason to. I'd `mat2gray` the crap out of your input image first and work on that instead. This even allows you to do `I.^2` which separates lower to higher intensities in one fell swoop. – Tasos Papastylianou Aug 22 '16 at 18:24
  • have a look at the related posts to the right. one of them links to http://blogs.mathworks.com/steve/2006/06/02/cell-segmentation/ which may also be of some help – Tasos Papastylianou Aug 22 '16 at 18:31