第六章 圖像聚類

6.1 K-Means聚類

6.1.1 SciPy聚類包

6.1.2 圖像聚類

6.1.3 在主成分上可視化圖像

6.1.4 像素聚類

6.2 層次聚類

6.2.1 圖像聚類

6.3 譜聚類

這一章會介紹幾種聚類方法，并就怎么使用它們對圖像進行聚類找出相似的圖像組進行說明。聚類可以用于識別，劃分圖像數據集、組織導航等。同時，我們也會用聚類相似的圖像進行可視化。

6.1 K-Means聚類

K-means是一種非常簡單的聚類算法，它能夠將輸入數據劃分成k個簇。關于K-means聚類算法的介紹可以參閱中譯本。

6.1.1 SciPy聚類包

盡管K-means聚類算法很容易實現，但我們沒必要自己去實現。SciPy矢量量化包sci.cluter.vq中有k-means的實現。這里我們演示怎樣使用它。

我們以2維示例樣本數據進行說明：

# coding=utf-8 """ Function: figure 6.1 An example of k-means clustering of 2D points """ from pylab import * from scipy.cluster.vq import *# 添加中文字體支持 from matplotlib.font_manager import FontProperties font = FontProperties(fname=r"c:\windows\fonts\SimSun.ttc", size=14)class1 = 1.5 * randn(100, 2) class2 = randn(100, 2) + array([5, 5]) features = vstack((class1, class2)) centroids, variance = kmeans(features, 2) code, distance = vq(features, centroids) figure() ndx = where(code == 0)[0] plot(features[ndx, 0], features[ndx, 1], '*') ndx = where(code == 1)[0] plot(features[ndx, 0], features[ndx, 1], 'r.') plot(centroids[:, 0], centroids[:, 1], 'go')title(u'2維數據點聚類', fontproperties=font) axis('off') show()

上面代碼中where()函數給出每類的索引。運行上面代碼，可得到原書P129頁圖6-1，即：

6.1.2 圖像聚類

現在我們用k-means對原書14頁的圖像進行聚類，文件selectedfontimages.zip包含了66張字體圖像。對于每一張圖像，我們用在前40個主成分上投影后的系數作為特征向量。下面為對其進行聚類的代碼：

# -*- coding: utf-8 -*- from PCV.tools import imtools import pickle from scipy import * from pylab import * from PIL import Image from scipy.cluster.vq import * from PCV.tools import pca# Uses sparse pca codepath. imlist = imtools.get_imlist('../data/selectedfontimages/a_selected_thumbs/')# 獲取圖像列表和他們的尺寸 im = array(Image.open(imlist[0])) # open one image to get the size m, n = im.shape[:2] # get the size of the images imnbr = len(imlist) # get the number of images print "The number of images is %d" % imnbr# Create matrix to store all flattened images immatrix = array([array(Image.open(imname)).flatten() for imname in imlist], 'f')# PCA降維 V, S, immean = pca.pca(immatrix)# 保存均值和主成分 #f = open('./a_pca_modes.pkl', 'wb') f = open('./a_pca_modes.pkl', 'wb') pickle.dump(immean,f) pickle.dump(V,f) f.close()# get list of images imlist = imtools.get_imlist('../data/selectedfontimages/a_selected_thumbs/') imnbr = len(imlist)# load model file with open('../data/selectedfontimages/a_pca_modes.pkl','rb') as f:immean = pickle.load(f)V = pickle.load(f) # create matrix to store all flattened images immatrix = array([array(Image.open(im)).flatten() for im in imlist],'f')# project on the 40 first PCs immean = immean.flatten() projected = array([dot(V[:40],immatrix[i]-immean) for i in range(imnbr)])# k-means projected = whiten(projected) centroids,distortion = kmeans(projected,4) code,distance = vq(projected,centroids)# plot clusters for k in range(4):ind = where(code==k)[0]figure()gray()for i in range(minimum(len(ind),40)):subplot(4,10,i+1)imshow(immatrix[ind[i]].reshape((25,25)))axis('off') show()

運行上面代碼，可得到下面的聚類結果：注：這里的結果譯者截的是原書上的結果，上面代碼實際運行出來的結果可能跟上面有出入。

6.1.3 在主成分上可視化圖像

