Perimeter used by regionprops

Hi all,

I am using the output of bwconncomp as an input to regionprops to calculate some values for blobs in a binary image. One of the values that I calculate is "solidity" which is a direct output of regionprops. However, I am also interested in calculating circularity, which I define as 4*pi*Area/Perimeter^2. However, I keep getting values that are greater than one, which is unexpected. From running bwconncomp and regionprops on some artificial binary images with some pixels arbitrarily set to 1 to represent particles I can tell that regionprops does not calculate perimeter in the way that I would normally think of it. Does anybody know if the perimeter calculation fails for very small clusters of "on" pixels, or is there some rationale behind what regionprops is doing, and I just don't understand it?

I'll post the code I use to troubleshoot this:

%make an empty matrix
a = zeros(5,5)

%put a small cross shaped particle onto it
a(2, 2:4) = 1    
a(1:3, 3) = 1

%call bwconncomp
cc = bwconncomp(a)

%call regionprops
blobFacts = regionprops(cc, 'Perimeter', 'Solidity', 'Area')

%report values
p = blobFacts.Perimeter
s = blobFacts.Solidity
a = blobFacts.Area

%manually calculate circularity
c = (4*pi*a)/p^2

This gives a circularity of 1.9865, which is higher than 1 (what I intuitively think of as the max possible circularity). The area of the particle is 5, which is expected, but the perimeter is 5.6240, and I can't seem to figure out how regionprops is getting that value. Thanks in advance,

Clay

Answer 1

[nsunjaya]

The 'perimeter' calculation in "regionprops" is done based on a formula explained in the following paper:

A.M. Vossepoel, and A.W.M. Smeulders, "Vector code probability and metrication error in the representation of straight lines of finite length", Computer Graphics and Image Processing, Volume 20, Issue 4, December 1982, Pages 347–364.

You can also refer to the following blog post which gives a general overview of the algorithm used (please note that this is not an official documentation):

http://www.crisluengo.net/index.php/archives/310

You can see the perimeter is computed using a local function called "computePerimeterFromBoundary" in the "regionprops" function which corresponds to the formula discussed in the paper and blog post above:

perimeter = sum(isEven)*0.980 + sum(~isEven)*1.406 - sum(isCorner)*0.091;

Answer 2

[Matt J]

I think the example of a 2x2 square below is easier for discussion purposes. I assume you expect the calculated perimeter to be 8, since that is the theoretical perimeter of an ideal square of area 4. That would be the calculated result if regionprops traversed the boundary pixel edges. For most shapes, that will give a very poor approximation to the perimeter of the imaged object, so I think it's clear why regionprops doesn't do that.

Based on the documentation for regionprops, I would expect the calculated perimeter to be 4. In the doc, it sounds like the calculation traverses the boundary pixels center-to-center, accumulating the inter-pixel distance. Why it comes out to be 3.5560, I just don't know. However, even a computed perimeter of 4 will yield circularity=pi>1 in your formula. That is, unless you modify the area calculation so that only the area of the region bounded by the pixel centers (in this case 1) is computed. You could do that with polyarea().

Bottom line: the way regionprops computes area is not consistent with the way it computes perimeter, and so your formula for circularity will indeed break down for small pixel clusters.

    >> a

    a =

         0     0     0     0     0
         0     1     1     0     0
         0     1     1     0     0
         0     0     0     0     0
         0     0     0     0     0

    >> regionprops(a,'Area','Perimeter')

    ans = 

      struct with fields:

             Area: 4
        Perimeter: 3.5560

[Image Analyst]

" I assume you expect the calculated perimeter to be 8, since that is the theoretical perimeter of an ideal square of area 4" Not sure where you get 8. I would think 4. However if the square is at 45 degrees, then the distance between pixel centers is sqrt(2) so the perimeter would be 4*sqrt(2) instead of 4. Gotta admit, I'm not sure how it gets 3.556 instead of 4.

[Matt J]

Not sure where you get 8. I would think 4.

The perimeter of a length-2 square in continuous space is 8. Regionprops clearly considers this a length-2 square for the purpose of computing area (it works out to be 4), but not perimeter.

Gotta admit, I'm not sure how it gets 3.556 instead of 4.

I've put in a bug report.

[Image Analyst]

Can you explain, using your 2x2 example that you gave above (and replicated below), why the perimeter is 8, not 4?

[Matt J]

The formula for the perimeter of a square is 4*s where s is the length of the square. Because s= 2 pixels in this example, the perimeter is arguably 8 pixels.

In your figure illustration, you are not treating this as a 2x2 square, but rather a 1x1 square whose corners are defined by the pixel centers. As I told Clay, that's fine, and possibly what MATLAB intends, however, it is inconsistent with regionprops' area calculation. The area calculation looks at this as a 2x2 square, and gives a result of 4.

It looks like you had a similar conversation a few years ago in this thread.

[Image Analyst]

If you edit regionprops() you can see how it computes perimeter:

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function perimeter = computePerimeterFromBoundaryOld(B)

delta = diff(B).^2;
perimeter = sum(sqrt(sum(delta,2)));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function perimeter = computePerimeterFromBoundary(B)

delta = diff(B).^2;
if(size(delta,1) > 1)
    isCorner  = any(diff([delta;delta(1,:)]),2); % Count corners.
    isEven    = any(~delta,2);
    perimeter = sum(isEven)*0.980 + sum(~isEven)*1.406 - sum(isCorner)*0.091;
else
    perimeter = 0; % if the number of pixels is 1 or less.
end

[Matt J]

Bizarre!

[Image Analyst]

There is a similarly unusual computation for area in bwarea() where it's not strictly a pixel count, like in regionprops(), but a weighted pixel area depending on if that pixel is on an outside edge, outside corner, or inside corner.