[Chapter 26 medical impact segmentation system based on minimum error method and region growth -- matlab deep learning practical GUI project]

Chest film segmentation system based on minimum error method –matlab Deep learning practice GUI project

Minimum error threshold segmentation method

“ Minimum error method the so-called image segmentation is based on the gray level 、 colour 、 Spatial texture 、 Geometric features divide the image into several disjoint regions , Make these features in the same area , Showing consistency or similarity , There are obvious differences between different regions .
It is realized according to the probability distribution density of background and target pixels in the image , The idea is to find a threshold , It is divided according to the threshold , Calculate the probability that the target point is misclassified into the background and the probability that the background point is misclassified into the target point , The total error partition probability is obtained . When the total error partition probability is the smallest , The desired optimal threshold is obtained .

Image preprocessing

Image enhancement results

Minimum error method segmentation results

The principle of region growing algorithm is to find the largest connected set of similar pixels in the image . Set pixels with similar gray values , Then, the pixels around the pixel are detected in turn 8 Neighborhood pixels , Add points with similar gray values to the previous set , If not, skip this point . Then browse other pixels in the image in turn , Until all pixels in the image are checked , The region thus obtained is the new region obtained according to the region growing algorithm .
Reference source :
Region growing image processing

Regional growth processing source code
The main function

clc; clear all; close all;
I = imread(fullfile(pwd, 'images/test.jpg'));
X = imadjust(I, [0.2 0.8], [0 1]);
%  Threshold segmentation 
bw = im2bw(X, graythresh(X));
[r, c] = find(bw);
rect = [min(c) min(r) max(c)-min(c) max(r)-min(r)];
Xt = imcrop(X, rect);
%  Automatically get seed points 
seed_point = round([size(Xt, 2)*0.15+rect(2) size(Xt, 1)*0.4+rect(1)]);
%  Region growth segmentation 
X = im2double(im2uint8(mat2gray(X)));
X(1:rect(2), :) = 0;
X(:, 1:rect(1)) = 0;
X(rect(2)+rect(4):end, :) = 0;
X(:, rect(1)+rect(3):end) = 0;
[J, seed_point, ts] = Regiongrowing(X, seed_point);
figure(1);
subplot(1, 2, 1); imshow(I, []);
hold on;
plot(seed_point(1), seed_point(2), 'ro', 'MarkerSize', 10, 'MarkerFaceColor', 'r');
title(' Seed points are automatically selected ');
hold off;
subplot(1, 2, 2); imshow(J, []); title(' Region growth image ');
%  Morphological post-processing 
bw = imfill(J, 'holes');
bw = imopen(bw, strel('disk', 5));
%  Edge extraction 
ed = bwboundaries(bw);
figure;
subplot(1, 2, 1); imshow(bw, []); title(' Morphological post-processing image ');
subplot(1, 2, 2); imshow(I); 
hold on;
for k = 1 : length(ed)
    %  edge 
    boundary = ed{
    k};
    plot(boundary(:,2), boundary(:,1), 'g', 'LineWidth', 2);
end
hold off;
title(' Edge marker image ');

Custom pick seed point location
The main function ：
main2.m

clc; clear all; close all;
I = imread(fullfile(pwd, 'images/test.jpg'));
X = imadjust(I, [0.2 0.8], [0 1]);
%  Region growth segmentation 
X = im2double(im2uint8(mat2gray(X)));
[J, seed_point, ts] = Regiongrowing(X);
figure(1);
subplot(1, 2, 1); imshow(I, []);
hold on;
plot(seed_point(1), seed_point(2), 'ro', 'MarkerSize', 10, 'MarkerFaceColor', 'r');
title(' Seed point selection ');
hold off;
subplot(1, 2, 2); imshow(J, []); title(' Region growth segmentation results ');

Morphological denoising .

Edge feature extraction

The main function :
main.m

clc; clear all; close all;
warning off all;
%  Read images 
filename = fullfile(pwd, 'images/test.jpg');
Img = imread(filename);
%  Graying 
if ndims(Img) == 3
    I = rgb2gray(Img);
else
    I = Img;
