MATLAB: Trying to do retinal vessel segmentation with mean-c thresholding

eyefundusMATLABmean-c thresholdingophthalmologyretinaretinal blood vesselsegmentationthresholding

Hello, I am trying to do a segmentation of retinal blood vessels with thresholding-based methods. I want to replicate this paper.
https://www.sciencedirect.com/science/article/pii/S2314728817300272
The pre-processing steps listed there are getting the green channel image, CLAHE, and median filter.
The steps of the segmentation is below.
a. A mean filter of window size N × N is selected.
b. The enhanced image is convolved with the mean.
c. A difference image is obtain by subtracting the convolved image from enhanced image.
d. The difference image is thresholded with the constant value C.
e. The complement of the thresholded image is calculated
Here is my code,
%objective : retinal blood vessel segmentation with thresholding
%%Pre-Processing%%

I=imread('r23_training.tif');
subplot(2,2,1),imshow(I);
title('a');
%green channel image
greenI=I(:,:,2);
subplot(2,2,2),imshow(greenI);
title('b');
%CLAHE
clahed = adapthisteq(greenI);
subplot(2,2,3),imshow(clahed);
title('c');
%median filter
medfiltered=medfilt2(clahed);
subplot(2,2,4),imshow(medfiltered);
title('d');
%%Pre-Processing%%
%%Segmentation%%
%create a 13x13 mean filter
meanfilter=[1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1]/169;
%convolve the enhanced image with mean filter
meanfiltered=uint8(conv2(medfiltered,meanfilter,'same'));
figure,imshow(meanfiltered);
title('meanfiltered');
diff=medfiltered-meanfiltered;
figure,imshow(diff);
title('difference');
%threshold the image with C = 0.042
thresholded = imbinarize(diff,0.042);
%complement the thresholded image
complemented = imcomplement(thresholded);
figure,imshow(complemented);
I wanted to have something like this,
But this is the result of my code.
Can anyone perhaps tell me what I did wrong? Thank you very much.
Also, below is the image (it's from DRIVE database).

Best Answer

Hi
You are converting a double image to uint8 format due to which you are unable to visualize the segmentation. Keep in double format. For example, run the below code:
%objective : retinal blood vessel segmentation with thresholding
%%Pre-Processing%%

I=imread('image.png');
subplot(2,2,1),imshow(I);
title('a');
%green channel image
greenI=I(:,:,2);
subplot(2,2,2),imshow(greenI);
title('b');
%CLAHE
clahed = adapthisteq(greenI);
subplot(2,2,3),imshow(clahed);
title('c');
%median filter
medfiltered=medfilt2(clahed);
subplot(2,2,4),imshow(medfiltered);
title('d');
%%Pre-Processing%%
%%Segmentation%%
%create a 13x13 mean filter
meanfilter=[1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1;1 1 1 1 1 1 1 1 1 1 1 1 1]/169;
%convolve the enhanced image with mean filter
meanfiltered=(conv2(medfiltered,meanfilter,'same'));
figure,imshow(meanfiltered);
title('meanfiltered');
diff=double(medfiltered)-meanfiltered;
figure,imshow(diff);
title('difference');
%threshold the image with C = 0.042
thresholded = imbinarize(diff,-15);  %%%%% check will other threshold
%complement the thresholded image
complemented = imcomplement(thresholded);
figure,imshow(complemented);
You will get the image as given below. After this, try to apply some morphological operation to remove the smaller component. Also, check with the algorithm if it is correct or not. Hope it will help!
Related Question