如何有效地为这个显微图像创建BW蒙版？-Java 学习之路

所以有些背景 . 我的任务是编写一个matlab程序来计算可见光显微图像中酵母细胞的数量 . 为此，我认为第一步将是细胞分割 . 在我得到真实的实验图像集之前，我开发了一种算法，使用 watershed 测试图像集 . 看起来像这样：

Original Image

分水岭的第一步是为细胞生成BW掩模 . 然后我将生成一个bwdist图像，该图像具有从BW蒙版生成的强制局部最小值 . 有了这个，我可以很容易地产生分水岭 .

BW mask

local minimum mask generated from BW mask

enter image description here

正如您所看到的，我的算法依赖于成功生成BW掩码 . 因为我需要从中生成bwdist图像和标记 . 最初，我按照以下步骤生成BW掩码：

生成图像的局部标准差sdImage = stdfilt（grayImage，ones（9））

std filter

使用BW阈值生成初始BW掩码binaryImage = sdImage <8;

initial BW filter

使用imclearborder来清除背景 . 使用其他一些代码在边框上添加单元格 .

final BW mask

Background finished. Here is my problem

但今天我收到了新的真实数据集 . 图像分辨率小得多，光照条件与测试图像集不同 . 颜色深度也小得多 . 这些使我的算法无用 . 就这个：

New image set

使用stdfilt无法生成良好的干净图像 . 而是生成这样的东西（注意：我已经调整了stdfilt函数和BW阈值的参数，以下是我能得到的最好结果）：

new stdfilt result

正如您所看到的，细胞中心有一些光像素，不需要比膜更暗 . 哪个引导bw阈值生成这样的东西：

new bw thresholding

bw阈值处理后的新bw图像具有不完整的膜或分段的细胞中心，并使它们不适合其他步骤 .

我最近才开始进行图像处理，不知道该怎么办 . 如果您有任何想法请帮助我！谢谢！

为了您的方便，我附上了一个来自Dropbox的链接subset of the images

1 回答

我认为你的方法存在根本问题 . 您的算法使用 stdfilt 以便对图像进行二值化 . 但这基本上意味着你假设背景和单元内存在低位"texture" . 这适用于您的第一张图片 . 但是，在第二个图像中，单元格内部存在"texture"，因此该假设被破坏 .

我认为一个更强的假设是每个单元格周围都有一个“环”（对于你发布的两个图像都有效） . 所以我采取了检测此环的方法 .

所以我的方法基本上是：

检测这些环（我使用'log'过滤器，然后根据正值进行二值化 . 但是，这导致了很多"chatter"
尝试通过过滤非常小和非常大的区域来删除一些"chatter"
现在，填写这些戒指 . 但是，在单元格之间仍然存在一些"chatter"和填充区域
同样，删除小区域和大区域，但由于单元格已填充，请增加可接受范围的边界 .
仍然存在一些不良区域，大多数不良区域将成为细胞之间的区域 . 通过观察区域边界周围的曲率可以检测细胞之间的区域 . 它们很多，在数学上表示为具有负曲率的边界的大部分 . 此外，要移除"chatter"的其余部分，这些区域的边界曲率将具有较大的标准偏差，因此也要删除具有较大标准偏差的边界 .

总的来说，最困难的部分是去除细胞之间的区域和“颤动”而不去除实际的细胞 .

无论如何，这里是代码（注意有很多启发式，而且它非常粗糙，并且基于旧项目，homeworks和stackoverflow答案的代码，所以它肯定远未完成）：

cell = im2double(imread('cell1.png'));
if (size(cell,3) == 3) 
    cell = rgb2gray(cell);
end

figure(1), subplot(3,2,1)
imshow(cell,[]);

% Detect edges
hw = 5;
cell_filt = imfilter(cell, fspecial('log',2*hw+1,1));

subplot(3,2,2)
imshow(cell_filt,[]);

% First remove hw and filter out noncell hws
mask = cell_filt > 0;
hw = 5;
mask = mask(hw:end-hw-1,hw:end-hw-1);

subplot(3,2,3)
imshow(mask,[]);

rp = regionprops(mask, 'PixelIdxList', 'Area');
rp = rp(vertcat(rp.Area) > 50 & vertcat(rp.Area) < 2000);

mask(:) = false;
mask(vertcat(rp.PixelIdxList)) = true;

subplot(3,2,4)
imshow(mask,[]);

% Now fill objects
mask1 = true(size(mask) + hw);
mask1(hw+1:end, hw+1:end) = mask;
mask1 = imfill(mask1,'holes');
mask1 = mask1(hw+1:end, hw+1:end);

mask2 = true(size(mask) + hw);
mask2(hw+1:end, 1:end-hw) = mask;
mask2 = imfill(mask2,'holes');
mask2 = mask2(hw+1:end, 1:end-hw);

mask3 = true(size(mask) + hw);
mask3(1:end-hw, 1:end-hw) = mask;
mask3 = imfill(mask3,'holes');
mask3 = mask3(1:end-hw, 1:end-hw);

mask4 = true(size(mask) + hw);
mask4(1:end-hw, hw+1:end) = mask;
mask4 = imfill(mask4,'holes');
mask4 = mask4(1:end-hw, hw+1:end);

mask = mask1 | mask2 | mask3 | mask4;

% Filter out large and small regions again
rp = regionprops(mask, 'PixelIdxList', 'Area');
rp = rp(vertcat(rp.Area) > 100 & vertcat(rp.Area) < 5000);

mask(:) = false;
mask(vertcat(rp.PixelIdxList)) = true;

subplot(3,2,5)
imshow(mask);

% Filter out regions with lots of positive concavity

% Get boundaries
[B,L] = bwboundaries(mask);

% Cycle over boundarys
for i = 1:length(B)
    b = B{i};

    % Filter boundary - use circular convolution
    b(:,1) = cconv(b(:,1),fspecial('gaussian',[1 7],1)',size(b,1));
    b(:,2) = cconv(b(:,2),fspecial('gaussian',[1 7],1)',size(b,1));

    % Find curvature
    curv_vec = zeros(size(b,1),1);
    for j = 1:size(b,1)
        p_b = b(mod(j-2,size(b,1))+1,:);  % p_b = point before
        p_m = b(mod(j,size(b,1))+1,:);    % p_m = point middle
        p_a = b(mod(j+2,size(b,1))+1,:);  % p_a = point after

        dx_ds = p_a(1)-p_m(1);              % First derivative
        dy_ds = p_a(2)-p_m(2);              % First derivative
        ddx_ds = p_a(1)-2*p_m(1)+p_b(1);    % Second derivative
        ddy_ds = p_a(2)-2*p_m(2)+p_b(2);    % Second derivative
        curv_vec(j+1) = dx_ds*ddy_ds-dy_ds*ddx_ds;
    end


    if (sum(curv_vec > 0)/length(curv_vec) > 0.4 || std(curv_vec) > 2.0)
        L(L == i) = 0;
    end
end

mask = L ~= 0;

subplot(3,2,6)
imshow(mask,[])

输出1：

enter image description here

输出2：

enter image description here

回复于 2024-05-05T18:57:49+08:00

如何有效地为这个显微图像创建BW蒙版？

1 回答

相关问题