# -*- coding: utf-8 -*- from PCV.tools import imtools, pca from PIL import Image, ImageDraw from pylab import * from PCV.clustering import hclusterimlist = imtools.get_imlist('../data/selectedfontimages/a_selected_thumbs') imnbr = len(imlist)# Load images, run PCA. immatrix = array([array(Image.open(im)).flatten() for im in imlist], 'f') V, S, immean = pca.pca(immatrix)# Project on 2 PCs. projected = array([dot(V[[0, 1]], immatrix[i] - immean) for i in range(imnbr)]) # P131 Fig6-3左圖 #projected = array([dot(V[[1, 2]], immatrix[i] - immean) for i in range(imnbr)]) # P131 Fig6-3右圖# height and width h, w = 1200, 1200# create a new image with a white background img = Image.new('RGB', (w, h), (255, 255, 255)) draw = ImageDraw.Draw(img)# draw axis draw.line((0, h/2, w, h/2), fill=(255, 0, 0)) draw.line((w/2, 0, w/2, h), fill=(255, 0, 0))# scale coordinates to fit scale = abs(projected).max(0) scaled = floor(array([(p/scale) * (w/2 - 20, h/2 - 20) + (w/2, h/2)for p in projected])).astype(int)# paste thumbnail of each image for i in range(imnbr):nodeim = Image.open(imlist[i])nodeim.thumbnail((25, 25))ns = nodeim.sizebox = (scaled[i][0] - ns[0] // 2, scaled[i][1] - ns[1] // 2,scaled[i][0] + ns[0] // 2 + 1, scaled[i][1] + ns[1] // 2 + 1)img.paste(nodeim, box)tree = hcluster.hcluster(projected) hcluster.draw_dendrogram(tree,imlist,filename='fonts.png')figure() imshow(img) axis('off') img.save('../images/ch06/pca_font.png') show()

運行上面代碼，可畫出原書P131圖6-3中的實例結果。

6.1.4 像素聚類

在結束這節前，我們看一個對像素進行聚類而不是對所有的圖像進行聚類的例子。將圖像區域歸并成“有意義的”組件稱為圖像分割。在第九章會將其單獨列為一個主題。在像素級水平進行聚類除了可以用在一些很簡單的圖像，在其他圖像上進行聚類是沒有意義的。這里，我們將k-means應用到RGB顏色值上，關于分割問題會在第九章第二節會給出分割的方法。下面是對兩幅圖像進行像素聚類的例子(注：譯者對原書中的代碼做了調整)：

# -*- coding: utf-8 -*- """ Function: figure 6.4 Clustering of pixels based on their color value using k-means. """ from scipy.cluster.vq import * from scipy.misc import imresize from pylab import * import Image# 添加中文字體支持 from matplotlib.font_manager import FontProperties font = FontProperties(fname=r"c:\windows\fonts\SimSun.ttc", size=14)def clusterpixels(infile, k, steps):im = array(Image.open(infile))dx = im.shape[0] / stepsdy = im.shape[1] / steps# compute color features for each regionfeatures = []for x in range(steps):for y in range(steps):R = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 0])G = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 1])B = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 2])features.append([R, G, B])features = array(features, 'f') # make into array# 聚類， k是聚類數目centroids, variance = kmeans(features, k)code, distance = vq(features, centroids)# create image with cluster labelscodeim = code.reshape(steps, steps)codeim = imresize(codeim, im.shape[:2], 'nearest')return codeimk=3 infile_empire = '../data/empire.jpg' im_empire = array(Image.open(infile_empire)) infile_boy_on_hill = '../data/boy_on_hill.jpg' im_boy_on_hill = array(Image.open(infile_boy_on_hill)) steps = (50, 100) # image is divided in steps*steps region print steps[0], steps[-1]#顯示原圖empire.jpg figure() subplot(231) title(u'原圖', fontproperties=font) axis('off') imshow(im_empire)# 用50*50的塊對empire.jpg的像素進行聚類 codeim= clusterpixels(infile_empire, k, steps[0]) subplot(232) title(u'k=3,steps=50', fontproperties=font) #ax1.set_title('Image') axis('off') imshow(codeim)# 用100*100的塊對empire.jpg的像素進行聚類 codeim= clusterpixels(infile_empire, k, steps[-1]) ax1 = subplot(233) title(u'k=3,steps=100', fontproperties=font) #ax1.set_title('Image') axis('off') imshow(codeim)#顯示原圖empire.jpg subplot(234) title(u'原圖', fontproperties=font) axis('off') imshow(im_boy_on_hill)# 用50*50的塊對empire.jpg的像素進行聚類 codeim= clusterpixels(infile_boy_on_hill, k, steps[0]) subplot(235) title(u'k=3,steps=50', fontproperties=font) #ax1.set_title('Image') axis('off') imshow(codeim)# 用100*100的塊對empire.jpg的像素進行聚類 codeim= clusterpixels(infile_boy_on_hill, k, steps[-1]) subplot(236) title(u'k=3，steps=100', fontproperties=font) axis('off') imshow(codeim)show()