end
%  Direct binarization 
bw_direct = im2bw(I);
figure; imshow(bw_direct); title(' Direct binary segmentation ');
%  Circle the air in the gastric area 
c = [1524 1390 1454 1548 1652 1738 1725 1673 1524];
r = [1756 1909 2037 2055 1997 1863 1824 1787 1756];
bw_poly = roipoly(bw_direct, c, r);
figure;
imshow(I, []);
hold on;
plot(c, r, 'r-', 'LineWidth', 2);
hold off;
title(' Gastric regional air selection ');
%  Set the intragastric air to 255
J = I;
J(bw_poly) = 255;
%  Image enhancement 
J = mat2gray(J);
J = imadjust(J, [0.532 0.72], [0 1]);
J = im2uint8(mat2gray(J));
figure; imshow(J, []); title(' Image enhancement processing ');
%  Histogram statistics 
[counts, gray_style] = imhist(J);
%  Brightness level 
gray_level = length(gray_style);
%  Calculate the gray probability 
gray_probability  = counts ./ sum(counts);
%  Statistical pixel mean 
gray_mean = gray_style' * gray_probability;
%  initialization 
gray_vector = zeros(gray_level, 1);
w = gray_probability(1);
mean_k = 0;
gray_vector(1) = realmax;
ks = gray_level-1;
for k = 1 : ks
    %  Iterative calculation 
    w = w + gray_probability(k+1);
    mean_k = mean_k + k * gray_probability(k+1);
    %  Judge whether convergence or not 
    if (w < eps) || (w > 1-eps)
        gray_vector(k+1) = realmax;
    else
        %  Calculate the mean 
        mean_k1 = mean_k / w;
        mean_k2 = (gray_mean-mean_k) / (1-w);
        %  Calculate variance 
        var_k1 = (((0 : k)'-mean_k1).^2)' * gray_probability(1 : k+1);
        var_k1 = var_k1 / w;
        var_k2 = (((k+1 : ks)'-mean_k2).^2)' * gray_probability(k+2 : ks+1);
        var_k2 = var_k2 / (1-w);
        %  Calculate the objective function 
        if var_k1 > eps && var_k2 > eps
            gray_vector(k+1) = 1+w * log(var_k1)+(1-w) * log(var_k2)-2*w*log(w)-2*(1-w)*log(1-w);
        else
            gray_vector(k+1) = realmax;
        end
    end
end
%  Extreme value statistics 
min_gray_index = find(gray_vector == min(gray_vector));
min_gray_index = mean(min_gray_index);
%  Calculate threshold 
threshold_kittler = (min_gray_index-1)/ks;
%  Threshold segmentation 
bw__kittler = im2bw(J, threshold_kittler);
%  Show 
figure; imshow(bw__kittler, []); title(' Minimum error segmentation ');
%  Morphological post-processing 
bw_temp = bw__kittler;
%  Reverse color 
bw_temp = ~bw_temp;
%  Fill the hole 
bw_temp = imfill(bw_temp, 'holes');
%  Denoise 
bw_temp = imclose(bw_temp, strel('disk', 5));
bw_temp = imclearborder(bw_temp);
%  Area markers 
[L, num] = bwlabel(bw_temp);
%  Area properties 
stats = regionprops(L);
Ar = cat(1, stats.Area);
%  Extract the target and clean it up 
[Ar, ind] = sort(Ar, 'descend');
bw_temp(L ~= ind(1) & L ~= ind(2)) = 0;
%  Denoise 
bw_temp = imclose(bw_temp, strel('disk',20));
bw_temp = imfill(bw_temp, 'holes');
figure;
subplot(1, 2, 1); imshow(bw__kittler, []); title(' Binary image to be processed ');
subplot(1, 2, 2); imshow(bw_temp, []); title(' Morphological post-processing image ');
%  Extract lung margin 
ed = bwboundaries(bw_temp);
%  Show 
figure;
subplot(2, 2, 1); imshow(I, []); title(' Original image ');
subplot(2, 2, 2); imshow(J, []); title(' Enhance the image ');
subplot(2, 2, 3); imshow(bw_temp, []); title(' Binary image ');
subplot(2, 2, 4); imshow(I, []); hold on;
for k = 1 : length(ed)
    %  edge 
    boundary = ed{
    k};
    plot(boundary(:,2), boundary(:,1), 'g', 'LineWidth', 2);
end
title(' The lung margin shows signs ');
figure;
subplot(1, 2, 1); imshow(bw_temp, []); title(' Binary image ');
subplot(1, 2, 2); imshow(I, []); hold on;
for k = 1 : length(ed)
    %  edge 
    boundary = ed{
    k};
    plot(boundary(:,2), boundary(:,1), 'g', 'LineWidth', 2);
end
title(' The lung margin shows signs ');

References

% MAINFORM MATLAB code for MainForm.fig
% MAINFORM, by itself, creates a new MAINFORM or raises the existing
% singleton*.
%
% H = MAINFORM returns the handle to a new MAINFORM or the handle to
% the existing singleton*.
%
% MAINFORM(‘CALLBACK’,hObject,eventData,handles,…) calls the local
% function named CALLBACK in MAINFORM.M with the given input arguments.
%
% MAINFORM(‘Property’,‘Value’,…) creates a new MAINFORM or raises the
% existing singleton*. Starting from the left, property value pairs are
% applied to the GUI before MainForm_OpeningFcn gets called. An
% unrecognized property name or invalid value makes property application
% stop. All inputs are passed to MainForm_OpeningFcn via varargin.
%
% *See GUI Options on GUIDE’s Tools menu. Choose “GUI allows only one
% instance to run (singleton)”.
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help MainForm

Threshold segmentation test function ：

clc; clear all; close all;
I = imread('IMG00016.jpg');
bw = im2bw(I, graythresh(I));
figure; 
subplot(1, 2, 1); imshow(I, []); title(' Original image ');
subplot(1, 2, 2); imshow(bw, []); title(' Threshold segmentation image ');

Custom region growth function :

function [J, seed_point, ts] = Regiongrowing(I, seed_point)
%  Statistics take time 
t1 = cputime;
%  Parameter detection 
if nargin < 2
    %  Display and select the seed point 
    figure; imshow(I,[]);  hold on;
    seed_point = ginput(1);   
    plot(seed_point(1), seed_point(2), 'ro', 'MarkerSize', 10, 'MarkerFaceColor', 'r');
    title(' Seed point selection ');
    hold off;
end
%  Variable initialization 
seed_point = round(seed_point);
x = seed_point(2);
y = seed_point(1);
I = double(I);
rc = size(I);
J = zeros(rc(1), rc(2));
%  Parameter initialization 
seed_pixel = I(x,y);
seed_count = 1;
pixel_free = rc(1)*rc(2); 
pixel_index = 0;
pixel_list = zeros(pixel_free, 3);
pixel_similarity_min = 0;
pixel_similarity_limit = 0.1;
%  Neighborhood 
neighbor_index = [-1 0;
        1 0;
        0 -1;
        0 1];
%  Loop processing 
while pixel_similarity_min < pixel_similarity_limit && seed_count < rc(1)*rc(2)   
    %  Add neighborhood points 
    for k = 1 : size(neighbor_index, 1)
        %  Calculate adjacent positions 
        xk = x + neighbor_index(k, 1); 
        yk = y + neighbor_index(k, 2);
        %  Region growing 
        if xk>=1 && yk>=1 && xk<=rc(1) && yk<=rc(2) && J(xk,yk) == 0
            %  Meet the conditions 
            pixel_index = pixel_index+1;
            pixel_list(pixel_index,:) = [xk yk I(xk,yk)]; 
            %  Update status 
            J(xk, yk) = 1;
        end
    end
    %  Update space 
    if pixel_index+10 > pixel_free 
        pixel_free = pixel_free+pixel_free;
        pixel_list(pixel_index+1:pixel_free,:) = 0;
    end
    %  Statistical iteration 
    pixel_similarity = abs(pixel_list(1:pixel_index,3) - seed_pixel);
    [pixel_similarity_min, index] = min(pixel_similarity);
    %  Update status 
    J(x,y) = 1; 
    seed_count = seed_count+1;
    seed_pixel = (seed_pixel*seed_count + pixel_list(index,3))/(seed_count+1);
    %  Storage location 
    x = pixel_list(index,1); 
    y = pixel_list(index,2);
    pixel_list(index,:) = pixel_list(pixel_index,:);
    pixel_index = pixel_index-1;
end
%  Return results 
J = mat2gray(J);
J = im2bw(J, graythresh(J));
%  Statistics take time 
t2 = cputime;
ts = t2 - t1;

The source code of this article corresponds to

Chest film segmentation system based on minimum error method
Download link –> Portal 1

Liver image segmentation system based on region growing
Download link –> Portal