[Scikit-learn] 2.3 Clustering - Spectral clustering

From: 2.3.5 Clustering - Spectral clustering

From: 漫談 Clustering (4): Spectral Clustering

From: 漫談 Clustering (番外篇): Dimensionality Reduction

---

傳統方式的弊端

事實上，Laplacian Eigenmap （拉普拉斯特徵映射） 假設數據分佈在一個嵌套在高維空間中的低維流形上， Laplacian Matrix  則是流形的 Laplace Beltrami operator 的一個離散近似。關於流形以及 Laplacian Eigenmap 這個模型的理論知識就不在這裏作更多地展開了，下面看一個比較直觀的例子 Swiss Roll 。

Swiss Roll 是一個像麪包圈同樣的結構，能夠看做一個嵌套在三維空間中的二維流形，以下圖所示，左邊是 Swiss Roll ，右邊是從 Swiss Roll 中隨機選出來的一些點，爲了標明點在流形上的位置，給它們都標上了顏色。

【顯然，這個東西無法用傳統方式作聚類分析，直觀的看上去，除了距離，還應該考慮結構因素】

而下圖是 Laplacian Eigenmap 和 PCA 分別對 Swiss Roll 降維的結果，能夠看到 LE 成功地把這個流形展開把在流形上屬於不一樣位置的點映射到了不一樣的地方，而 PCA 的結果則很糟糕，幾種顏色都混雜在了一塊兒。

 

因此，Spectral Clustering便登場了！

能夠看到，只有三個比較人性化：SpectralClustering, AgglomerativeClustering, DBSCAN

 

SpectralClustering

1. 根據數據構造一個 Graph ，Graph 的每個節點對應一個數據點，將類似的點鏈接起來，而且邊的權重用於表示數據之間的類似度。把這個 Graph 用鄰接矩陣的形式表示出來，記爲  。一個最偷懶的辦法就是：直接用咱們前面在 K-medoids 中用的類似度矩陣做爲  。
2. 把  的每一列元素加起來獲得  個數，把它們放在對角線上（其餘地方都是零），組成一個  的矩陣，記爲  。並令  。
3. 求出  的前  個特徵值（在本文中，除非特殊說明，不然「前  個」指按照特徵值的大小從小到大的順序） 以及對應的特徵向量  。
4. 把這  個特徵（列）向量排列在一塊兒組成一個  的矩陣，將其中每一行看做  維空間中的一個向量，並使用 K-means 算法進行聚類。聚類的結果中每一行所屬的類別就是原來 Graph 中的節點亦即最初的  個數據點分別所屬的類別。

 

其實，若是你熟悉 Dimensional Reduction （降維）的話，大概已經看出來了，Spectral Clustering 其實就是：

　　經過 Laplacian Eigenmap 的降維方式降維以後再作 K-means 的一個過程。

不過，爲何要恰好降到  維呢？其實整個模型還能夠從另外一個角度導出來，因此，讓咱們不妨先跑一下題。

在 Image Processing （我好像以前有據說我對這個領域深惡痛絕？）裏有一個問題就是對圖像進行 Segmentation （區域分割），也就是讓類似的像素組成一個區域，好比，咱們通常但願一張照片裏面的人（前景）和背景被分割到不一樣的區域中。在 Image Processing 領域裏已經有許多自動或半自動的算法來解決這個問題，而且有很多方法和 Clustering 有密切聯繫。好比咱們在談 Vector Quantization 的時候就曾經用 K-means 來把顏色類似的像素聚類到一塊兒，不過那還不是真正的 Segmentation ，由於若是僅僅是考慮顏色類似的話，圖片上位置離得很遠的像素也有可能被聚到同一類中，咱們一般並不會把這樣一些「遊離」的像素構成的東西稱爲一個「區域」，但這個問題其實也很好解決：只要在聚類用的 feature 中加入位置信息（例如，原來是使用 R、G、B 三個值來表示一個像素，如今加入 x、y 兩個新的值）便可。

 

Ref: http://blog.csdn.net/abcjennifer/article/details/8170687

將邊權賦值爲兩點之間的similarity，作聚類的目標就是最小化類間connection的weight。

好比對於下面這幅圖，分割以下

可是這樣有可能會有問題，好比：

因爲Graph cut criteria 只考慮了類間差小，而沒考慮internal cluster density.因此會有上面分割的問題。這裏引入Normalised-cut（Shi & Malik, 97'）。

 

Graph Cut

另外一方面，還有一個常常被研究的問題就是 Graph Cut ，簡單地說就是把一個 Graph 的一些邊切斷，讓他被打散成一些獨立聯通的 sub-Graph ，而這些被切斷的邊的權值的總和就被稱爲 Cut 值。

