大津の二値化法
大津の二値化法(おおつのにちかほう、Otsu's method、大津の方法などとも呼ばれる)とは、コンピュータビジョンやにおいて、自動画像しきい値処理を行うために使用される手法。その名は大津展之にちなむ。最も単純な形式では、アルゴリズムはピクセルをforegroundとbackgroundの2つのクラスに分ける1つの強度しきい値を返す。このしきい値はクラス内の強度分散を最小化することにより、または同等にクラス間の分散を最大化することにより決定される。
大津の二値化はフィッシャーの判別分析の1次元の離散に相当するものであり、ジェンクス最適化法に関連しており、強度ヒストグラムで行われる大域最適なK平均法と同等である。マルチレベルのしきい値処理へ拡大することは最初の論文で説明されており、以降計算効率の高い実装が提案されている。
アルゴリズムは、2つのクラスの分散の加重和として定義されるクラス内分散を最小化するしきい値を探す。
\sigma^2_w(t)=\omega_0(t)\sigma^2_0(t)+\omega_1(t)\sigma^2_1(t)
加重\omega_{0}と\omega_{1}はしきい値tにより分けられた2つのクラスの確率であり、\sigma^2_{0}と\sigma^2_{1}はこれら2つのクラスの分散である。
クラスの確率\omega_{0,1}(t)はヒストグラムのL個のビンから計算される。
\begin{align}
\omega_0(t) & =\sum_{i=0}^{t-1} p(i)\\[4pt]
\omega_1(t) & =\sum_{i=t}^{L-1} p(i)
\end{align}
2つのクラスでは、クラス内分散を最小化することはクラス間分散を最大化することと同じである。
\begin{align}
\sigma^2_b(t) & =\sigma^2-\sigma^2_w(t)=\omega_0(t)(\mu_0-\mu_T)^2+\omega_1(t)(\mu_1-\mu_T)^2 \\
& =\omega_0(t) \omega_1(t) \left[\mu_0(t)-\mu_1(t)\right]^2
\end{align}
これはクラス確率\omegaとクラス平均\muで表され、クラス平均\mu_0(t)、\mu_1(t)および\mu_Tは、
\begin{align}
\mu_0(t) & = \frac{\sum_{i=0}^{t-1} i p(i)}{\omega_0(t)} \\[4pt]
\mu_1(t) & = \frac{\sum_{i=t}^{L-1} i p(i)}{\omega_1(t)} \\
\mu_T & = \sum_{i=0}^{L-1} ip(i)
\end{align}
である。以上の関係は以下のように簡単に確かめることができる。
\begin{align}
\omega_0\mu_0+\omega_1\mu_1 & = \mu_T \\
\omega_0+\omega_1 & =1
\end{align}
クラス確率とクラス平均は、繰り返し計算することができる。このアイデアは効果的なアルゴリズムを生む。
各強度レベルのヒストグラムと確率を計算する。
\omega_i(0)と\mu_i(0)の初期値を設定する。
t = 1, \ldots  最大強度までの全ての考えうるしきい値で次のステップを行う。
\omega_iと\mu_iを更新する。
\sigma^2_b(t)を計算する。
望ましいしきい値は\sigma^2_b(t)が最大となる値である。
histogramCountsはさまざまなグレーレベルのグレースケール画像(通常は8ビット画像)の256画素のヒストグラムである。levelは画像のしきい値 (double) である。
MatlabにはImage Processing Toolboxに組み込み関数graythresh()とmultithresh()があり、それぞれ大津の手法とマルチ大津の手法(Multi Otsu's method)で実装されている。
この実装にはNumPyライブラリが必要であるOpenCVやScikit-imageなどの画像処理専用のPythonライブラリは、アルゴリズムの組み込み実装を提案する。
大津の方法は、ヒストグラムが2つのピークの間に深く鋭い谷を持つ二峰性分布を有するときにうまく機能する。
他の全ての大域的なしきい値処理と同様に、ノイズが大きく、対象のサイズが小さく、明暗が不均一で、クラス内の分散がクラス間の分散よりも大きい場合にパフォーマンスが低下する。このような場合、大津の方法を局所的に適用することが開発されている。
さらに、大津の方法の数学的根拠は、画像のヒストグラムを分散が等しくサイズが等しい2つの正規分布の混合としてモデル化する。しかし、大津のしきい値処理はこれらの仮定を満たさない場合でも満足のいく結果をもたらす可能性がある。同様に統計検定(大津の方法が密接に関連する)は、動作の仮定が完全に満たされていない場合でも正しく実行される。しかし、これらの仮定から生じる深刻な逸脱を無くすために、Kittler-Illingworthの方法などの大津の方法のいくつかのバリエーションが提案されている。
大津の方法の限界に対処するためにさまざまな拡張が開発されてきた。人気のある拡張の1つは、ノイズの多い画像のオブジェクトセグメンテーションのタスクに適した2次元大津の方法である。ここでは、特定のピクセルの強度値をそのすぐ近傍の平均強度と比較することでセグメンテーション結果を改良する。
各ピクセルで、近傍の平均グレーレベル値が計算される。与えられたピクセルのグレーレベルをL個の離散値に分割し、平均グレーレベルも同じL個の値に分割する。その後、ピクセルのグレーレベルと近傍の平均のペア(i, j)が形成される。各ペアはL\times L個の考えうる2次元ビンの1つに属する。ペア(i, j)の出現の総数(頻度)f_{ij}を画像のピクセルの総数Nで割ったものは、2次元ヒストグラムで同時確率密度関数を定義する。
P_{ij} = \frac{f_{ij} N, \qquad \sum_{i=0}^{L-1}\sum_{j=0}^{L-1} P_{ij}=1
そして、2次元大津の方法は、2次元ヒストグラムに基づいて次のように開発されている。
2つのクラスの確率は、次のように表せる。
\begin{align}
\omega_0 & = \sum_{i=0}^{s-1} \sum_{j=0}^{t-1} P_{ij} \\
\omega_1 & = \sum_{i=s}^{L-1} \sum_{j=t}^{L-1} P_{ij}
\end{align}
2つのクラスの強度平均ベクトルと合計平均ベクトルは次のように表せる。
\begin{align}
\mu_0 & =[\mu_{0i}, \mu_{0j}]^T = \left[\sum_{i=0}^{s-1} \sum_{j=0}^{t-1} i \frac{P_{ij}{\omega_0}, \sum_{i=0}^{s-1}\sum_{j=0}^{t-1} j \frac{P_{ij}{ \omega_0} \right]^T \\
\mu_1 & =[\mu_{1i}, \mu_{1j}]^T = \left[\sum_{i=s}^{L-1}\sum_{j=t}^{L-1} i \frac{P_{ij}{\omega_1}, \sum_{i=s}^{L-1}\sum_{j=t}^{L-1} j \frac{P_{ij}{\omega_1} \right]^T \\
\mu_T & =[\mu_{Ti}, \mu_{Tj}]^T = \left[\sum_{i=0}^{L-1} \sum_{j=0}^{L-1} i P_{ij}, \sum_{i=0}^{L-1}\sum_{j=0}^{L-1} j P_{ij}\right]^T
\end{align}
ほとんどの場合、非対角の確率はわずかであるため、次のことを簡単に確認することができる。
\omega_0+\omega_1 \cong 1
\omega_0\mu_0+\omega_1\mu_1 \cong \mu_T
クラス間離散行列は次のように定義される。
S_b = \sum_{k=0}^1\omega_k[(\mu_k-\mu_T)(\mu_k-\mu_T)^T]
離散行列のトレースは次のように表される。
\begin{align}
& \operatorname{tr}(S_b) \\[4pt]
= {} & \omega_0[(\mu_{0i}-\mu_{Ti})^2 + (\mu_{0j}-\mu_{Tj})^2] + \omega_1[(\mu_{1i}-\mu_{Ti})^2 + (\mu_{1j}-\mu_{Tj})^2] \\[4pt]
= {} & \frac{(\mu_{Ti}\omega_0-\mu_i)^2 + (\mu_{Tj}\omega_0-\mu_j)^2}{\omega_0(1-\omega_0)}
\end{align}
\mu_i = \sum_{i=0}^{s-1}\sum_{j=0}^{t-1}iP_{ij}
\mu_j = \sum_{i=0}^{s-1}\sum_{j=0}^{t-1}jP_{ij}
1次元の大津の方法と同様に、最適なしきい値(s,t)は、\operatorname{tr}(S_b)を最大化することで取得される。
sとtは、1次元の大津の方法と同様に繰り返し取得される。sとtの値は、\operatorname{tr}(S_b)の最大値を取得するまで変更される。つまり、
\operatorname{tr}(S_b)を評価するために、高速再起動的プログラミングアルゴリズムを使用して時間パフォーマンスを向上させることができることに留意してください。しかし、動的プログラミングのアプローチを使用しても、2次元大津の方法は依然として時間計算量が非常に複雑である。したがって、計算コストを削減するために多くの研究が行われてきた。
合計領域テーブルを使用して3つのテーブルを作製する場合は、P_{ij}を合計し、i * P_{ij}を合計し、j * P_{ij}を合計し、その時のランタイムの複雑性は(O(N_pixels), O(N_bins*N_bins))の最大値である。しきい値に関して粗い解像度のみが必要な場合はN_binsを減らすことができることに注意してください。
Summed-area table参照
histsは、グレースケール値と近傍平均グレースケール値のペアの256\times 256の2次元ヒストグラムである。
totalは、所与の画像のペアの数である。これは各方向の2次元ヒストグラムのビンの数により決まる。
thresholdは取得されたしきい値である。
画像のクラスのグレーレベルは正規分布とみなすことができるが、サイズや分散が等しくない場合、大津のアルゴリズムの仮定は満たされない。Kittler-Illingworthのアルゴリズム(最小エラーしきい値処理とも呼ばれる)はこのような場合を処理するための大津の方法のバリエーションである。このアルゴリズムを数学的に説明する方法はいくつかある。それらの1つは、試験される各しきい値について結果の2値画像の正規分布のパラメータがデータが与えられた最尤推定により推定されたことを考慮することである。
このアルゴリズムは大津の方法よりも優れているように見えることがあるが、推定される新しいパラメータが導入されるためアルゴリズムが過剰パラメータ化されて不安定になる可能性がある。大津の方法からの仮定が少なくとも部分的に有効であると考えられる多くの場合において、オッカムの剃刀に従ってKittler-Illingworthのアルゴリズムよりも大津の方法を優先する方が好ましい場合がある。