Normalized Cuts and Image Segmentation

文章目录

Shi J. and Malik J. Normalized cuts and image segmentation. In IEEE Transactions on Pattern Analysis and Machine Intelligence.

在Digital Image Preprocessing的书上看到了这个算法, 对于其公式结果的推出不是很理解, 于是下载下来看了看. 本文主要讲的是一种利用图结构进行图像分割的算法.

主要内容

假设 f ( x , y ) , x = 1 , 2 , ⋯ M , y = 1 , 2 , ⋯ N f(x, y), x=1,2,\cdots M, y=1,2,\cdots N f(x,y),x=1,2,⋯M,y=1,2,⋯N为一张图片, 我们想要对其进行分割. 给定某一个距离函数, 可以用于衡量任意两点 i , j i, j i,j的相似度:
w i j = w ( i , j ) . w_{ij} = w(i, j). wij​=w(i,j).

把图片的每一个pixel看成一个节点, pixel和pixel之间的边为一条无向边, 则整体构成了一个无向的图 G = ( V , E ) G = (V, E) G=(V,E), 每条边的权重如上所述是 w i j w_{ij} wij​, 故易知 w i j = w j i w_{ij} = w_{ji} wij​=wji​. 我们的目标是将图分成相斥的两块 A , B A, B A,B, 即满足:
A ⋃ B = V , A ⋂ B = ∅ . A \bigcup B = V, A \bigcap B = \empty. A⋃B=V,A⋂B=∅.

以往的做法是, 找到一个分割, 使得下列指标最小:
c u t ( A , B ) = ∑ i ∈ A , j ∈ B w i j , cut(A, B) = \sum_{i \in A, j \in B} w_{ij}, cut(A,B)=i∈A,j∈B∑​wij​,
但是这种策略往往会导致不均匀的分割, 即最角落里的元素被单独分割出来:

Normalized Cuts and Image Segmentation

于是作者提出了一种新的指标:
N c u t ( A , B ) = c u t ( A , B ) a s s o c ( A , V ) + c u t ( A , B ) a s s o c ( B , V ) , Ncut(A, B) = \frac{cut(A, B)}{assoc(A, V)} + \frac{cut(A, B)}{assoc(B, V)}, Ncut(A,B)=assoc(A,V)cut(A,B)​+assoc(B,V)cut(A,B)​,
其中 a s s o c ( A , V ) = ∑ i ∈ A , j ∈ V w i j assoc(A, V) = \sum_{i \in A, j \in V} w_{ij} assoc(A,V)=∑i∈A,j∈V​wij​.
注意到:
N c u t ( A , B ) = c u t ( A , B ) c u t ( A , B ) + a s s o c ( A , A ) + c u t ( A , B ) c u t ( A , B ) + a s s o c ( B , B ) , Ncut(A, B) = \frac{cut(A, B)}{cut(A, B) + assoc(A, A)} + \frac{cut(A, B)}{cut(A, B) + assoc(B, B)}, Ncut(A,B)=cut(A,B)+assoc(A,A)cut(A,B)​+cut(A,B)+assoc(B,B)cut(A,B)​,
所以只有到 a s s o c ( A , A ) , a s s o c ( B , B ) assoc(A, A), assoc(B, B) assoc(A,A),assoc(B,B)都足够大的时候Ncut才会足够小, 这说明该指标更关注了内部的一种紧密性.

求解


x i = + 1 ,  if  i ∈ A , x i = − 1  if  i ∈ B d i = ∑ j w i j . x_i = +1, \text{ if } i \in A, \quad x_i = -1 \text{ if } i \in B \\ d_i = \sum_{j}w_{ij}. xi​=+1, if i∈A,xi​=−1 if i∈Bdi​=j∑​wij​.


N c u t ( A , B ) = ∑ x i > 0 , x j < 0 − w i j x i x j ∑ x i > 0 d i + ∑ x i < 0 , x i > 0 − w i j x i x j ∑ x i < 0 d i . Ncut(A, B) = \frac{\sum_{x_i > 0, x_j < 0} -w_{ij}x_i x_j}{\sum_{x_i > 0}d_i} +\frac{\sum_{x_i < 0, x_i > 0} -w_{ij}x_i x_j}{\sum_{x_i < 0}d_i}. Ncut(A,B)=∑xi​>0​di​∑xi​>0,xj​<0​−wij​xi​xj​​+∑xi​<0​di​∑xi​<0,xi​>0​−wij​xi​xj​​.

