matlab - How do I efficiently created a BW mask for this microscopic image? -
so background. tasked write matlab program count number yeast cells inside visible-light microscopic images. think first step cell segmentation. before got real experiment image set, developed algorithm use test image set utilizing watershed. looks this:
the first step of watershed generating bw mask cells. generate bwdist image imposed local minimums generated bw mask. can generate watershed easily.
as can see algorithm rely on successful generation of bw mask. because need generate bwdist image , markers it. originally, generate bw mask following following steps:
- generate local standard deviation of image sdimage = stdfilt(grayimage, ones(9))
- use bw thresholding generate initial bw mask binaryimage = sdimage < 8;
- use imclearborder clear background. use other code add cells on border back.
background finished. here problem
but today received new real data sets. image resolution smaller , light condition different test image set. color depth smaller. these make algorithm useless. here it:
using stdfilt failed generate clean images. instead generate stuff (note: have adjusted parameters stdfilt function , bw threshold value, following best result can get) :
as can see there light pixels in center of cells not necessary darker membrane. lead bw thresholding generate stuff this:
the new bw image after bw thresholding have either incomplete membrane or segmented cell centers , make them unsuitable other steps.
i start image processing , have no idea how should proceed. if have idea please me! thanks!
for convience, have attached link dropbox subset of images
i think there's fundamental problem in approach. algorithm uses stdfilt
in order binarize image. means you're assuming there there low "texture" in background and within cell. works first image. however, in second image there "texture" within cell, assumption broken.
i think stronger assumption there "ring" around each cell (valid both images posted). took approach of detecting ring instead.
so approach essentially:
- detect these rings (i use 'log' filter , binarize based on positive values. however, results in lot of "chatter"
- try remove of "chatter" filtering out small , large regions
- now, fill in these rings. however, there still "chatter" , filled regions between cells left
- again, remove small , large regions, since cells filled, increase bounds acceptable.
- there still bad regions, of bad areas going regions between cells. regions between cells detectable observing curvature around boundary of region. "bend inwards" lot, expressed mathematically large portion of boundary having negative curvature. also, remove rest of "chatter", these regions have large standard deviation in curvature of boundary, remove boundaries large standard deviation well.
overall, difficult part removing regions between cells , "chatter" without removing actual cells.
anyway, here's code (note there lot of heuristics , it's rough , based on code older projects, homeworks, , stackoverflow answers it's far finished):
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 , 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,[]); % 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 , 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 lots of positive concavity % boundaries [b,l] = bwboundaries(mask); % cycle on boundarys = 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); 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,[])
output1:
output2:
Comments
Post a Comment