Is it possible to speed up regionprops function?

Hello!

I'm writing code for counting particles passing under a camera that uses some filters and the function regionprops to get the area and the centroid (for tracking purposes). The issue is that I want to reduce the processing time. I'm using parallel computing to do it and it improved. Using the profiler tool I see that the line that take the most time to run is the one using regionprops. I would like to know if there is ay way to speed up this function (?) Is it possible to vectorize this funtion?

Thanks in advance for any help!

[DGM]

It depends how exactly you're using it, though full disclosure, I have zero familiarity with trying to approach this for parallelization.

If you're passing a label array (e.g. from bwlabel()) to regionprops(), that's slower than simply passing the binary image.

How to speed up regionprops() itself? Last I checked, many of regionprops() tasks are done in m-code (maybe not still in newer versions), and so one could certainly carve out a purpose-built subset of the code to cut out a little bit of the overhead. That said, the way to make any short m-code loops faster would probably be to rewrite it and compile it as mex.

[Clemente Gotelli]

Thank you @DGM for your answer! The code takes a several grayscale images, applies some filters to them (where the last function is imbinarize()) and after that applies regionprops() to get the data. So, I think the only option is to make a short m-code, as you suggest. The problem is that I'm not really expert with Matlab (first time I hear about mex compilation), so I don't think I'll be able to do it </3. Do you think Python could have better performance for this kind of functions?

Thanks again for your time!

[DGM]

I have no familiarity with python, and writing for mex is over my head.

I've browsed some parts of the tools regionprops() uses, but I haven't exactly wrapped my head around any of the routines yet.

I can think of naive ways to do the centroid and area calculations, but they all seem really inefficient. This might be one of the things that's being done in mex. It's been about a year since I looked. EDIT: at least in R2019b, it's all in m-code in the base file. It's a bunch of looping, so I don't know how much faster it can get in m-code. I doubt the overhead is a significant contributor to the time if that's the case.

Maybe someone else has a better idea.

[Clemente Gotelli]

Ahhh, ok, I understand. Seems difficult to speed up :( Anyway, thank you very much for your response! ;)

Answer 1

[DGM]

Consider the following demonstration of an alternative example.

For an 800x1000 logical image with roughly 450 objects, this is significantly faster than regionprops(). It's not really doing anything differently than regionprops; it's just stripping out all the overhead.

% 804x1004 logical image
A = ~imread('dotgrid.png');

% original method
a = timeit(@() testA(A))

a = 0.0110

% using simplified mcode
b = timeit(@() testB(A))

b = 0.0035

% time ratio
a/b

ans = 3.1097

function testA(A)
    S = regionprops(A,'area','centroid');
end

function testB(A)
    Acc = bwconncomp(A);
    nobj = Acc.NumObjects;
    imsz = Acc.ImageSize;
    objarea = zeros(nobj,1);
    objcentroid = zeros(nobj,2);
    for k = 1:nobj
        thispl = Acc.PixelIdxList{k};
        objarea(k) = numel(thispl);
        [r c] = ind2sub(imsz,thispl);
        objcentroid(k,:) = [mean(c) mean(r)];
    end
end

[Clemente Gotelli]

I implemented this short m-code and it is faster! Thank you very much for your help!

Answer 2

[Matt J]

What are the dimensions of the images? Perhaps you could downsample them.

Additionally, you might try reading the images into the odd slices of a 3D binary array,

BW=false(M,N,P);

BW(:,:,1)=Image1
BW(:,:,3)=Image2
BW(:,:,5)=Image3
...

Then you could process them all in a single call to regionprops(BW,___) or regionprops3(BW,___). I'm not sure how much of the difference that would make though.

[Clemente Gotelli]

Hello! Thanks for your answer. The dimensions of the images are 2700 x 1000.

What would be the advantage of reading the odd slices of a 3D array?

[Matt J]

Hello! Thanks for your answer. The dimensions of the images are 2700 x 1000.

I would try binning the image down to 675x250.

What would be the advantage of reading the odd slices of a 3D array?

To vectorize the call to regionprops. Instead of calling regionprops many times in a loop, call it once on the entire 3D stack. I'm not sure it will make much difference.

[Clemente Gotelli]

I tried this and it wasn't a big improvement. Sometimes it was actually slower... I realized some functions worked better in for loops compared with one call to process the entire 3D array. Crazy...

Answer 3

[Walter Roberson]

There is more than one way you can handle finding centroid and area.

In one of the approaches, you start by labeling the binary image. Having done that, you proceed object number by object number, comparing the array of labels to the current label, find() the coordinates, mean() along each dimension to find the centroid, count the coordinates to find the area.

In another approach, you scan the matrix looking for the next unexamined True location, and start walking through the adjacent locations to find connected locations, adding their coordinates to running totals and counting as you go, and marking each connected location as having already been examined, stopping when you no longer have adjacent locations to visit. You might use a breadth-first search strategy, perhaps. And you use the running total of coordinates to calculate the centroid and the count is the area.

Now, in theory, running a labelling pass can be done in parallel, breaking the array up into pieces... but you have to do fix-ups to account for particles that cross the boundaries. The fix-ups can be non-trivial. Imagine for example

   x  x  x  x
   ———————————
     x  x  x

Each of the x appears to be a separate object when considered only in its own section, but really it is all one object. And this kind of pattern could turn and follow vertical edges and then turn again and follow the next horizontal boundary of examination, so the isolated seemingly distinct particles at the top edge of the section could turn out to be the same object as ones on the bottom edge, joined by information in adjacent sections. So the fix-ups are not necessarily "nice" and can involve a bunch of work that is not easy to run in parallel.

There is another approach to labelling involving "equivalence classes" that involves scanning the entire matrix twice, once doing provisional labeling and then finalizing. (I think... guess I would have to implement it to be sure that I have not overlooked anything.) At the moment I do not see how it could be run in parallel sections, but perhaps there is a way that I am overlooking at the moment.

The alternative, non-labeling algorithm, can only be run in parallel if you can find a dividing line between sections. But I think my argument about the alternating blocks means that you can't reliably do that using only local knowledge.

Perhaps there are additional strategies that can be found by researching "flood filling".

So... best case might work nicely in parallel, but fixing up for the worst case could potentially add considerably to the cost of labeling?

Once you have the label matrix then the mass-examination algorithm I described above can be run in parallel.

There are perhaps more efficient ways to proceed once you have the label matrix. Hmmm, ah yes I think I just thought of one that could be run in parallel... mumble mumble... how much would the per-particle cost add up to...