容易证明(但是不容易想到):
[ 1 + x 2 ] i = x i = + 1 ,   if  i ∈ A , [ 1 + x 2 ] i = 0 ,   if  i ∈ B . [\frac{1+x}{2}]_i = x_i = +1, \: \text{if } i \in A, \quad [\frac{1+x}{2}]_i = 0, \: \text{if } i \in B. [21+x​]i​=xi​=+1,if i∈A,[21+x​]i​=0,if i∈B.
[ 1 − x 2 ] i = − x i = + 1 ,   if  i ∈ B , [ 1 − x 2 ] i = 0 ,   if  i ∈ A . [\frac{1-x}{2}]_i = -x_i = +1, \: \text{if } i \in B, \quad [\frac{1-x}{2}]_i = 0, \: \text{if } i \in A. [21−x​]i​=−xi​=+1,if i∈B,[21−x​]i​=0,if i∈A.

[ W ] i j = w i j , D i i = d i , [W]_{ij} = w_ij, \\ D_{ii} = d_i, [W]ij​=wi​j,Dii​=di​,
且 D i i D_{ii} Dii​为对角矩阵.
所以我们能够证明以下事实:
4 ⋅ c u t ( A , B ) = ( 1 + x ) T W ( 1 − x ) 4 ⋅ a s s o c ( A , V ) = 2 ⋅ ( 1 + x ) T D 1 = ( 1 + x ) T D ( 1 + x ) 4 ⋅ a s s o c ( B , V ) = 2 ⋅ ( 1 − x ) T D 1 = ( 1 − x ) T D ( 1 − x ) a s s o c ( V , V ) = ∑ i d i = 1 T D 1 ( 1 + x ) T D ( 1 − x ) = 0. 4\cdot cut(A, B) = (1+x)^T W (1 - x) \\ 4 \cdot assoc(A, V) = 2\cdot (1 + x)^T D 1 = (1 + x)^T D (1 + x) \\ 4 \cdot assoc(B, V) = 2\cdot (1 - x)^T D 1 = (1 - x)^T D (1 - x) \\ assoc(V, V) = \sum_i d_i = 1^T D 1 \\ (1 + x)^T D (1 - x) = 0. 4⋅cut(A,B)=(1+x)TW(1−x)4⋅assoc(A,V)=2⋅(1+x)TD1=(1+x)TD(1+x)4⋅assoc(B,V)=2⋅(1−x)TD1=(1−x)TD(1−x)assoc(V,V)=i∑​di​=1TD1(1+x)TD(1−x)=0.

又注意到:
∑ x i > 0 , x j < 0 − w i j x i x j = ∑ x i > 0 [ d i − ∑ x j > 0 w i j ] = 1 4 ( 1 + x ) T ( D − W ) ( 1 + x ) , \sum_{x_i > 0, x_j < 0} -w_{ij} x_i x_j = \sum_{x_i > 0} [d_i - \sum_{x_j >0} w_{ij}] = \frac{1}{4}(1 + x)^T (D - W) (1 + x), xi​>0,xj​<0∑​−wij​xi​xj​=xi​>0∑​[di​−xj​>0∑​wij​]=41​(1+x)T(D−W)(1+x),
于是同理可证:
( 1 + x ) T W ( 1 − x ) = ( 1 + x ) T ( D − W ) ( 1 + x ) = ( 1 − x ) T ( D − W ) ( 1 − x ) . (1 + x)^TW(1 - x) = (1 + x)^T (D - W)(1 +x) = (1 - x)^T (D - W)(1 -x) . (1+x)TW(1−x)=(1+x)T(D−W)(1+x)=(1−x)T(D−W)(1−x).


k = a s s o c ( A , V ) a s s o c ( V , V ) , k = \frac{assoc(A, V)}{assoc(V, V)}, k=assoc(V,V)assoc(A,V)​,

1 − k = a s s o c ( B , V ) a s s o c ( V , V ) . 1 - k = \frac{assoc(B, V)}{assoc(V, V)}. 1−k=assoc(V,V)assoc(B,V)​.

