Counting cells and extensions

The image below has 3 channels - blue (DAPI nuclei), red (Ph3 cell death), green (Ki67 cell proliferation)

In [14]:
library(SimpleITK)
source("setup.R")
cntrl <- ReadImage("images/Control.tif")
# it is a 3D tif, so need to turn it into colour to display
red <- cntrl[,,1]
green <- cntrl[,,2]
blue <- cntrl[,,3]
# first use of a filter
cntrl.colour <- ComposeImageFilter()$Execute(red, green, blue)
cntrl.colour

Standard microscopy analysis approach

  1. Segment the dapi using a threshold.
  2. Split touching cells using a watershed.
  3. Count them.
  4. Threshold the dying cells (don't need to split).
  5. Count them.
  6. Threshold the proliferating cells.
  7. Count them.
  8. Fine tuning - use DAPI channel to refine segmentation of the other channels.
In [15]:
segChannel
Out[15]:
function (dapi, dtsmooth = 1, osmooth = 0.5) 
{
    dapi.smooth <- SmoothingRecursiveGaussian(dapi, osmooth)
    th <- LiThresholdImageFilter()
    th$SetOutsideValue(1)
    th$SetInsideValue(0)
    B <- th$Execute(dapi.smooth)
    g <- splitBlobs(B, dtsmooth)
    return(list(thresh = B, labels = g$labels, dist = g$dist, 
        peaks = g$peaks))
}
In [16]:
splitBlobs
Out[16]:
function (mask, smooth = 1) 
{
    DT <- DanielssonDistanceMapImageFilter()
    DT <- DT$UseImageSpacingOn()
    distim <- DT$Execute(!mask)
    distimS <- SmoothingRecursiveGaussian(distim, smooth, TRUE)
    distim <- distimS * Cast(distim > 0, "sitkFloat32")
    peakF <- RegionalMaximaImageFilter()$SetForegroundValue(1)
    peakF <- peakF$FullyConnectedOn()
    peaks <- peakF$Execute(distim)
    markers <- ConnectedComponent(peaks, TRUE)
    WS <- MorphologicalWatershedFromMarkers(-1 * distim, markers, 
        TRUE, TRUE)
    WS <- WS * Cast(distim > 0, WS$GetPixelID())
    return(list(labels = WS, dist = distimS, peaks = peaks))
}

Apply thresholding and splitting to each channel and count the cells

In [20]:
dapi.cells <- segChannel(blue, 1)
ph3.cells <- segChannel(red, 1)
Ki67.cells <- segChannel(green, 1)

## stats
dapi.stats <- LabelShapeStatisticsImageFilter()
lmap <- dapi.stats$Execute(dapi.cells$labels)
dapi.cells.count <- dapi.stats$GetNumberOfLabels()

ph3.stats <- LabelShapeStatisticsImageFilter()
lmap <- ph3.stats$Execute(ph3.cells$labels)
ph3.cells.count <- ph3.stats$GetNumberOfLabels()

Ki67.stats <- LabelShapeStatisticsImageFilter()
lmap <- Ki67.stats$Execute(Ki67.cells$labels)
Ki67.cells.count <- Ki67.stats$GetNumberOfLabels()

c(dapi.cells.count, ph3.cells.count, Ki67.cells.count)
Out[20]:
  1. 439
  2. 101
  3. 2
In [19]:
dapi.cells$dist
In [18]:
dapi.cells$labels > 0
In [21]:
ph3.cells$labels > 0