The threshold: Otsu and the histogram.

The idea

In a photo of a tetrazolium plate there are two kinds of pixel: seed and blue paper. If each pixel becomes a number, the histogram of those numbers usually has two humps, one for the paper and one for the seeds. Choosing the threshold means choosing where to cut the valley between them.

In 1979, Nobuyuki Otsu proposed a criterion with no parameters at all [1]. For each possible cut, the histogram is split into two classes. The best cut is the one that makes the two classes differ most in their means, weighted by how many pixels each has. The animation Live segmentation runs this calculation on one seed.

Before the histogram, a channel

The photo has three numbers per pixel (red, green and blue) and the threshold needs one. The choice of that number decides much of the result. Three candidates make sense for a light seed on blue paper.

Gray is the video luma, \(Y'=0{,}299R+0{,}587G+0{,}114B\): the seed is lighter than the paper. The \((R+G)/2-B\) channel measures how much more red and green than blue a pixel has: high on the cream seed coat and the red embryo, strongly negative on the paper. The third is b*, the yellow-to-blue axis of the CIELAB color space, standardized in 1976 [5]:

where \(Y\) and \(Z\) come from the linear R, G and B values through a fixed matrix, and \(Y_n\), \(Z_n\) are those of the reference white. Blue paper has a strongly negative b*. Cream or red seed sits well above it.

photoOtsu on grayOtsu on (R+G)/2 − BOtsu on b*

What defeats gray is uneven light. In this crop, the paper in the shadow strip has a gray value of 43.0 and the paper below has 59.6. The typical seed has 87.6. The light moves the paper by 16.6 levels out of a margin of 28 to the seed, more than half of it, and the Otsu threshold (66) ends up a few levels from the light paper. On b*, the same shadow moves the paper from −50.6 to −38.5, or 12.1 units, out of a margin of 36.7 to the seed (−13.9): one third.

The reason lies in formula (1) itself. If the light on the paper is multiplied by \(k\), the linear values X, Y and Z are too, and \(b^*\) ends up multiplied by \(k^{1/3}\). The paper stays on the blue side, closer to zero but still far from the seed. Shadow mainly changes lightness, and b* measures color.

This does not make b* a champion. In the Where it fails section, the \((R+G)/2-B\) channel holds up better on an image with few seeds, because on this plate the paper is smoother in it. The right channel is the one that separates the seed from the paper in units of the paper's noise and changes little when the light changes.

Otsu's calculation

Once the channel is chosen, each pixel falls into one of \(L\) histogram bins (here \(L=256\)). Let \(p_i\) be the fraction of pixels in bin \(i\). A cut at \(k\) separates the bins \(0,\dots,k\) (background) from the bins \(k+1,\dots,L-1\) (seed).

The Otsu threshold is \(k^*=\arg\max_k\sigma_B^2(k)\). The two pieces of the histogram end up as tight as possible around their means, and the means as far apart as possible. Otsu himself proposed measuring how good the separation is by \(\eta=\sigma_B^2(k^*)/\sigma_T^2\), a number between 0 and 1 that does not change if the channel is multiplied by a constant or shifted by another [1].

Two values give \(\eta\) its scale. A single Gaussian, cut at the mean, gives \(\eta=2/\pi\approx 0{,}637\): that is what the method finds when there are no two classes. Two narrow, distant classes take \(\eta\) close to 1. In the crops on this page, with seeds, \(\eta\) falls between 0.65 and 0.78.

A single pass

Recomputing means and variances for each of the 256 cuts would be wasteful. With the cumulative sums \(\omega(k)=\sum_{i\le k}p_i\) and \(\mu(k)=\sum_{i\le k}i\,p_i\), formula (2) becomes

because \(\mu_0=\mu/\omega\), \(\mu_1=(\mu_T-\mu)/(1-\omega)\) and \(\mu_0-\mu_1=(\mu-\omega\mu_T)/\big(\omega(1-\omega)\big)\). Each cut costs a single calculation. The whole method is one pass over the pixels to build the histogram and one pass over the bins. This is what runs in the instruments below:

// h: histogram with 256 bins; returns the cut k*
let n = sum(h), mT = Σ i·h[i]/n, w = 0, m = 0, best = -1, kStar = 0;
for (let i = 0; i < 256; i++) {
  w += h[i]/n;  m += i*h[i]/n;              // accumulated ω(k) and μ(k)
  const sB = (mT*w - m)**2 / (w*(1 - w));   // equation (3)
  if (sB > best) { best = sB; kStar = i; }
}

The threshold sits midway between the means

On the instrument's clean plate, in b*, the two classes of the Otsu cut have means −41.52 and −12.42. The midpoint is −26.97, and the Otsu threshold is −27.0.

This is the same fixed point as the iterative method of Ridler and Calvard [4]: start from any threshold, compute the two means, put the threshold midway between them and repeat. Among the fixed points, Otsu takes the one with the largest \(\sigma_B^2\). The result also explains the failures. The threshold only sees the means. The width of each class and its size matter only through how they shift the means.

On the real plate

The instrument below loads five real crops, computes the three channels in your browser and shows the mask, the histogram and the \(\sigma_B^2(t)\) curve. The Otsu button jumps to the maximum of the curve. The Minimum error button uses the criterion from the Alternatives section.

On the full plate, the three channels agree: Otsu marks 26.8% of the image on b*, 28.8% on \((R+G)/2-B\) and 28.8% on gray. On the shadow strip, \(\eta\) is practically the same in all three channels (0.654, 0.659 and 0.654) and the masks are very different, as in Figure 1. \(\eta\) measures the histogram, and the histogram does not know where each pixel is.

On the plate glare, all three channels get it wrong (30.0% to 37.1% of the image), because the glare strip is lighter than the paper and, in part, more yellow. On the paper-only crop, Otsu still returns a threshold and marks between 41.5% and 49.6% of the image as seed, with \(\eta\) between 0.60 and 0.66, close to the \(2/\pi\) of a single Gaussian. The method always answers, and the answer does not say whether there was anything to find.

Where it fails

Few seeds on a large background

To isolate the effect of class size, Figure 2 builds histograms from real pixels of one plate: 8,817 pixels from inside the annotated seeds and 223,729 pixels from a quadrant with no seed. The mixture is exact, the weighted sum of the two histograms, for each seed fraction. Only the proportion changes.

The jump has a calculation behind it. Take a small seed fraction \(f\), the difference \(\Delta\) between the seed and paper means, and the standard deviation \(\sigma\) of the paper. The right cut gives \(\sigma_B^2\approx f\,\Delta^2\) by equation (2). Cutting the paper in half, if it is close to Gaussian, gives \(\sigma_B^2\approx\frac{2}{\pi}\sigma^2\). Otsu switches cuts when the second passes the first:

The estimate ignores the width of the seed class and gets the breakdown point right within 13%. What matters is the ratio \(\Delta/\sigma\), the distance between seed and paper in units of the paper's noise. On this plate each seed takes up about 0.4% of the 946 × 946 px crop: on gray, a crop with a single seed is already below the breakdown point.

Uneven light and glare

A global threshold assumes the paper has the same value across the whole image. The shadow and glare in Figure 1 and in the instrument break that assumption, and no method that looks only at the histogram fixes it, because the histogram has lost the position of the pixels.

A histogram without two humps

With a single class, as in the paper crop, Otsu cuts the noise in half. With weak contrast, the two humps merge and the valley disappears. In both cases the method returns a number with the same confidence as always.

Alternatives

Minimum error

Kittler and Illingworth modeled each class as a Gaussian with its own mean, variance and weight and chose the cut that best explains the histogram under that model [2]. The criterion is

and the threshold is the \(k\) that minimizes \(J\). Because the width of each class enters the calculation, the small class does not lose to the large one. On the mixture in Figure 2, minimum error does not jump: on b*, with 1% seed, it marks 1.02%. On the real crops, though, it fails in another way. On the shadow strip, on gray, it marks 0.4% of the image and misses almost all the seeds. On the paper crop, on b*, it marks 89%. When the histogram does not have the shape of two Gaussians, the model has nothing to hold on to.

Gaussian mixture

The smooth version of the same model fits the two Gaussians with the EM algorithm and cuts where the two weighted densities cross. The same problem shows up with another histogram in SeedCounter: that of a shape descriptor among the seeds in a scene, where isolated seeds and clusters form two humps. A two-component mixture would find that cut without a hand-picked constant. This is a proposal from the project, not yet tested.

Local threshold

Instead of one \(t\) for the image, one \(t(x)\) for each neighborhood. Sauvola and Pietikäinen use the mean \(m(x)\) and the standard deviation \(s(x)\) in a window: \(t(x)=m(x)\,[1+\kappa\,(s(x)/R-1)]\) [6]. Another way out is to flatten the light first: fit a smooth surface to the paper and subtract it, and only then apply a global threshold. Sezgin and Sankur cataloged and compared 40 thresholding methods, grouped by the information they use: histogram shape, clustering, entropy, object attributes, spatial correlation and local surface [3].

In SeedCounter

References

1. Otsu N. (1979). A threshold selection method from gray-level histograms. IEEE Transactions on Systems, Man, and Cybernetics 9(1), 62–66. doi:10.1109/TSMC.1979.4310076
2. Kittler J., Illingworth J. (1986). Minimum error thresholding. Pattern Recognition 19(1), 41–47. doi:10.1016/0031-3203(86)90030-0
3. Sezgin M., Sankur B. (2004). Survey over image thresholding techniques and quantitative performance evaluation. Journal of Electronic Imaging 13(1), 146–165. doi:10.1117/1.1631315
4. Ridler T. W., Calvard S. (1978). Picture thresholding using an iterative selection method. IEEE Transactions on Systems, Man, and Cybernetics 8(8), 630–632. doi:10.1109/TSMC.1978.4310039
5. Robertson A. R. (1977). The CIE 1976 color-difference formulae. Color Research & Application 2(1), 7–11. doi:10.1002/j.1520-6378.1977.tb00104.x
6. Sauvola J., Pietikäinen M. (2000). Adaptive document image binarization. Pattern Recognition 33(2), 225–236. doi:10.1016/S0031-3203(99)00055-2