上面代碼中，先載入一幅圖像，然后用一個steps*steps的方塊在原圖中滑動，對窗口中的圖像值求和取平均，將它下采樣到一個較低的分辨率，然后對這些區域用k-means進行聚類。運行上面代碼，即可得出原書P133頁圖6-4中的圖。

6.2 層次聚類

層次聚類(或稱凝聚聚類)是另一種簡單但有效的聚類算法。下面我們我們通過一個簡單的實例看看層次聚類是怎樣進行的。

from pylab import * from PCV.clustering import hclusterclass1 = 1.5 * randn(100,2) class2 = randn(100,2) + array([5,5]) features = vstack((class1,class2))tree = hcluster.hcluster(features) clusters = tree.extract_clusters(5) print 'number of clusters', len(clusters) for c in clusters:print c.get_cluster_elements()

上面代碼首先創建一些2維數據點，然后對這些數據點聚類，用一些閾值提取列表中的聚類后的簇群，并將它們打印出來，譯者在自己的筆記本上打印出的結果為：

number of clusters 2 [197, 107, 176, 123, 173, 189, 154, 136, 183, 113, 109, 199, 178, 129, 163, 100, 148, 111, 143, 118, 162, 169, 138, 182, 193, 116, 134, 198, 184, 181, 131, 166, 127, 185, 161, 171, 152, 157, 112, 186, 128, 156, 108, 158, 120, 174, 102, 137, 117, 194, 159, 105, 155, 132, 188, 125, 180, 151, 192, 164, 195, 126, 103, 196, 179, 146, 147, 135, 139, 110, 140, 106, 104, 115, 149, 190, 170, 172, 121, 145, 114, 150, 119, 142, 122, 144, 160, 187, 153, 167, 130, 133, 165, 191, 175, 177, 101, 141, 124, 168] [0, 39, 32, 87, 40, 48, 28, 8, 26, 12, 94, 5, 1, 61, 24, 59, 83, 10, 99, 50, 23, 58, 51, 16, 71, 25, 11, 37, 22, 46, 60, 86, 65, 2, 21, 4, 41, 72, 80, 84, 33, 56, 75, 77, 29, 85, 93, 7, 73, 6, 82, 36, 49, 98, 79, 43, 91, 14, 47, 63, 3, 97, 35, 18, 44, 30, 13, 67, 62, 20, 57, 89, 88, 9, 54, 19, 15, 92, 38, 64, 45, 70, 52, 95, 69, 96, 42, 53, 27, 66, 90, 81, 31, 34, 74, 76, 17, 78, 55, 68]

6.2.1 圖像聚類

# -*- coding: utf-8 -*- import os import Image from PCV.clustering import hcluster from matplotlib.pyplot import * from numpy import *# create a list of images path = '../data/sunsets/flickr-sunsets-small/' imlist = [os.path.join(path, f) for f in os.listdir(path) if f.endswith('.jpg')] # extract feature vector (8 bins per color channel) features = zeros([len(imlist), 512]) for i, f in enumerate(imlist):im = array(Image.open(f))# multi-dimensional histogramh, edges = histogramdd(im.reshape(-1, 3), 8, normed=True, range=[(0, 255), (0, 255), (0, 255)])features[i] = h.flatten() tree = hcluster.hcluster(features)# visualize clusters with some (arbitrary) threshold clusters = tree.extract_clusters(0.23 * tree.distance) # plot images for clusters with more than 3 elements for c in clusters:elements = c.get_cluster_elements()nbr_elements = len(elements)if nbr_elements > 3:figure()for p in range(minimum(nbr_elements,20)):subplot(4, 5, p + 1)im = array(Image.open(imlist[elements[p]]))imshow(im)axis('off') show()hcluster.draw_dendrogram(tree,imlist,filename='sunset.pdf')

運行上面代碼，可得原書P140圖6-6。同時會在上面腳本文件所在的文件夾下生成層次聚類后的簇群樹：我們對前面字體圖像同樣創建一個樹，正如前面在主成分可視化圖像中，我們添加了下面代碼：

tree = hcluster.hcluster(projected) hcluster.draw_dendrogram(tree,imlist,filename='fonts.png')

運行添加上面兩行代碼后前面的例子，可得對字體進行層次聚類后的簇群樹：

6.3 譜聚類

譜聚類是另一種不同于k-means和層次聚類的聚類算法。關于譜聚類的原理，可以參閱中譯本。這里，我們用原來k-means實例中用到的字體圖像。

