perimeter from computeFeatures.shape is too low

[daboe01]

you can see this e.g. when you calculate cirularity from a perfect circle blob (2 x pi x (area/perimeter^2))

[aoles]

Thank you for reporting this. The observed discrepancy is an intrinsic artifact of the discrete nature of image data. In particular, the area and perimeter are calculated by counting pixels representing a given object, which is just an approximation of a continues shape. Let the following examples serve as an illustration.

library(EBImage)

## create empty background image
bg <- Image(0, c(99, 99))

## ideal circle with r=45
circle <- drawCircle(bg, x=50, y=50, r=45, col=1, fill=TRUE)
display(circle)

## shape characteristics
s <- computeFeatures.shape(circle)
s
#   s.area s.perimeter s.radius.mean s.radius.sd s.radius.min s.radius.max
# 1   6493         256      45.00955   0.2495668     44.55334     45.39824

In this case the area of the object is the total number of pixels constituting it,

table(circle)
# circle
# 0    1 
# 3308 6493 

while the perimeter is calculated by summing over the border pixels. This can be seen when we draw a ring which is similar to the circle above, but not filled. Then the perimeter is the same as for the filled circle. The ring's area is equal to its perimeter, as the object is constituted only from border pixels.

ring <- drawCircle(bg, x=50, y=50, r=45, col=1, fill=FALSE)
display(ring)

computeFeatures.shape(ring)
#   s.area s.perimeter s.radius.mean s.radius.sd s.radius.min s.radius.max
# 1    256         256      45.00955   0.2495668     44.55334     45.39824

As you noticed correctly, this approach underestimates the perimeter, which in turn causes the circularity of a perfect circle, defined as 4pi * (area/perimeter^2), to be higher than 1.

4 * pi * (s[1,"s.area"] / s[1,"s.perimeter"]^2)
# [1] 1.245017

I will look into possible improvements of calculating the perimeter.

[daboe01]

thank you for looking into this. maybe this resource is helpful: http://stackoverflow.com/questions/15912817/computing-the-circularity-of-a-binary-image best whishes, daniel

Am 12.08.2016 um 12:09 schrieb Andrzej Oleś notifications@github.com:

Thank you for reporting this. The observed discrepancy is an intrinsic artifact of the discrete nature of image data. In particular, the area and perimeter are calculated by counting pixels representing a given object, which is just an approximation of a continues shape. Let the following examples serve as an illustration.

library(EBImage)

create empty background image

bg <- Image(0, c(99, 99))

ideal circle with r=45

circle <- drawCircle(bg, x=50, y=50, r=45, col=1, fill=TRUE) display(circle)

shape characteristics

s <- computeFeatures.shape(circle) s

s.area s.perimeter s.radius.mean s.radius.sd s.radius.min s.radius.max

1 6493 256 45.00955 0.2495668 44.55334 45.39824

In this case the area of the object is the total number of pixels constituting it,

table(circle)

circle

0 1

3308 6493

while the perimeter is calculated by summing over the border pixels. This can be seen when we draw a ring which is similar to the circle above, but not filled. Then the perimeter is the same as for the filled circle. The ring's area is equal to its perimeter, as the object is constituted only from border pixels.

ring <- drawCircle(bg, x=50, y=50, r=45, col=1, fill=FALSE) display(ring)

computeFeatures.shape(ring)

s.area s.perimeter s.radius.mean s.radius.sd s.radius.min s.radius.max

1 256 256 45.00955 0.2495668 44.55334 45.39824

As you noticed correctly, this approach underestimates the perimeter, which in turn causes the circularity of a perfect circle, defined as 4pi * (area/perimeter^2), to be higher than 1.

4 * pi * (s[1,"s.area"] / s[1,"s.perimeter"]^2)

[1] 1.245017

I will look into possible improvements of calculating the perimeter.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or mute the thread.

[aoles]

Dear Daniel, Thank you for getting back to me, and for the link to the discussion on SO.

Indeed, in order to improve the perimeter estimation I will need to re-implement the algorithm so that it draws a chain of segments connecting the boundary pixels, and counts its length. It might take some time until I get back to this, though, as I'm busy with other stuff right now.