综上可得:
N c u t ( A , B ) = c u t ( A , B ) k 1 T D 1 + c u t ( A , B ) ( 1 − k ) 1 T D 1 = c u t ( A , B ) k ( 1 − k ) 1 T D 1 . Ncut(A, B) = \frac{cut(A, B)}{k1^T D1} + \frac{cut(A, B)}{(1-k)1^TD1} = \frac{cut(A, B)}{k(1-k)1^TD1}. Ncut(A,B)=k1TD1cut(A,B)​+(1−k)1TD1cut(A,B)​=k(1−k)1TD1cut(A,B)​.


[ ( 1 + x ) − b ( 1 − x ) ] T ( D − W ) [ ( 1 + x ) − b ( 1 − x ) ] = ( 1 + x ) T ( D − W ) ( 1 + x ) + b 2 ( 1 − x ) T ( D − W ) − 2 b ( 1 + x ) T ( D − W ) ( 1 − x ) = 4 ( 1 + b 2 ) c u t ( A , B ) − 2 b ( 1 + x ) T D ( 1 − x ) + 2 b ( 1 + x ) T W ( 1 − x ) = 4 ( 1 + b 2 ) c u t ( A , B ) − 0 + 8 b c u t ( A , B ) = 4 ( 1 + b ) 2 c u t ( A , B ) . \begin{array}{ll} &[(1 + x) - b(1-x)]^T (D-W)[(1+x) - b(1-x)] \\ =& (1+x)^T(D-W)(1+x) + b^2 (1-x)^T(D-W) \\ &- 2b (1+x)^T(D-W)(1-x) \\ =&4(1+b^2)cut(A, B) - 2b (1 + x)^TD(1-x) + 2b(1 + x)^T W(1-x) \\ =&4(1+b^2)cut(A, B) - 0 + 8b cut(A, B) \\ =&4(1 + b)^2 cut(A, B). \end{array} ====​[(1+x)−b(1−x)]T(D−W)[(1+x)−b(1−x)](1+x)T(D−W)(1+x)+b2(1−x)T(D−W)−2b(1+x)T(D−W)(1−x)4(1+b2)cut(A,B)−2b(1+x)TD(1−x)+2b(1+x)TW(1−x)4(1+b2)cut(A,B)−0+8bcut(A,B)4(1+b)2cut(A,B).​

( 1 + k 1 − k ) 2 = 1 ( 1 − k ) 2 , (1 + \frac{k}{1-k})^2 = \frac{1}{(1-k)^2}, (1+1−kk​)2=(1−k)21​,

4 ⋅ N c u t ( A , B ) = 4 ( 1 + b ) 2 b 1 T D 1 = [ ( 1 + x ) − b ( 1 − x ) ] T ( D − W ) [ ( 1 + x ) − b ( 1 − x ) ] b 1 T D 1 , b = k 1 − k . 4\cdot Ncut(A,B) = \frac{4(1+b)^2}{b1^TD1} = \frac{[(1 + x) - b(1-x)]^T (D-W)[(1+x) - b(1-x)]}{b1^TD1}, \\ b = \frac{k}{1-k}. 4⋅Ncut(A,B)=b1TD14(1+b)2​=b1TD1[(1+x)−b(1−x)]T(D−W)[(1+x)−b(1−x)]​,b=1−kk​.
令 y = ( 1 + x ) − b ( 1 − x ) y = (1 + x) - b(1 - x) y=(1+x)−b(1−x), 且
y T D y = ∑ x i > 0 d i + b 2 ∑ x i < 0 d i = b ( ∑ x i < 0 d i + b ∑ x i < 0 d i ) = b 1 T D 1. y^TDy = \sum_{x_i > 0}d_i + b^2 \sum_{x_i < 0}d_i = b( \sum_{x_i < 0}d_i + b \sum_{x_i < 0}d_i) = b1^TD1. yTDy=xi​>0∑​di​+b2xi​<0∑​di​=b(xi​<0∑​di​+bxi​<0∑​di​)=b1TD1.
4 ⋅ N c u t ( A , B ) = y T ( D − W ) y y T D y . 4 \cdot Ncut(A, B) = \frac{y^T(D-W)y}{y^TDy}. 4⋅Ncut(A,B)=yTDyyT(D−W)y​.

