How to identify the nth largest components in binary image stacks

Hello,

I have some code that is programmed to identify the largest blob within a stack of binary images and delete all surrounding blobs, exporting the new image files.

I want to edit to code to identify nth largest blobs and delete everything else. I will probably go as high as 3rd to 5th largest blobs.

Currently the code looks like this:

for i=0:num_of_files
CC(:,i+1) = bwconncomp(BW(:,:,i+1),26);
end

n1 =zeros(1023,997,11);
for i = 0:num_of_files   
newimage=zeros(1023,997);
[~, maxcell(i+1)] = max(cellfun(@numel, CC(i+1).PixelIdxList));
newimage(CC(i+1).PixelIdxList{1,maxcell(i+1)}) = 255;
n1 (:,:,i+1)= newimage;
imwrite(newimage,sprintf('%s%d.%s', name_of_output, i, file_type));

How can I edit to take highest 3 blobs into account in each image slice?

Regards, R.

[Guillaume]

Can you explain why you have:

CC(:, i+1) = ...

Instead of

CC(i+1) = ...

?

Can you also explain why you're using a connectivity of 26 when you're passing 2D slices to bwconncomp?

Also, it looks like the number of files/slices is num_of_files+1, so the num_of_files variable name is misleading.

Why isn't n1 declared as

n1 = zeros(1023, 997, num_of_files+1);

or even better

n1 = zeros(size(BW));

?

Note that if you're putting 255 in the image, then n1 should be declared of type uint8:

n1 = zeros(size(BW), 'uint8');

As for your question, can't you use maxk instead of max, assuming you're on R2017b or later?

[RachelCahalane]

Thanks for your comments.

I am currently on matlab version R2015b. Is there another function to use for older matlab versions?

Answer 1

[Guillaume]

First, I would encourage you to address the questions I've raised as they're all logic gaps in your code.

Rather than looping twice over CC let's do it all in one go:

%inputs:
%   num_of_files:    variable with a misleading name that is the number of files minus 1, scalar non-negative integer.
%   BW:              binary image stack, a 3D array.
%   name_of_output:  base name of file to be written for each layer, char vector.
%   file_type:       image type to be written by imwrite, char vector.
%   blobs_to_keep:   how many of the largest blobs to keep in each layer, scalar positive integer.
blobs_to_keep = 3;
 
%outputs:
%   outimage:        uint8 image stack where only the first blobs_to_keep largest blobs of each layer in BW is kept.
   
outimage = zeros(size(BW), 'uint8');
for layer = 1:num_of_files + 1
   CC = bwconncomp(BW(:, :, layer), 8);  %more than 8 is pointless on 2D images
   blobsizes = cellfun(@numel, CC.pixelIdxList);
   [~, order] = sort(blobsizes, 'descend');
   pixelstokeep = vertcat(CC.PixelIdxList{order(1:blobstokeep)});  
   outimage(pixelstokeep + prod(CC.ImageSize) * (layer - 1)) = 255;  
   %adding prod(CC.ImageSize) * (layer - 1) ensures the pixels in the correct layer are sets
   %this avoids using an intermediate 2D image to set the pixels
   imwrite(outimage(:, :, layer), sprintf('%s%d.%s', name_of_output, layer-1, file_type));
end

[RachelCahalane]

Thanks @Guillaume It works perfectly except I do get an error at the end as follows:

>> Isolate_Nth_Largest_Components Index exceeds matrix dimensions.

Error in Isolate_Nth_Largest_Components (line 25) pixelstokeep = vertcat(CC.PixelIdxList{order(1:blobs_to_keep)});

So for an original image cohort of 1992 I am only getting 1076 image results?

[Guillaume]

Index exceeds matrix dimensions

That would be because that particular layer has less than blobs_to_keep blobs.

   pixelstokeep = vertcat(CC.PixelIdxList{order(1:min(blobstokeep, numel(CC.PixelIdxList)))});

would get rid of the error. I'd also add

if blobstokeep > numel(CC.pixelIdxList)
    warning('layer %d has only %n blobs. Keeping all the blobs.', layer, numel(CC.pixelidxList));
end

so you know that it happened.