# -*- coding: utf-8 -*- from PCV.tools import imtools, pca from PIL import Image, ImageDraw from pylab import * from scipy.cluster.vq import *imlist = imtools.get_imlist('../data/selectedfontimages/a_selected_thumbs') imnbr = len(imlist)# Load images, run PCA. immatrix = array([array(Image.open(im)).flatten() for im in imlist], 'f') V, S, immean = pca.pca(immatrix)# Project on 2 PCs. projected = array([dot(V[[0, 1]], immatrix[i] - immean) for i in range(imnbr)]) # P131 Fig6-3左圖 #projected = array([dot(V[[1, 2]], immatrix[i] - immean) for i in range(imnbr)]) # P131 Fig6-3右圖n = len(projected) # compute distance matrix S = array([[ sqrt(sum((projected[i]-projected[j])**2)) for i in range(n) ] for j in range(n)], 'f') # create Laplacian matrix rowsum = sum(S,axis=0) D = diag(1 / sqrt(rowsum)) I = identity(n) L = I - dot(D,dot(S,D)) # compute eigenvectors of L U,sigma,V = linalg.svd(L) k = 5 # create feature vector from k first eigenvectors # by stacking eigenvectors as columns features = array(V[:k]).T # k-means features = whiten(features) centroids,distortion = kmeans(features,k) code,distance = vq(features,centroids) # plot clusters for c in range(k):ind = where(code==c)[0]figure()gray()for i in range(minimum(len(ind),39)):im = Image.open(imlist[ind[i]])subplot(4,10,i+1)imshow(array(im))axis('equal')axis('off') show()

上面我們在前個特征向量上計算標準的k-means。下面是運行上面代碼的結果：注意，由于在k-means階段會給出不同的聚類結果，所以你運行上面代碼出來的結果可能跟譯者的是不一樣的。

同樣，我們可以在不知道特征向量或是沒有嚴格相似性定義的情況下進行譜聚類。原書44頁的位置地理圖像是通過它們之間有多少局部描述子匹配相連接的。48頁的相似性矩陣中的元素是為規范化的匹配特征點數。我們同樣可以對其進行譜聚類，完整的代碼如下：

# -*- coding: utf-8 -*- from PCV.tools import imtools, pca from PIL import Image, ImageDraw from PCV.localdescriptors import sift from pylab import * import glob from scipy.cluster.vq import *#download_path = "panoimages" # set this to the path where you downloaded the panoramio images #path = "/FULLPATH/panoimages/" # path to save thumbnails (pydot needs the full system path)download_path = "F:/dropbox/Dropbox/translation/pcv-notebook/data/panoimages" # set this to the path where you downloaded the panoramio images path = "F:/dropbox/Dropbox/translation/pcv-notebook/data/panoimages/" # path to save thumbnails (pydot needs the full system path)# list of downloaded filenames imlist = imtools.get_imlist('../data/panoimages/') nbr_images = len(imlist)# extract features #featlist = [imname[:-3] + 'sift' for imname in imlist] #for i, imname in enumerate(imlist): # sift.process_image(imname, featlist[i])featlist = glob.glob('../data/panoimages/*.sift')matchscores = zeros((nbr_images, nbr_images))for i in range(nbr_images):for j in range(i, nbr_images): # only compute upper triangleprint 'comparing ', imlist[i], imlist[j]l1, d1 = sift.read_features_from_file(featlist[i])l2, d2 = sift.read_features_from_file(featlist[j])matches = sift.match_twosided(d1, d2)nbr_matches = sum(matches > 0)print 'number of matches = ', nbr_matchesmatchscores[i, j] = nbr_matches print "The match scores is: \n", matchscores# copy values for i in range(nbr_images):for j in range(i + 1, nbr_images): # no need to copy diagonalmatchscores[j, i] = matchscores[i, j]n = len(imlist) # load the similarity matrix and reformat S = matchscores S = 1 / (S + 1e-6) # create Laplacian matrix rowsum = sum(S,axis=0) D = diag(1 / sqrt(rowsum)) I = identity(n) L = I - dot(D,dot(S,D)) # compute eigenvectors of L U,sigma,V = linalg.svd(L) k = 2 # create feature vector from k first eigenvectors # by stacking eigenvectors as columns features = array(V[:k]).T # k-means features = whiten(features) centroids,distortion = kmeans(features,k) code,distance = vq(features,centroids) # plot clusters for c in range(k):ind = where(code==c)[0]figure()gray()for i in range(minimum(len(ind),39)):im = Image.open(imlist[ind[i]])subplot(5,4,i+1)imshow(array(im))axis('equal')axis('off') show()

改變聚類數目k，可以得到不同的結果。譯者分別測試了原書中k=2和k=10的情況，運行結果如下：?k=2k=10注：對于聚類后，圖像小于或等于1的類，在上面沒有顯示。

from:?http://yongyuan.name/pcvwithpython/chapter6.html

總結

以上是生活随笔為你收集整理的Python计算机视觉：第六章 图像聚类的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。