min ⁡ x N c u t ( A , B ) = min ⁡ y 1 4 y T ( D − W ) y y T D y , s . t . y i ∈ { 1 , 1 − b } . \min_x Ncut(A, B) = \min_y \frac{1}{4} \frac{y^T(D-W)y}{y^TDy}, \\ \mathrm{s.t.} \quad y_i \in \{1, 1 - b\}. xmin​Ncut(A,B)=ymin​41​yTDyyT(D−W)y​,s.t.yi​∈{1,1−b}.
倘若我们能放松条件至实数域中, 此时只需要通过求解下列系统:
( D − W ) y = λ D y ⇔ D − 1 2 ( D − W ) D − 1 2 z = λ z , z = D 1 2 y . (D-W)y = \lambda Dy \Leftrightarrow D^{-\frac{1}{2}}(D-W)D^{-\frac{1}{2}} z = \lambda z, z = D^{\frac{1}{2}}y. (D−W)y=λDy⇔D−21​(D−W)D−21​z=λz,z=D21​y.
需要注意的是:
( D − W ) 1 = 0 , (D-W)1 = 0, (D−W)1=0,
此时 z 0 = D 1 2 1 z_0 = D^{\frac{1}{2}}1 z0​=D21​1,
故 1 1 1实际上上述系统的一个解, 且对应最小的特征值, 但其不是我们所要的解. 因为 y y y必须要还满足:
y T D 1 = ∑ x i > 0 d i − b ∑ x i < 0 d i = 0 , y^T D 1 = \sum_{x_i > 0}d_i - b \sum_{x_i < 0} d_i = 0, yTD1=xi​>0∑​di​−bxi​<0∑​di​=0,
这意味着, 我们要的恰恰是
D − 1 2 ( D − W ) D − 1 2 z = λ z , z = D 1 2 y D^{-\frac{1}{2}}(D-W)D^{-\frac{1}{2}} z = \lambda z, z = D^{\frac{1}{2}}y D−21​(D−W)D−21​z=λz,z=D21​y
倒数第二小的特征值对应的特征向量 z 1 z_1 z1​, 于是 y 1 = D − 1 2 z 1 y_1 = D^{-\frac{1}{2}}z_1 y1​=D−21​z1​.

相似度

文中采用如下的计算方式:

w i j = { e − ∥ F i − F j ∥ 2 / σ I 2 ⋅ e − ∥ X i − X j ∥ 2 / σ X 2 if  ∥ X i − X j ∥ < r 0 else . w_{ij} = \left \{ \begin{array}{ll} e^{-\|F_i - F_j\|^2 / \sigma^2_I} \cdot e^{-\|X_i - X_j\|^2 / \sigma^2_X} & \text{if } \|X_i - X_j\| < r \\ 0 & \text{else}. \end{array} \right. wij​={e−∥Fi​−Fj​∥2/σI2​⋅e−∥Xi​−Xj​∥2/σX2​0​if ∥Xi​−Xj​∥<relse.​
其中 F F F对应颜色之类的距离, 如直接取密度值, 而 X X X对应空间距离, r r r限定了搜索范围, 同样会导致 W W W变成系数矩阵, 对应特征求解加速有帮助.

总的算法流程

  1. 计算权重矩阵 W W W以及 D D D;
  2. 通过
    D − 1 2 ( D − W ) D − 1 2 z = λ z D^{-\frac{1}{2}}(D-W)D^{-\frac{1}{2}} z = \lambda z D−21​(D−W)D−21​z=λz
    计算得到倒数第二小的特征值所对应的特征向量 z 1 z_1 z1​并令 y 1 = z 1 y_1=z_1 y1​=z1​;
  3. 通过某种方法(如网格搜索)找到一个阈值 t t t:
    x i = 1 ,   if  y i > t ,   else  − 1. x_i = 1, \: \text{if }y_i > t, \: \text{else } -1. xi​=1,if yi​>t,else −1.
    且 x x x的划分下
    N c u t ( A , B ) Ncut(A, B) Ncut(A,B)
    较小.
  4. 对于 A , B A, B A,B可以重复上述分割过程, 直到满足区域数目或者其它某种条件(比如文中说的特征向量的分布过于均匀时停止).

skimage.future.graph.cut

skimage.future.graph.cut

上一篇:使用pd.cut进行分箱操作


下一篇:<数据库> if 条件语句的使用 SQL26 计算25岁以上和以下的用户数量