Best, Andrzej

[daboe01]

thank you, keep up your great work!

Am 12.08.2016 um 14:33 schrieb Andrzej Oleś notifications@github.com:

Dear Daniel, Thank you for getting back to me, and for the link to the discussion on SO.

Indeed, in order to improve the perimeter estimation I will need to re-implement the algorithm so that it draws a chain of segments connecting the boundary pixels, and counts its length. It might take some time until I get back to this, though, as I'm busy with other stuff right now.

Best, Andrzej

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or mute the thread.

[dcornet]

Any hope for a solution on this issue? I'm using EBImage (thanks !) for blob shape characterization but the varying underestimation of the perimeter introduce some bias. Best,

denis

[dcornet]

For any purpose, here is a hacky estimate of the perimeter using above mentionned links and the trajr library. Actually, the less the resolution the higher the error. But I guess the smoothing function could be better use to manage this. For objects with holes inside it may be necessary to use floodFill before.

library(EBImage) library(trajr) library(ggplot2)

Usefull functions from https://stackoverflow.com/questions/15912817

s1 <- function(z) cbind(rep(0,nrow(z)), z[,-ncol(z)] ) # left value s2 <- function(z) cbind(z[,-1], rep(0,nrow(z)) ) # right value s3 <- function(z) rbind(rep(0,ncol(z)), z[-nrow(z),] ) # upper value s4 <- function(z) rbind(z[-1,], rep(0,ncol(z)) ) # lower value edge <- function(z) z & !(s1(z)&s2(z)&s3(z)&s4(z)) area <- function(z) sum(z) perimeter <- function(z) sum( z != 0 & s1(z) s2(z) s3(z) s4(z) == 0) circularity <- function(z) 4 pi * area(z) / perimeter(z)^2

Create empty background image

bg <- Image(0, c(400, 400))

Ideal circle

circle <- drawCircle(bg, x=150, y=150, r=80, col=1, fill=TRUE) display(circle) area(circle) perimeter(circle) circularity(circle)

Get edge pixels coordinates

e <- edge(circle) Dim<-dim(e) e<-as.data.frame(array(e, dim=prod(Dim))) e$x<-rep(1:Dim[1], each=Dim[2]) e$y<-rep(Dim[2]:1, Dim[1]) colnames(e)<-c("times", "x", "y") e<-subset(e, times==TRUE)

Get the pixels order to produce a path

e$times<-0 for (i in 1 : (nrow(e)-1)) { if (i==1) { e$times[i]<-1 ei<-subset(e, x %in% c(e$x[i], e$x[i]+1, e$x[i]-1) & y %in% c(e$y[i], e$y[i]+1, e$y[i]-1) & times==0) ei<-ei[base::which.min(abs(ei$x-e[i,"x"])+abs(ei$y-e[i,"y"])), ] e[which(e$x==ei$x & e$y==ei$y), "times"]<-2 } else { epos<-e[which(e$times==i), ] ei<-subset(e, x %in% c(epos$x, epos$x+1, epos$x-1) & y %in% c(epos$y, epos$y+1, epos$y-1) & times==0) ei<-ei[base::which.min(abs(ei$x-epos$x)+abs(ei$y-epos$y)), ] e[which(e$x==ei$x & e$y==ei$y), "times"]<-i+1 } e <- e[order(e$times),] e<-e[,c(2,3,1)] e[nrow(e)+1, ]<-c(e[1,1], e[1,2], nrow(e)+1) ggplot(e, aes(x, y))+geom_path()+coord_fixed()

Define a traject (path), smooth it and measure its length

trj <- TrajFromCoords(e, "x", "y", "times", spatialUnits="pixels", timeUnits="s") plot(trj) nn<-(floor(nrow(trj)/64) + floor(nrow(trj)/64) %% 2)+1 trjs<-TrajSmoothSG(trj, p = 4, n=nn, 0) plot(trjs) perim<-TrajLength(trjs)

Check circularity

4 pi area(circle) / perim^2