若是用一張圖片中的全部像素來組成一個 Graph ，並把（好比，顏色和位置上）類似的節點鏈接起來，邊上的權值表示類似程度，那麼把圖片分割爲幾個區域的問題實際上等價於把 Graph 分割爲幾個 sub-Graph 的問題，而且咱們能夠要求分割所得的 Cut 值最小，亦即：那些被切斷的邊的權值之和最小，直觀上咱們能夠知道，權重比較大的邊沒有被切斷，表示比較類似的點被保留在了同一個 sub-Graph 中，而彼此之間聯繫不大的點則被分割開來。咱們能夠認爲這樣一種分割方式是比較好的。

實際上，拋開圖像分割的問題不談，在 Graph Cut 相關的一系列問題中，Minimum cut （最小割）自己就是一個被普遍研究的問題，而且有成熟的算法來求解。只是單純的最小割在這裏一般並非特別適用，不少時候只是簡單地把和其餘像素聯繫最弱的那一個像素給分割出去了，相反，咱們一般更但願分割出來的區域（的大小）要相對均勻一些，而不是一些很大的區塊和一些幾乎是孤立的點。爲此，又有許多替代的算法提出來，如 Ratio Cut 、Normalized Cut 等。

 

Laplacian Eigenmap

From: http://www.cnblogs.com/lmsj918/p/4075228.html

主要內容是Laplacian Eigenmap中的核心推導過程。

有空仍是多點向這位師兄請教，每次都會撿到很多金子。

Reference : 《Laplacian Eigenmaps for Dimensionality Reduction and Data Representation》,2003,MIT

 

From: Spectral clustering for image segmentation

我很是喜歡的一個示例代碼

print(__doc__) # Authors: Emmanuelle Gouillart <emmanuelle.gouillart@normalesup.org> # Gael Varoquaux <gael.varoquaux@normalesup.org> # License: BSD 3 clause

import numpy as np import matplotlib.pyplot as plt from sklearn.feature_extraction import image from sklearn.cluster import spectral_clustering l = 100 x, y = np.indices((l, l)) center1 = (28, 24) center2 = (40, 50) center3 = (67, 58) center4 = (24, 70) radius1, radius2, radius3, radius4 = 16, 14, 15, 14 circle1 = (x - center1[0]) ** 2 + (y - center1[1]) ** 2 < radius1 ** 2 circle2 = (x - center2[0]) ** 2 + (y - center2[1]) ** 2 < radius2 ** 2 circle3 = (x - center3[0]) ** 2 + (y - center3[1]) ** 2 < radius3 ** 2 circle4 = (x - center4[0]) ** 2 + (y - center4[1]) ** 2 < radius4 ** 2

# ############################################################################# # 4 circles
img = circle1 + circle2 + circle3 + circle4 # We use a mask that limits to the foreground: the problem that we are # interested in here is not separating the objects from the background, # but separating them one from the other.
mask = img.astype(bool) img = img.astype(float) img += 1 + 0.2 * np.random.randn(*img.shape) # Convert the image into a graph with the value of the gradient on the # edges. 攜帶了梯度值，有什麼用麼？
graph = image.img_to_graph(img, mask=mask) # Take a decreasing function of the gradient: we take it weakly # dependent from the gradient the segmentation is close to a voronoi
graph.data = np.exp(-graph.data / graph.data.std()) # Force the solver to be arpack, since amg is numerically # unstable on this example
labels = spectral_clustering(graph, n_clusters=4, eigen_solver='arpack') label_im = -np.ones(mask.shape) label_im[mask] = labels plt.matshow(img) plt.matshow(label_im) # ############################################################################# # 2 circles
img = circle1 + circle2 mask = img.astype(bool) img = img.astype(float) img += 1 + 0.2 * np.random.randn(*img.shape) graph = image.img_to_graph(img, mask=mask) graph.data = np.exp(-graph.data / graph.data.std()) labels = spectral_clustering(graph, n_clusters=2, eigen_solver='arpack') label_im = -np.ones(mask.shape) label_im[mask] = labels plt.matshow(img) plt.matshow(label_im) plt.show()

Result:

 

算法理解：

# Convert the image into a graph with the value of the gradient on the # edges. 攜帶了梯度值，有什麼用麼？
graph = image.img_to_graph(img, mask=mask) # Take a decreasing function of the gradient: we take it weakly # dependent from the gradient the segmentation is close to a voronoi
graph.data = np.exp(-graph.data / graph.data.std()) # Force the solver to be arpack, since amg is numerically # unstable on this example
labels = spectral_clustering(graph, n_clusters=4, eigen_solver='arpack') label_im = -np.ones(mask.shape) label_im[mask] = labels

 

泰森多邊形 voronoi

 

 

參考：Dimensionality Reduction: Why we take Eigenvectors of the Similarity Matrix?