Conclusion

steps

  • blur the image
  • apply a threshold to turn cell into ‘blobs’
  • count the ‘blobs’ ## comments
  • watershed function can not resolve two linked cells.
im <- readImage("papSmear_bg.jpg")
par(mfrow=c(1,2))
plot(im)
EBImage::display(im, method = "raster") # raster method means within R

make the image brighter by mutiplying all the values by 2

# convert image de gray scale
# im <- channel(im.c, "gray")
im2 <- im * 2
par(mfrow=c(1,2))
plot(im)
plot(im2)

# gaussian blur
im.blur <- gblur(im, sigma= 0.01)
im2.blur <- gblur(im2, sigma = 1)
par(mfrow=c(1,2))
plot(im.blur)
plot(im2.blur)

threshold using Otsu’s method used to reduce of graylevel image to a binary image.

# only gray images are accepted
# gives a threshold value using Otsu algorithm
otsu(im.blur)
## [1] 0.2988281
#otsu(im2.blur)

#im.c2.g.blur.otsu <- im.c2.g.blur > otsu(im.c2.g.blur)
nuc = readImage(system.file("images", "nuclei.tif", package="EBImage"))
nuc_gblur = gblur(nuc, sigma = 5)
#display(nuc_gblur, method = "raster", all = TRUE)
threshold = otsu(nuc)
threshold_gblur <- otsu(nuc_gblur) 
threshold
## [1] 0.3535156 0.4082031 0.3808594 0.4121094
threshold_gblur
## [1] 0.2519531 0.3144531 0.2871094 0.3183594
#display(nuc_gblur > threshold, all =FALSE)
par(mfrow=c(1,2))
display(nuc > threshold, method = "raster", all = TRUE)
display(nuc_gblur > threshold, method = "raster", all = TRUE)

max(bwlabel(nuc > threshold))
## [1] 181
max(bwlabel(nuc_gblur > threshold_gblur))
## [1] 1066
display(colorLabels(bwlabel(nuc > threshold)), method = "raster", all = TRUE)
display(colorLabels(bwlabel(nuc_gblur > threshold_gblur)), method = "raster", all = TRUE)

im.blur.thr.cnt <- bwlabel(im.blur > otsu(im.blur))
N <- max(im.blur.thr.cnt)
N
## [1] 195
display(colorLabels(im.blur.thr.cnt), method= "raster")

R function to get label of cells

getCelLab <- function(img, bright, sigma ){

im <- readImage(img)
im <- EBImage::channel(im, mode= "gray")
im2 <- im * bright
im.blur <- gblur(im, sigma)
plot(im.blur)
## filtering option 

im.blur.thr.cnt <- bwlabel(im.blur > otsu(im.blur))
#im.blur.the.cnt <- bwlabel(watershed( distmap(im.blur), 2 ))
#im.blur.thr.cnt <- bwlabel(adaptative_thre(31, "disc"))
#im.blur.thr.cnt <- propagate(seeds = im.blur, x = im.blur, lambda = 100)


N <- max(im.blur.thr.cnt)
display(colorLabels(im.blur.thr.cnt), method= "raster") 
return(N)
}
#par(mfrow=c(2,2))
getCelLab("papSmear.jpg", 2, 5)

## [1] 16
getCelLab("papSmear_g.jpg", 2, 2)

## [1] 22
getCelLab("papSmear_bg.jpg", 2, 0.1)

## [1] 195
getCelLab("basale_inverse_gray.jpg", 2, 1)

## [1] 43
adaptative_thre <- function(size, shape){
disc = makeBrush(size, shape)
#disc <- fillHull(disc)
disc = disc / sum(disc)
offset = 0.05
im.blur_bg = filter2( im.blur, disc )
im.blur_th = im.blur > im.blur_bg + offset
return(im.blur_th)
}
plot(adaptative_thre(31, "disc"))

plot(adaptative_thre(9, "disc"))

getCelLab with adaptative_thre

getCelLab <- function(img, bright, sigma ){

im <- readImage("papSmear_bg.jpg")
im2 <- im * bright
im.blur <- gblur(im, sigma)
#im.blur.thr.cnt <- bwlabel(im.blur > otsu(im.blur))
im.blur.thr.cnt <- bwlabel(adaptative_thre(31, "disc"))
N <- max(im.blur.thr.cnt)
display(colorLabels(im.blur.thr.cnt), method= "raster") 
return(N)
}
#par(mfrow=c(2,2))
getCelLab("papSmear.jpg", 2, 5)

## [1] 780
getCelLab("papSmear_g.jpg", 2, 2)

## [1] 780
getCelLab("papSmear_bg.jpg", 2, 0.1)

## [1] 780

get graded magnetude image

im <- load.image("papSmear_g.jpg")
dx <- imager::imgradient(im,"x")
dy <- imager::imgradient(im,"y")
grad.mag <- sqrt(dx^2+dy^2)
plot(grad.mag)